基于保序回歸的X波段天氣雷達(dá)衰減訂正算法

郭春輝，袁微，王旭，白愛(ài)娟

（1.肇慶市氣象局，廣東肇慶 526040；2.成都信息工程大學(xué)，四川成都 610225）

X波段雷達(dá)在探測(cè)降水過(guò)程中能提供更高的時(shí)間和空間分辨率，而且具有天線體積小、易于移動(dòng)等優(yōu)點(diǎn)，并在定量估測(cè)降水和粒子的相態(tài)識(shí)別中，有著非常明顯的優(yōu)勢(shì)［1-3］。但由于X波段波長(zhǎng)較短，回波信號(hào)衰減問(wèn)題比起S波段和C波段更為嚴(yán)重，所以有必要對(duì)其回波進(jìn)行衰減訂正。Testud等［4］提出了降水廓線訂正算法（ZPHI）利用差分傳播相移（φDP）進(jìn)行訂正，ZPHI算法是分段進(jìn)行衰減訂正的，但是ZH-Kdp（ZH為水平反射率、Kdp為差分傳播相移率）訂正法整段仍然使用固定的系數(shù)；伍魏等［5］針對(duì)X波段雙偏振雷達(dá)資料先進(jìn)行φDP退折疊后重算Kdp，采用ZH-Kdp訂正法與自適應(yīng)約束算法對(duì)反射率進(jìn)行衰減訂正，結(jié)果表明ZH-Kdp訂正法效果更好；王晗等［6］采用ZH-Kdp和多點(diǎn)ZPHI的方法多反射率進(jìn)行衰減訂正，結(jié)果表明多點(diǎn)ZPHI法優(yōu)于ZH-Kdp法，ZPHI法的誤差主要來(lái)源于φDP質(zhì)量不高；胡志群等［7］比較了φDP訂正法和ZH訂正法，指出Kdp訂正法要優(yōu)于ZH訂正法，但是當(dāng)Kdp很小時(shí)會(huì)產(chǎn)生較大的誤差，最后提出ZH-Kdp綜合訂正法，但因ZH-Kdp參數(shù)受溫度等影響較大，難以使用固定的系數(shù)完成衰減訂正；畢永恒等［8］在ZPHI等方法的基礎(chǔ)上，進(jìn)行不斷調(diào)整參數(shù)本身，提出了自適應(yīng)約束方法，該方法能有效改善雷達(dá)的衰減，但在遠(yuǎn)距離處會(huì)出現(xiàn)過(guò)訂正的情況。為了克服單部雷達(dá)反射率衰減訂正難的缺點(diǎn)，Testud等［9］引入了高分辨率天氣雷達(dá)（HRWR）重疊網(wǎng)絡(luò)的概念，可以確定和訂正衰減；Katharina［10］提出結(jié)合C波段雷達(dá)來(lái)訂正X波段雷達(dá)的衰減，結(jié)果表明在層狀云降水、對(duì)流型降水和混合型降水中均得到較好的訂正效果。

本研究綜合分析了多種訂正方法的優(yōu)缺點(diǎn)，認(rèn)為僅使用X波段雷達(dá)本身進(jìn)行衰減訂正是不夠的，擬聯(lián)合S波段雷達(dá)的探測(cè)，采用保序回歸算法對(duì)X波段雷達(dá)進(jìn)行衰減訂正。將S波段雷達(dá)數(shù)據(jù)作為標(biāo)準(zhǔn)來(lái)評(píng)估衰減訂正效果，通過(guò)平均誤差、相關(guān)性等參數(shù)檢驗(yàn)訂正算法的性能，并綜合對(duì)比保序回歸算法、ZH-Kdp訂正法及自適應(yīng)約束方法的訂正效果。

1 雷達(dá)網(wǎng)絡(luò)及性能指標(biāo)

肇慶地區(qū)安裝了兩部X波段全固態(tài)雙偏振天氣雷達(dá)，一部S波段天氣雷達(dá)，地理位置及相對(duì)距離如圖1所示。本研究選擇德慶站X波段雷達(dá)和S波段雷達(dá)作為研究對(duì)象，主要參數(shù)如表1所示。

2 訂正算法

圖1 肇慶雷達(dá)網(wǎng)絡(luò)分布圖

本研究提出了一種結(jié)合X波段和S波段雷達(dá)觀測(cè)數(shù)據(jù)進(jìn)行衰減訂正的新方法，S波段雷達(dá)受衰減影響很小，而且能覆蓋X波段雷達(dá)網(wǎng)絡(luò)區(qū)域，覆蓋X波段雷達(dá)網(wǎng)絡(luò)區(qū)域的S波段雷達(dá)觀測(cè)到比X波段雷達(dá)更粗糙的時(shí)間和空間分辨率的反射率。因此，需要對(duì)S波段反射率場(chǎng)進(jìn)行數(shù)據(jù)轉(zhuǎn)換［11-15］以匹配X波段雷達(dá)極坐標(biāo)網(wǎng)格，每個(gè)X波段雷達(dá)像素被分配其中心所在的S波段雷達(dá)像素的反射率值。

