長沙機場陣列天氣雷達風場驗證

李 渝 馬舒慶 楊 玲* 甄小瓊4) 喬 丹

1)(成都信息工程大學電子工程學院，成都 610225) 2)(中國氣象局大氣探測重點開放實驗室，成都 610225) 3)(中國氣象局氣象探測中心，北京 100081) 4)(中國科學院大氣物理研究所，北京 100029) 5)(雷象科技(北京)有限公司，北京 100089)

引 言

多普勒天氣雷達能獲得降水粒子的反射率因子、徑向速度和譜寬等，其中徑向速度是有效照射體積內(nèi)降水粒子沿波束方向上的平均速度，但徑向速度只是實際風場中的一個分量，不能直接代表實際風。利用徑向速度合成或反演三維風場，是雷達氣象學中的重要研究方向。精細的三維風場可以研究中小尺度天氣系統(tǒng)結(jié)構(gòu)和運動特性。

國內(nèi)外關(guān)于雙(多)多普勒天氣雷達聯(lián)合探測風場的研究工作和外場試驗已取得一定成果，主要技術(shù)方法分為以下3類：①3部或以上多普勒天氣雷達直接使用公式合成共同觀測區(qū)域內(nèi)的風場，Armijo[1]和Ray等[2]在笛卡爾坐標系下給出求解方程組，韓頌雨等[3]在動態(tài)地球坐標系下給出求解方程組；②對于雙多普勒天氣雷達風場反演方法，Lhermitte等[4]、Chong等[5]、Fankhaus等[6]、劉黎平等[7-8]和羅昌榮等[9]通過雙多普勒天氣雷達同步觀測獲得散射體兩個方向的徑向速度資料,利用這兩個方向的資料與其他約束條件(如質(zhì)量連續(xù)方程)聯(lián)立求解大氣風場；③雙(多)多普勒天氣雷達三維變分(three-dimensional variational data assimilation，3DVAR)風場反演方法。

第3類方法中通常以是否考慮天氣系統(tǒng)平移和內(nèi)演變化情況又分為兩類：①未考慮天氣系統(tǒng)平移和內(nèi)演變化。Scialom等[10]對多部多普勒天氣雷達進行分析，以標準正交函數(shù)的形式建立徑向速度的分析模型，將地面邊界層約束和質(zhì)量連續(xù)方程作為變分約束加入分析提高垂直速度的反演能力。Gao等[11-12]對雙多普勒天氣雷達變分反演方法進行優(yōu)化，其中以質(zhì)量連續(xù)方程為弱約束條件，采用準牛頓共軛梯度算法，通過有限的內(nèi)存使代價函數(shù)最小化提高運算效率。Liou等[13]改善多部多普勒天氣雷達變分反演方法的效果，解決傳統(tǒng)方法中遇到的迭代不穩(wěn)定性問題，進一步提高風場反演的精度。王艷春等[14-15]分析雙多普勒天氣雷達三維變分方法反演風場的能力，誤差分析結(jié)果表明:水平風速和水平風向的反演結(jié)果均比較可靠。但此類方法未考慮雷達體掃時真實天氣系統(tǒng)也在快速移動的事實，未進行平流校正，資料時空不匹配會帶來一定誤差。②考慮天氣系統(tǒng)平移和內(nèi)演變化。Liu等[16]在3DVAR反演算法中加入時間線性插值，它能減少傳統(tǒng)方法忽略觀測天氣系統(tǒng)時間演化的誤差，反演結(jié)果能捕獲到風場的快速演化過程。Shapiro等[17]在3DVAR反演算法中引入泰勒冰凍湍流假設(shè)，通過計算兩次體掃資料的反射率因子可以大致推導實際數(shù)據(jù)的時空位置，從而減少反演誤差，在不太強的中小尺度天氣系統(tǒng)中效果較好。在此基礎(chǔ)上，Potvin等[18]進一步優(yōu)化算法，進行平流校正和內(nèi)演變化改善，水平風場的平移和演化對于提高垂直風速的反演精度也起到巨大作用，并且得出雷達體掃時間越短，反演效果越好的試驗結(jié)論。

