謝千里,李百鋒
(中國民用航空中南地區(qū)空中交通管理局,廣東 廣州510403)
在我國經(jīng)濟(jì)的高速發(fā)展的同時(shí),世界氣候也在產(chǎn)生非常規(guī)的變化,極端天氣頻繁出現(xiàn),對社會(huì)和人民的生活和經(jīng)濟(jì)活動(dòng)產(chǎn)生嚴(yán)重的影響,其中對民航飛行的安全也構(gòu)成嚴(yán)重威脅。對天氣的監(jiān)測可謂多管齊下,其中天氣雷達(dá)有不可或缺的重要作用。隨著我國雷達(dá)網(wǎng)的建設(shè),各地基本得到覆蓋,而簡單的水平拼圖已不能滿足日益繁忙的民航空管的多高度層,垂直剖面天氣狀況的需求,在研究與開發(fā)相關(guān)軟件的過程中,提出并實(shí)現(xiàn)了一種天氣雷達(dá)柵格化三維拼圖算法,以解決獲取任意高度層和任意水平方向垂直剖面天氣回波圖的需求,對航空氣象提供了更的便利性,對通用氣象相關(guān)領(lǐng)域也有一定的借鑒意義。
天氣雷達(dá)體掃范圍是由十個(gè)左右的不同仰角掃描出來的離散圓椎面組成的一個(gè)很扁的底面為球面的倒圓椎體,如圖1所示。
圖1 雷達(dá)體掃數(shù)據(jù)范圍
進(jìn)行多雷達(dá)拼圖,各部雷達(dá)的數(shù)據(jù)的掃描角分辨率參數(shù)往往不一樣,例如有的雷達(dá)掃描一圈徑向角度有360個(gè),有的有512個(gè),有的是600個(gè),故需要先進(jìn)行角度的歸一化,例如,歸一化為360個(gè),即每度一個(gè)徑向數(shù)據(jù)。
把雷達(dá)掃描到的錐體裝進(jìn)一個(gè)立方體中,該立方體底部為海拔0 m,高為統(tǒng)一固定值,例如15 000 m,立方體由細(xì)小的柵格組成,每個(gè)柵格是一個(gè)小正方體,邊長為統(tǒng)一固定值,例如300 m,柵格是變換后回波的最小尺寸,如圖2所示。
體掃數(shù)據(jù)裝進(jìn)立體柵格化過程中存在很多空隙,在這些空隙處需要進(jìn)行插值,有分水平方向插值和垂直方向插值。
處在相鄰兩個(gè)方位中的某一柵格點(diǎn)C,計(jì)算出它與方位m的夾角α,與方位m+1的夾角β.如果α≤ β,則C點(diǎn)在徑向m上取值;如果α>β,則C點(diǎn)在徑向m+1上取值。見圖3.
圖3 水平方向徑向定位
氣象雷達(dá)在垂直方向上,掃描仰角只有10個(gè)左右,為了獲得盡可能大的垂直分辨率,垂直方向的插值算法如圖4所示。
圖4 垂直方向柵格插值
通過水平方向徑向定位,已將C點(diǎn)的取值確定在某一個(gè)方位的徑向上,假定為徑向m.
氣象雷達(dá)在垂直方向上,掃描仰角只有10個(gè)左右,每個(gè)仰角層都在同一方位上存在徑向m,為了獲得盡可能大的垂直分辨率,垂直方向的插值算法如圖4所示。
處在相鄰兩個(gè)仰角夾角中的某一柵格點(diǎn)C,以雷達(dá)天線點(diǎn)O為圓心,以O(shè)C為半徑畫弧,與仰角n的徑向m相交于點(diǎn)A,與仰角n+1的徑向m相交于點(diǎn)B,計(jì)算出夾角α,夾角β.
情況1:C點(diǎn)處于最小仰角與最大仰角之間。
如果α<β,則C點(diǎn)回波強(qiáng)度值=徑向1中A點(diǎn)的值;
如果α>β,則C點(diǎn)回波強(qiáng)度值=徑向2中B點(diǎn)的值。
如果α=β,則C點(diǎn)回波強(qiáng)度值=(徑向1中A點(diǎn)的值+徑向2中B點(diǎn)的值)/2.
情況2:C點(diǎn)處于最小仰角之下,則OA為水平線,射線OB為最小仰角徑向。
如果β≤ 0.5°,則C點(diǎn)回波強(qiáng)度值=徑向2中B點(diǎn)的值;
如果β>0.5°,則C點(diǎn)回波強(qiáng)度值=0.
情況3:C點(diǎn)處于最大仰角之上,則射線OA為最大仰角徑向,OB為垂直于地面的垂線。
如果α≤ 0.5°,則C點(diǎn)回波強(qiáng)度值=徑向1中A點(diǎn)的值;
如果α>0.5°,則C點(diǎn)回波強(qiáng)度值=0.
在氣象雷達(dá)拼圖應(yīng)用中,單部雷達(dá)覆蓋范圍較小,半徑為150 km~460 km,分辨率為百米級別,故單部雷達(dá)覆蓋范圍內(nèi)的各點(diǎn)位置誤差受地球不規(guī)則形狀的影響不大,所以在計(jì)算單部雷達(dá)回波位置時(shí),可以把地球看作標(biāo)準(zhǔn)球形。
如圖5所示,O點(diǎn)為雷達(dá)天線所在位置,h0為雷達(dá)天線所處的海拔高度,A點(diǎn)為探測范圍空中某點(diǎn),B點(diǎn)為A點(diǎn)與地心的連線和海平面的交點(diǎn)。r為A點(diǎn)距離雷達(dá)天線的直線距離,R為地球半徑,線段AB即為A點(diǎn)的海拔高度。
圖5 柵格海拔高度計(jì)算
即A點(diǎn)所處位置的海拔高度AB由兩個(gè)部分組成:
AB=h1+h2
上式中,h1為A點(diǎn)距雷達(dá)天線所處的海平面延伸平面的垂直高度;h2為因地球曲率而增加的高度。
A點(diǎn)距雷達(dá)站的水平距離為弧BC的長度B(C.
由幾何關(guān)系推導(dǎo)出
由以上4式可以計(jì)算出某點(diǎn)的海拔高度和離雷達(dá)站的水平距離,通過該算法進(jìn)行三維柵格化。
各地的雷達(dá)站分布在地球曲面上,而本算法拼圖的結(jié)果是一個(gè)平整的立方體,每一層都是平面,故需要將各個(gè)雷達(dá)站從橢球表面投影到平面上來,這與地圖投影的目的和使用的方法是一樣的。目前使用最廣泛的地圖投影算法有高斯-克呂格投影、蘭勃特投影、墨卡托投影。氣象行業(yè)使用最多的是蘭勃特投影中的等角割圓錐投影。
圖6 三維拼圖
經(jīng)過前面的步驟,單路雷達(dá)的柵格立方體分辨率和海拔高度已全部統(tǒng)一,故拼圖只需按投影后的相對位置拼在一起,形成一個(gè)大的柵格立方體,如上圖6所示。重疊部分的處理可根據(jù)需要采用不同的算法,例如常用的取最大值法、取平均值法、加權(quán)法[1]等等。
拼圖柵格立方體以3維數(shù)組形式存放于計(jì)算機(jī)內(nèi)存中,而展示給用戶的,是其所關(guān)注的平面產(chǎn)品,可根據(jù)需要在此立方體中計(jì)算或切割出來并生成圖片呈現(xiàn)給用戶,在此基礎(chǔ)上輸出產(chǎn)品只需要極其簡單的運(yùn)算即可。例如通過簡單的比較運(yùn)算,可生成組合反射率圖,水平切割立方體,得到對應(yīng)高度層的CAPPI圖,任意方向垂直切割立方體,得到該處的垂直剖面圖。
本算法在實(shí)際應(yīng)用中的效果如以下圖7~10所示。無論在水平方向的組合反射率、CAPPI,還是垂直方向任意角度和距離剖面,拼圖均達(dá)到了良好的視覺效果。
圖7 西南4省25部氣象雷達(dá)拼圖組合反射率
圖8 云南省9部氣象雷達(dá)拼圖組合反射率
圖9 云南省氣象局拼圖海拔4500米CAPPI
圖10 任意垂直剖面
本算法簡潔高效,所獲得的回波立體柵格,用途廣泛,可以根據(jù)需要生成各種類型的產(chǎn)品,或者進(jìn)行數(shù)據(jù)統(tǒng)計(jì)分析,或特征研究,對于氣象行業(yè)有重要意義。