得到同一網(wǎng)格的反射率值Zs（r）和Zx（r）后（其中Zs（r）與Zx（r），單位均為 mm6/m3，r為與雷達(dá)中心的距離（m），下同），定義最大比衰減Kmax（dB）為

表1 X波段與S波段雷達(dá)參數(shù)

派生的Kmax具有與X波段雷達(dá)反射率場(chǎng)相同的數(shù)據(jù)間隙。此外，由于兩個(gè)雷達(dá)系統(tǒng)的誤差，空間和掃描幾何結(jié)構(gòu)的時(shí)間不匹配，Kmax不會(huì)沿著雷達(dá)波束單調(diào)增加中心到外邊界（如圖2黑色十字）。理論上，由于X波段衰減的影響比S波段的影響更大，因此與雷達(dá)的距離越大，Kmax必須增加（或在沒(méi)有降水的情況下保持恒定）。為了得到沿每個(gè)雷達(dá)徑向的Kmax的連續(xù)單調(diào)增加函數(shù)，就需要對(duì)Kmax進(jìn)行保序回歸［10］。保序回歸是統(tǒng)計(jì)上最優(yōu)的方法，以強(qiáng)制單調(diào)增加物理約束，與觀察的偏差最小。保序回歸的概念是計(jì)算給定的單調(diào)遞增函數(shù)y值a1，...，an和它們的正權(quán)重w1，...，wn，以雷達(dá)徑向數(shù)據(jù)作為有限的實(shí)數(shù)集合，訓(xùn)練一個(gè)模型來(lái)最小化方程

其中，y1≤y2≤...≤yn。在這項(xiàng)研究中，本研究使用簡(jiǎn)單的未加權(quán)線性序列保序回歸，即所謂的相鄰最優(yōu)化迭代算法（PAVA）。在這種情況下，所有權(quán)重都等于1。PAVA迭代解決了最小化問(wèn)題，并且基于以下定理來(lái)獲得最優(yōu)解：

如果在迭代的第一步ai＞ai-1，違反了此約束，則設(shè)置新索引j=i。然后aj與aj-1合并為新值a′j-1，根據(jù)以下公式計(jì)算新值a′j-1和新的權(quán)重w′j-1：

使用a′j-1、aj、w′j-1和wj再次測(cè)試式（3）的約束。當(dāng)滿足約束時(shí)，保持a′的值及其位置i、i自動(dòng)加1進(jìn)入下一個(gè)迭代步驟。重復(fù)該過(guò)程直到不再違反式（3）中的約束。最后，將a′的值分配給函數(shù)y中它們對(duì)應(yīng)的位置i。

沿著每個(gè)雷達(dá)徑向和每個(gè)時(shí)間步長(zhǎng)對(duì)Kmax進(jìn)行保序回歸。沿著雷達(dá)徑向的Kmax保序回歸后，KIR呈現(xiàn)為單調(diào)遞增函數(shù)（圖2中的藍(lán)線所示）。為了獲得衰減校正后的反射率dBZXcorr，將該KIR函數(shù)添加到實(shí)測(cè)的徑向X波段反射率dBZX中，即

由于在高分辨率X波段實(shí)測(cè)量中簡(jiǎn)單地增加了校正功能，因此與S波段測(cè)量相比，保持了小尺寸的反射率結(jié)構(gòu)。只有在由于非常強(qiáng)的衰減而導(dǎo)致信號(hào)消失的情況下，才需要使用來(lái)自S波段雷達(dá)的較粗分辨率數(shù)據(jù)作為填充剩余數(shù)據(jù)間隙的最佳可能近似值。

圖2 沿雷達(dá)徑向的K max的保序回歸

3 訂正結(jié)果分析

2019年4月14日肇慶地區(qū)經(jīng)歷了一次較強(qiáng)的降雨過(guò)程，觀測(cè)資料 14：30—14：36（北京時(shí)，下同）時(shí)段的雨區(qū)位于S波段雷達(dá)與X波段雷達(dá)對(duì)比的最佳觀測(cè)范圍，選取該時(shí)段的雷達(dá)數(shù)據(jù)進(jìn)行衰減訂正分析。