目前，國內(nèi)對雙(多)多普勒天氣雷達3DVAR風場反演方法的外場試驗驗證研究較少。本文使用部署在湖南省長沙黃花機場的陣列天氣雷達資料，由于雷達探測范圍較小，且在嚴格的協(xié)同探測模式下3個收發(fā)子陣共同覆蓋的三維精細探測區(qū)域資料時差在2 s以下，保證各個收發(fā)子陣同時探測到大氣中各點徑向速度的空間和時間一致性。這能夠減小最終風場合成和反演的誤差[19]。目前降水云中真實的三維風場信息很難獲得，現(xiàn)有風場觀測資料的時空分辨率較低，二者難以匹配，風場的“真值”很難衡量。本文為了驗證陣列反演風場的有效性，使用機場1部L波段邊界層風廓線雷達產(chǎn)品和陣列合成風場結(jié)果作為參考值，與陣列反演風場結(jié)果進行相互比較。在3個收發(fā)子陣共同覆蓋的三維精細探測區(qū)域(空間同一點存在3個徑向速度時)，將陣列合成風場作為參考值，對陣列反演風場進行驗證評估。

1 陣列天氣雷達簡介和研究個例

本文使用的長沙機場陣列天氣雷達由3個相控陣接收發(fā)射子陣(簡稱收發(fā)子陣)為一組進行同步探測，在水平方向上每個收發(fā)子陣采用機械掃描方式，覆蓋0°～360°方位，垂直方向具有4個發(fā)射波束和64個接收波束，覆蓋0°～90°仰角。每個收發(fā)子陣的最大探測距離為20.28 km，徑向分辨率為30 m[20]。張沛源等[21]研究3部多普勒天氣雷達聯(lián)合測量大氣風場的誤差分布及最佳布局，從理論上證明3部雷達以等邊三角形布局為最佳。機場陣列天氣雷達布局近似于等邊三角形(圖1)，在3個收發(fā)子陣測速精度相同的條件下，能減小風場合成和反演的誤差。外場試驗中收發(fā)子陣天線旋轉(zhuǎn)速度為30°·s-1，1次體掃時間為12 s。單個收發(fā)子陣在60°范圍內(nèi)掃描時間為2 s，120°范圍內(nèi)掃描時間為4 s。故在嚴格的協(xié)同探測模式下，3個收發(fā)子陣在共同覆蓋區(qū)域(簡稱三維精細探測區(qū))的探測資料時差在2 s以下，兩個收發(fā)子陣共同覆蓋區(qū)域的探測資料時差在4 s以下。圖1中三維精細探測區(qū)外也能獲取探測資料，這些區(qū)域被稱為普通探測區(qū)。在進行風場合成和反演前，對雷達基數(shù)據(jù)進行地物雜波濾除、徑向速度退模糊等質(zhì)量控制處理[22]，并且將3個收發(fā)子陣的體掃資料通過Cressman方法插值到水平和垂直分辨率均為100 m的直角坐標系網(wǎng)格場中。

圖1 3個收發(fā)子陣陣列天氣雷達布局及探測區(qū)示意圖(紅色矩形為3DVAR風場反演范圍，紅色圓點為機場風廓線雷達站點)

2019年4—9月長沙機場地區(qū)降水較多，也是強對流天氣發(fā)生最多的月份，本文選取機場外場試驗中10次降水過程，按照降水類型分為兩類,第1類為穩(wěn)定性降水，第2類為對流性降水。具體降水個例及時間(北京時，下同)、與風廓線雷達產(chǎn)品對比時段及陣列合成風場與陣列反演風場選用的誤差對比分析時刻見表1。

表1 降水個例及描述

2 風場算法簡介

2.1 陣列合成風場算法

陣列合成風場計算原理與文獻[2]中3部多普勒天氣雷達分析方法近似，在插值后的直角坐標系網(wǎng)格場中，每個收發(fā)子陣坐標為(xi,yi,zi)，空間任意一點坐標在網(wǎng)格場中表示為(x,y,z)。u，v和w(w=w+wt)代表在x,y,z方向的風分量(wt為降水粒子的下落速度[23])，Ri為網(wǎng)格點到每個收發(fā)子陣的距離。每個收發(fā)子陣測得的徑向速度Vi在直角坐標系中與風分量的公式為

