今天介紹進階的觀測資料格式-Rainbow5
雷達觀測原始資料儲存常使用Rainbow5格式，除了Rainbow5外，也有其他資料儲存格是可以儲存雷達資料，如netcdf及NEXRAD等。而原始資料代表的意思是資料皆未處理過，如果要得到比較可信的資料品質，必須要解出來之後進行資料後處理。
氣象局開放資料提供的雷達回波資料(進第七天)，就是經由雷達原始資料處理後再以json或xml資料格式提供給民眾使用。

通常資料檔名會是 yyyymmddHHMMSS00變數.vol或yyyymmddHHMMSS00變數.azi，以我現在要示範的檔名為例，檔名為2022052607405300dBZ.vol，也就是說這一份資料為2022年5月26日7點40分53秒的dBZ(回波)觀測值。那結尾的vol及azi分別代表的是整體資料及單層資料。雷達會依照不同仰角進行360度掃描，所以每一個仰角都會有一份360度的資料，即單層資料，而把每一份單層整合起來即即整體資料。
那原始資料長什麼樣子呢？讓我們先用瀏覽器打開此檔案看一下。

上圖是資料的一開始，沒錯，就是xml檔案，也就是說這一份rainbow5資料的部分資料是由xml檔案格式儲存。

再繼續看

上圖是顯示xml(END XML)檔案格式後的資料狀態，完全看不懂(XD)，因為是以binary的方式儲存。而這些binary就是實際的回波「值」囉，而xml檔案格式部分是提供一些雷達及雷達掃描策略的基本資訊，譬如每一條掃描波束的單位步及長度範圍。

那現在若要取得雷達原始資料必須向氣象局申購，氣象局的開放資料是沒有提供的。

現在開始讀取rainbow5的資料吧。
主要使用的有pyart 及 wradlib，這兩個套件都是處理雷達資料的python套件。
如果要解rainbow5資料的話，建議使用wradlib，那pyart可以取其提供的colorbar使用。
先安裝套件吧
conda install -c conda-forge arm_pyart wradlib
或
pip3 install arm_pyart wradlib

import wradlib as wrl
rvolread = wrl.io.read_rainbow("day24/2022052607405300dBZ.vol")

#用xml的方式解資料吧
'''
第一個標記是'volume'，volume下有6個標記，分別為@version', '@datetime', '@type', '@owner', 'scan', 'sensorinfo'，掃描相關的資料全數在scan標記。
scan標記下有10個標記，分別為'@datetimehighaccuracy', '@name', '@time', '@date', 'unitid', 'advancedchanged', 'detailedchanged', 'scantime', 'pargroup', 'slice'，而要取得的變數(這一份示範資料是雷達回波)，就在slice標記。

取得slice標記後，會得到一個list物件，此list內每1個元素即代表一個仰角的360度掃描資訊。那這一份有20個仰角，所以list的長度為20。
'''

elevation_angle = len(rvolread["volume"]["scan"]['slice'])
print(elevation_angle)  #會得到20

可以取得你想要的仰角層有關資料，譬如

ele11 = rvolread["volume"]["scan"]["slice"][12]
#取回波值所有資訊
dbzall = ele11["slicedata"]["rawdata"]
'''
用dbzall.keys()可以得到8個key，分別為
'@blobid', '@rays', '@type', '@bins', 
'@min', '@max', '@depth', 'data'

data顧名思義就是原始值，但這個原始值必須要經過轉換才會是真的回波值，轉換方式如下
realdata = min + data *(max - min) /2 ** depth
'''
rawd = dbzall["data"]
mind,maxd,depthd = float(dbzall["@min"]), float(dbzall["@max"]),\
                   float(dbzall["@depth"])
realdata = mind + rawd*(maxd-mind)/2**depthd

上述就可以得到原始的值囉

那取得值，值的維度格點是掃描角度乘以波束掃瞄範圍(徑向範圍)，以這一份資料來說是360乘以1600。
波束掃描範圍如下

'''
若一份資料有好幾個雷達掃描仰角資料，那波束掃瞄範圍會存在最低仰角的資料訊息中。
rangestep是單位掃描距離，stoprange是停止掃描的最大範圍，
也就是rangestep*n會等於stoprange
'''
rangstep = float(rvolread["volume"]["scan"]["slice"][0]["rangestep"])
stoprange = float(rvolread["volume"]["scan"]["slice"][0]["stoprange"])
r = np.arange(0,stoprange,rangstep)

現在用直角坐標系視覺化

from pathlib import Path
import matplotlib.pyplot as plt
import pyart

fpath = Path("font/msjh.ttf")
anglerange = np.arange(1,361) #設定360度
rs, angles = np.meshgrid(r,anglerange) #轉換網格角度對應波束掃瞄距離

cmap = pyart.graph.cm_colorblind.HomeyerRainbow #借用pyart的colormap
pcl = plt.pcolormesh(rs,angles,realdata,cmap=cmap)
cbar = plt.colorbar(pcl)
cbar.set_label("雷達回波",rotation=270,font=fpath)

我們可以看到在這個仰角下，當雷達掃描到方位角大約40至60度時(北方為零度，順時針開始計算角度)，完全沒有雷達回波資訊，代表有遮蔽物。這顆雷達是在南屯，在仰角6的狀況下呈現的回波。

那如果把仰角14.6視覺化

就可以看到遮蔽效果不見了。
不過現在用的資料都是未經過資料控管，看到的視覺化結果很像是真的量測到回波，但必須要進一步確認，這就是雷達專家們的工作囉。

明天應該是說明在地圖上繪製雷達回波，畢竟有誰會想看直角坐標系呢XD