本研究選擇X波段雷達(dá)3.4°仰角數(shù)據(jù)進(jìn)行分析（圖3a為反射率產(chǎn)品，圖3b為φDP產(chǎn)品），同時(shí)將S波段雷達(dá)反射率數(shù)據(jù)匹配至X波段雷達(dá)3.4°仰角的極坐標(biāo)網(wǎng)格中（圖3c）。圖3a與3c反射率數(shù)據(jù)進(jìn)行對(duì)比，發(fā)現(xiàn)在方位120°E—240°E范圍內(nèi)，X波段雷達(dá)發(fā)生了強(qiáng)烈的衰減，這一點(diǎn)也可以從φDP值在此范圍內(nèi)出現(xiàn)較大突變得到驗(yàn)證，如圖3b紅色橢圓框內(nèi)的紅色區(qū)域。在兩部雷達(dá)共同探測(cè)區(qū)內(nèi)，靠近X波段雷達(dá)，兩部雷達(dá)都觀察到相似的反射率值，但隨著距離的增加，X波段雷達(dá)出現(xiàn)了衰減直至回波信號(hào)完全消失。

圖3 2019年4月14日14：30德慶X波段雷達(dá)3.4°仰角反射率產(chǎn)品（a）；3.4°仰角φDP產(chǎn)品（b）以及轉(zhuǎn)換后的S波段雷達(dá)3.4°仰角反射率（c）

圖4給出了由S波段雷達(dá)和X波段雷達(dá)反射率數(shù)據(jù)沿X波段雷達(dá)反射率場(chǎng)的155°方位角的分布。在30 km以內(nèi)，X波段雷達(dá)受衰減影響較大，觀測(cè)數(shù)據(jù)與S波段雷達(dá)近似，超出30 km后，X波段雷達(dá)受衰減影響越來(lái)越大，與S波段雷達(dá)探測(cè)數(shù)據(jù)差距也越來(lái)越大，最大比衰減量達(dá)32 dBz。

圖4 2019年4月14日14：30方位角155°的不同訂正結(jié)果的反射率因子分布

新提出的衰減訂正方法結(jié)合了X波段和S波段雷達(dá)反射率數(shù)據(jù)，并且與X波段雷達(dá)的校準(zhǔn)誤差引起的不確定性無(wú)關(guān)。圖5a中的最大比衰減Kmax（由式（1）計(jì)算）的場(chǎng)在雷達(dá)覆蓋區(qū)域的外部顯示高達(dá)20 dB。由于不完美校準(zhǔn)的雷達(dá)系統(tǒng)和雜波消除，Kmax不是單調(diào)遞增的函數(shù)，隨著與雷達(dá)的距離增加或甚至可能產(chǎn)生負(fù)值（圖5a中的數(shù)據(jù)間隙）。因此，執(zhí)行式（2）—式（4）的保序回歸算法后的Kmax場(chǎng)如圖5b所示，其中Kmax隨著雷達(dá)距離的遞增，單調(diào)遞增至覆蓋區(qū)域的最高值。然后將Kmax場(chǎng)添加到X波段雷達(dá)的dBZ場(chǎng)完成反射率數(shù)據(jù)訂正（圖5c）。

圖5 以S和X波段雷達(dá)反射率的對(duì)數(shù)比計(jì)算的最大比衰減K max（a）以及保序回歸后的最大比衰減

獲得的衰減訂正后的dBZ場(chǎng)與圖3c中關(guān)于S波段雷達(dá)觀測(cè)非常一致，如方位120°E—240°E范圍。沿保序回歸訂正的155°方位角的反射率分布在圖4中用藍(lán)色虛線表示，它與受衰減影響區(qū)域的S波段雷達(dá)獲得的反射率分布比較吻合，而且保留了X波段雷達(dá)觀測(cè)到的小尺度變化。即使在距雷達(dá)較遠(yuǎn)的高衰減區(qū)，也能與S波段雷達(dá)觀測(cè)dBZ近似一致。兩個(gè)系統(tǒng)之間的平均偏差僅為1.1 dBz，標(biāo)準(zhǔn)偏差為4.2 dB，相關(guān)系數(shù)為0.88。此外，訂正后的X波段反射率場(chǎng)顯示了比S波段雷達(dá)更多的降水事件細(xì)節(jié)。