(1)

(2)

聯(lián)立式(1)和式(2)可在網(wǎng)格場存在3個徑向速度時計算出風分量u，v和w。

2.2 陣列反演風場算法

陣列反演風場使用Shapiro等[17]、Potvin等[18]和North等[24]提出的3DVAR風場反演算法，具體算法思想如下：記直角坐標系網(wǎng)格場中每個格點上3個方向的風分量為u，v，w，每個收發(fā)子陣坐標為(xi,yi,zi)，空間任意一點坐標在網(wǎng)格場中表示為(x,y,z)。利用徑向速度觀測約束JO，質(zhì)量守恒約束JM，地面邊界層約束JP和空間平滑約束JS對3個收發(fā)子陣的徑向速度觀測資料進行迭代求解，風場反演結(jié)果在J的全局極小值處取得最優(yōu)解。定義代價函數(shù)J：

J(u，v，w)=JO+JM+JP+JS。

(3)

等式右邊第1項JO定義為

(4)

(z-zi)[wi-|wti|]}，

(5)

(6)

第2項JM定義為

i=1,…,n。

(7)

基態(tài)大氣參考密度為

(8)

其中，ρs單位為kg·m-3，z單位為m。

第3項JP定義為

(9)

其中，非地面的網(wǎng)格點權(quán)重系數(shù)設(shè)置為0。

第4項JS定義為

+λSw(2w)2,

(10)

(11)

上述公式中，每一項代價函數(shù)的權(quán)重系數(shù)采用North等[24]在外場試驗中通過靈敏度分析試驗得到的最優(yōu)解結(jié)論。權(quán)重系數(shù)分別設(shè)置為λO=1，λM=500，λP=1，λSu=1，λSv=1，λSw=0.1。

3 三維精細探測區(qū)域風場驗證

進行三維精細探測區(qū)域風場驗證前，檢查選用個例徑向速度資料是否存在較大偏差。具體做法如下：首先，找到兩個收發(fā)子陣連線上的中點位置A，B，C(子陣1與子陣2、子陣1與子陣3、子陣2與子陣3)；然后，取出該點連續(xù)5個時刻不同高度上徑向速度的平均值；最后，比較兩個收發(fā)子陣在該點的值。該位置上兩個子陣的徑向速度在同一高度的大小應(yīng)接近且方向相反。雷達的徑向速度本身存在觀測誤差約為1 m·s-1，實際探測中還存在不確定的偏差[25]。因此，徑向速度會存在一定范圍的偏差。表2中給出個例3在16:00:00—16:00:48連續(xù)5個時刻各個收發(fā)子陣徑向速度的一致性分析結(jié)果，可以看出徑向速度觀測數(shù)據(jù)無較大偏差。

表2 不同子陣連線上中點位置的徑向速度一致性分析

3.1 陣列反演風場和陣列合成風場與風廓線雷達產(chǎn)品誤差

本節(jié)重點討論陣列反演風場和陣列合成風場與風廓線雷達產(chǎn)品誤差，風廓線雷達是風場探測的有利工具[26-28]，其位置見圖1紅色圓點。精確的探測位置是定點經(jīng)緯度(28.2°N，113.2°E)上空的風場，高度780 m以下垂直分辨率為60 m，高度780～5220 m 垂直分辨率為120 m。通過風廓線雷達站點經(jīng)緯度坐標匹配網(wǎng)格場中對應(yīng)位置，以位置為中心，計算半徑R=2500 m內(nèi)陣列反演和陣列合成風場每一層高度的平均水平風場。降水過程中回波經(jīng)過此位置上空時間較短，且陣列反演風場和陣列合成風場有時沒有風場解，故半徑取值較大。分析時間段見表1。選取風廓線不同高度隨時間變化的風場產(chǎn)品作為參考值與陣列反演風場和陣列合成風場進行相互比較。采用平均絕對偏差、均方根誤差、相對均方根誤差(相對參考值的均方根誤差百分比)3個評價參數(shù)分別評估分析時間段內(nèi)水平風速和水平風向兩個物理量，水平風向的平均絕對偏差和均方根誤差已經(jīng)能客觀表現(xiàn)風向誤差，不再計算水平風向的相對均方根誤差。

3.1.1 穩(wěn)定性降水

表1中有5次穩(wěn)定性降水過程，圖2是個例3在2019年5月12日15:30—17:30時段，風廓線雷達產(chǎn)品、陣列反演風場和陣列合成風場在風廓線站點上空高度隨時間變化的水平風場(風羽)。為了與風廓線雷達產(chǎn)品5 min 的時間分辨率匹配，本文選擇最接近風廓線雷達探測時間的陣列天氣雷達資料進行對比。圖2中水平軸時間間隔為5 min，垂直軸高度間隔為500 m，無風矢的位置表示無風。

由圖2可看到，風廓線雷達產(chǎn)品、陣列反演風場和陣列合成風場的水平風速在1 km和1.5 km高度存在差異，但風向隨高度變化趨勢一致。3種方法在2～5 km 高度表現(xiàn)為均勻的西南風，風速、風向相近。在風場均勻的降水過程中，風廓線雷達測風精度較準確[29]。將風廓線雷達產(chǎn)品的水平風作為參考值，分別計算陣列反演風場和陣列合成風場5個個例分析時間段內(nèi)水平風速和水平風向的平均絕對偏差和均方根誤差以及水平風速的相對均方根誤差。由表3可見，陣列反演風場風速和風向的平均絕對偏差分別為1.28～3.96 m·s-1,5.72°～10.81°，表明水平風在分析時段變化不大。表3中水平風速的均方根誤差為2.14～3.92 m·s-1，相對均方根誤差為19%～31%。陣列反演風場的水平風速存在一定誤差，但相對偏差在31%以下，水平風速比較合理。水平風向的均方根誤差為7.56°～17.69°，與風廓線雷達產(chǎn)品的風向差小于20°，風向表現(xiàn)較一致。陣列合成風場的誤差評價參數(shù)差別不大(表略)，說明陣列反演和陣列合成風場比較一致?？烧J為風廓線雷達產(chǎn)品，陣列反演風場和陣列合成風場3種方法的結(jié)果較一致。

圖2 2019年5月12日15:30—17:30不同高度的水平風場(a)風廓線雷達產(chǎn)品，(b)陣列反演風場，(c)陣列合成風場

表3 個例分析時段內(nèi)陣列反演風場與風廓線雷達產(chǎn)品水平風速和水平風向的平均絕對偏差、均方根誤差和相對均方根誤差

3.1.2 對流性降水

表1中有5次對流性降水過程,由于個例8降水回波未經(jīng)過風廓線雷達站點上空，所以只對其余個例進行對比分析。圖3為個例1在2019年4月26日18:00—20:00 時段風廓線雷達產(chǎn)品、陣列反演風場和陣列合成風場在風廓線站點上空高度隨時間變化的水平風場(風羽)?？梢钥吹?，在18:00—18:30風廓線雷達產(chǎn)品風場出現(xiàn)較大變化，陣列反演風場和陣列合成風場整體變化趨勢不大。與風廓線雷達產(chǎn)品相比，陣列反演風場和陣列合成風場(表略)的誤差評價參數(shù)差別不大，說明陣列反演和陣列合成風場比較一致。表4中風速和風向的平均絕對偏差分別為1.52～2.94 m·s-1和39.79°～55.89°，表明在分析時段內(nèi)水平風速變化不大，但水平風向變化顯著。表4中水平風速的均方根誤差為2.66～5.57 m·s-1，相對均方根誤差為44%～73%；水平風向的均方根誤差為33.88°～58.24°。表明水平風速和水平風向在誤差參數(shù)評價下表現(xiàn)不好，與風廓線雷達產(chǎn)品相比誤差較大。對流性降水出現(xiàn)時環(huán)境風場不均勻會造成風廓線雷達測得的水平風向、風速的測量誤差較大，測風精度下降[29]，故在對流性降水過程中，風廓線雷達產(chǎn)品與陣列反演風場和陣列合成風場差異較大。但陣列反演風場和合成風場在風廓線站點上水平風表現(xiàn)較一致，也驗證了在降水場很不均勻時風廓線反演的水平風會出現(xiàn)誤差的事實。