為了檢驗(yàn)保序回歸衰減訂正效果，比較了保序回歸訂正法、ZH-Kdp訂正法及自適應(yīng)約束算法的衰減訂正效果（圖6）。從圖3a與圖6a對(duì)比可以看出，反射率場(chǎng)得到一定程度的訂正，在方位180°附近尤為明顯，比較接近S雷達(dá)探測(cè)的數(shù)據(jù)（圖5c），但是大部分區(qū)域訂正后的回波很稀疏，強(qiáng)回波不夠集中，主要原因是雙偏振量Kdp數(shù)據(jù)難以控制，質(zhì)量不佳，導(dǎo)致訂正效果不理想。從圖6b與圖3a進(jìn)行對(duì)比，可以看出回波整體變化不明顯，框框內(nèi)回波還是有訂正的痕跡，反射率值有明顯變化。但整體變化不大，與S波段雷達(dá)探測(cè)結(jié)果差距較大（圖5c）。為了更加直觀的比較不同訂正方法的訂正效果，圖4增加了ZH-Kdp法（磚紅色虛線）和自適應(yīng)約束算法（紫色虛線）訂正后的反射率數(shù)據(jù)沿X波段雷達(dá)反射率場(chǎng)的155°方位角的分布。從圖4可以看出，在近X波段雷達(dá)處，衰減較小，與S波段觀測(cè)數(shù)據(jù)差異較小，幾種訂正法訂正后的數(shù)據(jù)也與原始數(shù)據(jù)近似，隨著距離的越來(lái)越大，X波段雷達(dá)受到衰減的強(qiáng)烈影響，在距離約為30 km處反射率衰減較大，與S波段雷達(dá)觀測(cè)的回波數(shù)據(jù)差距也越來(lái)越大。

圖6 2019年4月14日14：30衰減訂正后的反射率場(chǎng)（單位：dBz）

從3種不同訂正方法的相關(guān)系數(shù)及平均誤差（表2）可以看出，保序回歸算法相關(guān)系數(shù)為0.88，均高于其它兩種訂正法，平均誤差為1.1dBz，遠(yuǎn)低于其它兩種訂正法的結(jié)果，可以看出保序回歸訂正后的回波與S波段雷達(dá)具有較高的一致性，達(dá)到了預(yù)期的目標(biāo)。

表2 不同訂正方法的相關(guān)系數(shù)和平均誤差

為了進(jìn)一步分析衰減訂正的效果，選取方位120°E—240°E范圍內(nèi)降水回波訂正前后的反射率散點(diǎn)圖（圖7）。從圖7a中可看出，反射率因子大于35 dBz，絕大部分點(diǎn)分布在對(duì)角線的下側(cè)，說(shuō)明隨著回波值增大，衰減越來(lái)越嚴(yán)重。圖7c較圖7b點(diǎn)更集中于對(duì)角線，說(shuō)明ZH-Kdp訂正法優(yōu)于自適應(yīng)算法，散點(diǎn)圖效果不如曲線圖4效果明顯。圖7d是保序回歸算法后的散點(diǎn)圖，可以看出絕大部分點(diǎn)都集中在對(duì)角線兩側(cè)，說(shuō)明訂正后的X波段反射率數(shù)據(jù)與S波段雷達(dá)反射率數(shù)據(jù)保持了較高的一致性，訂正效果較好，同時(shí)存在一些凌亂、抖動(dòng)的數(shù)據(jù)，原因可能是受到雷達(dá)標(biāo)定及雜波等因素的影響。

圖7 訂正前后的反射率因子與S雷達(dá)反射率因子的散點(diǎn)圖

綜上所述，ZH-Kdp訂正法及自適應(yīng)約束算法在遠(yuǎn)雷達(dá)區(qū)衰減訂正改善很小，而且在大規(guī)模強(qiáng)降水或雜波污染的情況下，這種方法會(huì)變得不穩(wěn)定，因此在強(qiáng)雨帶衰減的背景下，基于保序回歸的訂正法更好，能較好訂正X波段雷達(dá)的衰減回波，與S雷達(dá)數(shù)據(jù)保持較高的相關(guān)性。

4 結(jié)論

1）X波段雷達(dá)受衰減影響大，使得經(jīng)過(guò)強(qiáng)雨帶后回波強(qiáng)度變?nèi)?，甚至探測(cè)不到回波，可以利用共同探測(cè)區(qū)域的S波段雷達(dá)數(shù)據(jù)對(duì)衰減區(qū)域進(jìn)行訂正，彌補(bǔ)單部雷達(dá)探測(cè)的不足。

2）從反射率訂正前后回波圖對(duì)比分析，基于保序回歸訂正后的回波與S波段雷達(dá)回波有較高的一致性。

3）從徑向數(shù)據(jù)曲線圖分析，在近距離處，3種訂正結(jié)果相差不大，說(shuō)明X波段雷達(dá)在降水強(qiáng)度不大的情況下衰減較?。辉谶h(yuǎn)距離處，保序回歸法的訂正效果要優(yōu)于其它兩種方法，訂正后的回波更接近S波段雷達(dá)回波；而ZH-Kdp訂正法優(yōu)于自適應(yīng)約束算法，但在遠(yuǎn)距離處均出現(xiàn)了訂正量不夠的情況，說(shuō)明在強(qiáng)雨帶衰減的背景下，僅使用X波段雷達(dá)本身進(jìn)行衰減校正是不夠的，體現(xiàn)了聯(lián)合S波段雷達(dá)進(jìn)行衰減訂正的優(yōu)勢(shì)。