表4 個例分析時段內(nèi)陣列反演風場與風廓線雷達產(chǎn)品水平風速和水平風向的平均絕對偏差、均方根誤差和相對均方根誤差

圖3 2019年4月26日18:00—20:00不同高度的水平風場(a)風廓線雷達產(chǎn)品，(b)陣列反演風場，(c)陣列合成風場

3.2 兩種風場算法的誤差分析

本節(jié)重點討論兩種風場算法對不同降水類型求解的水平風速與水平風向誤差。為了將陣列合成風場作為參考值，本節(jié)陣列反演風場只計算有3個徑向速度存在的網(wǎng)格點，保證誤差計算過程中兩種算法均在有風場解的格點進行對比。降水類型雖不同，但兩種風場算法在同一時刻求解結(jié)果應(yīng)接近。不同時刻的誤差分析結(jié)果差別不大，以10次降水過程中某一時刻的計算結(jié)果為例進行誤差分析，選用的對比體掃時刻見表1，誤差評價參數(shù)與3.1節(jié)中一致。

3.2.1 穩(wěn)定性降水

風場求解區(qū)域為圖1中紅色矩形區(qū)域，網(wǎng)格場水平范圍大小為25 km×25 km，高度為20 km，水平和垂直方向網(wǎng)格間距均為100 m。選取表1中5個穩(wěn)定性降水過程，圖4為各個高度層水平風速和水平風向的均方根誤差，以及水平風速的相對均方根誤差。由圖4可以看到，水平風速的均方根誤差為0.36～2.28 m·s-1，相對均方根誤差為1%～19%；水平風向的均方根誤差為1.51°～14.92°。水平風速的相對偏差低于20%，水平風向差小于15°，即在穩(wěn)定性降水條件下水平風速和水平風向在各個高度層較為一致，誤差較小。

圖4 5個穩(wěn)定性降水個例誤差(a)水平風速均方根誤差，(b)水平風速相對均方根誤差，(c)水平風向均方根誤差

2019年8月25日16:50—18:20陣列天氣雷達探測到降水回波自東移入三維精細探測區(qū)，向西南方向發(fā)展。風場與回波移動方向幾乎一致，低層為均勻的東北風，中高層轉(zhuǎn)為偏東風。圖5a、圖5b為8月25日17:05:12陣列合成風場和陣列反演風場在3 km，5 km高度上的水平風羽圖，選取紅色矩形區(qū)域內(nèi)10 km × 10 km范圍進行展示。由圖5可知，3 km高度為均勻的東北風，5 km高度為比較均勻的偏東風。兩種算法的水平風場空間分布和大小方向非常接近，在視覺上看不出明顯差異。

圖5 2019年8月25日17:05:12陣列合成風場和陣列反演風場不同高度的水平風矢(填色為反射率因子)(a)陣列合成風場,3 km高度，(b)陣列反演風場,3 km高度，(c)陣列合成風場,5 km高度,(d)陣列反演風場,5 km高度

3.2.2 對流性降水

與3.2.1節(jié)一致，將陣列合成風場作為參考值，選取表1中5個對流性降水過程求解的陣列反演風場進行誤差分析。對流性降水過程中回波強度大、持續(xù)時間短、隨時間變化迅速，且風場復(fù)雜多變。圖6為各個高度層水平風速和水平風向的均方根誤差，以及水平風速的相對均方根誤差。由圖6可知，水平風速的均方根誤差為0.19～3.73 m·s-1，相對均方根誤差為1%～29%；水平風向的均方根誤差為2.04°～26.35°。水平風速的相對偏差低于30%，水平風向差小于30°。與穩(wěn)定性降水過程相比，水平風向差出現(xiàn)較大差異是由于陣列合成風場是計算網(wǎng)格場中每一個點的風分量，而陣列反演風場是網(wǎng)格場中全局風分量的最優(yōu)解，式(10)中還對結(jié)果進行平滑處理。故在對流性降水過程中，陣列反演風場的每個格點風向經(jīng)過平滑后比陣列合成風場整體偏向更為一致。結(jié)果表明：兩種算法在對流性降水條件下，水平風速和水平風向在各個高度層較為合理。與穩(wěn)定性降水過程相比，對流性降水水平風速和水平風向誤差略大。

圖6 5個對流性降水個例誤差 (a)水平風速均方根誤差，(b)水平風速相對均方根誤差，(c)水平風向均方根誤差

2019年8月21日14:00—16:00中等強度的分散對流云團在三維精細探測區(qū)內(nèi)以40 km·h-1的速度向南偏西方向移動，其間存在短時強降水并伴有雷暴天氣。圖7為8月21日14:52:00時刻陣列合成風場和陣列反演風場在3 km和5 km高度上的水平風羽圖。由圖7可知，局地強回波位于三維精細探測區(qū)之間，周圍為層云區(qū)，中心最強回波超過45 dBZ。從風場結(jié)構(gòu)看，兩種算法結(jié)果比較一致。3 km和5 km高度分別以東北風和偏東風為主，而且出現(xiàn)明顯的渦旋結(jié)構(gòu)(圖7中紅色圓圈)；風速、風向變化大，且在強回波區(qū)域風速較大?？傮w上，從視覺上看風場細節(jié)存在差異，但空間分布較一致。

圖7 2019年8月21日14:52:00陣列合成風場和陣列反演風場不同高度的水平風矢(填色為反射率因子)(a)陣列合成風場,3 km高度，(b)陣列反演風場,3 km高度，(c)陣列合成風場,5 km高度,(d)陣列反演風場,5 km高度

4 結(jié)論與討論

利用2019年4—9月高時空分辨率的陣列天氣雷達資料開展3DVAR風場反演。為探討其結(jié)果在不同降水類型下的反演能力，選用10個降水個例，在三維精細探測區(qū)采用1部L波段邊界層風廓線雷達產(chǎn)品和陣列合成風場作為參考值對其進行驗證評估。得到以下結(jié)論：

1)在穩(wěn)定性降水條件下，反演和合成結(jié)果與風廓線雷達產(chǎn)品較為一致，比較合理。在對流性降水條件下，得到的結(jié)果誤差較大。由于在對流性降水中環(huán)境風場不均勻性會造成風廓線雷達測風精度下降，與陣列反演和陣列合成風場相比出現(xiàn)較大差異。

2)對比不同降水條件下的陣列反演風場與陣列合成風場，兩種算法得到的風場結(jié)構(gòu)符合各類天氣系統(tǒng)的基本特征，水平風場的空間分布和風速、風向非常相近。

3)在穩(wěn)定性降水條件下，兩種算法得到的水平風速相對偏差低于19%,水平風向差小于14.92°；在對流性降水條件下，兩種算法得到的水平風速相對偏差低于29%,水平風向差小于26.35°?？傮w上，兩種算法效果較一致，誤差在可接受范圍內(nèi)，陣列反演風場在穩(wěn)定性降水條件下優(yōu)于對流性降水。

盡管本文已對陣列天氣雷達反演風場的結(jié)果進行不同的驗證評估，但由于外場觀測資料時空分辨率不匹配和分析方法有限，仍存在不足：①由于計算原理不同，利用陣列合成風場對陣列反演風場進行定量評估，無法對兩個子陣共同覆蓋的區(qū)域進行誤差分析；②風廓線雷達與陣列天氣雷達的測風原理不同，與風廓線雷達產(chǎn)品對比分析時未考慮穩(wěn)定性降水和對流性降水過程的不同偏差，后期還需完善資料對比分析的相對偏差統(tǒng)計，從而進一步提高陣列天氣雷達反演風場的準確性和實際應(yīng)用的可靠性；③將陣列天氣雷達資料從極坐標系插值到直角坐標系中存在一定誤差，導致求解的風場結(jié)果也有誤差；④目前，陣列天氣雷達反演風場研究中的降水粒子下落末速度是以雷達反射率因子參數(shù)化進行估計，降水粒子下落末速度計算如果存在較大誤差會影響最終垂直速度的結(jié)果。本文未對垂直速度進行驗證分析，未來將開展對降水粒子下落末速度的估算，進一步討論反演垂直速度的準確性和合理性。