王 昊，丁 楠，張文淵，馮遵德，趙長勝，閆香蓉

1.江蘇師范大學地理測繪與城鄉(xiāng)規(guī)劃學院，江蘇 徐州 221116；2.中國礦業(yè)大學環(huán)境與測繪學院，江蘇 徐州 221116

大氣水汽是影響氣候變化的關(guān)鍵氣象要素，對大氣能量輸送和對流天氣的演變具有重要影響作用[1-2]。傳統(tǒng)水汽探測手段，如無線電探空技術(shù)、水汽微波輻射計等，由于具有諸多局限性，已不能滿足目前實時、快捷、精確地獲取水汽要素的要求。隨著GNSS技術(shù)在氣象學中的應(yīng)用發(fā)展，利用GNSS技術(shù)進行對流層三維水汽場的反演已快速成為一種新興水汽探測手段。GNSS水汽層析技術(shù)具有精確、快速、全天候和低成本等諸多優(yōu)勢[3]，已在大氣水汽監(jiān)測中得到越來越多的應(yīng)用。

文獻[4]率先發(fā)展了GNSS水汽層析技術(shù)，研究了垂直分層對水汽層析結(jié)果精度的影響規(guī)律，基于夏威夷Kilauea GPS網(wǎng)進行了水汽層析試驗，通過對比歐洲中期天氣預(yù)報中心(European Centre for Medium-Range Weather Forecasts，ECMWF)的數(shù)據(jù)，得出了垂直分層的最小分辨率需大于300 m以削弱噪聲影響并保證層析結(jié)果精度的結(jié)論。垂直分層是GNSS水汽層析的關(guān)鍵技術(shù)之一，基于文獻[4]的水汽層析思路，諸多學者對這一關(guān)鍵技術(shù)進行了研究，目前主要分為兩類：一類為采用恒定高度間隔的均勻分層；另一類為采用不同高度間隔的非均勻分層。在均勻分層的研究中，文獻[5]采用3種不同高度間隔(600、800、1000 m)的等距均勻垂直分層進行層析試驗，發(fā)現(xiàn)采用均勻分層時以800 m為間隔最優(yōu)。在非均勻分層的研究中，文獻[6]以四川省為例對非均勻分層與均勻分層進行對比試驗，研究并分析了采用非均勻分層對層析結(jié)果的影響。文獻[7]基于探空數(shù)據(jù)研究了層析區(qū)域的頂層邊界高度，提出了兩種不同非均勻分層模型，并對比了兩種模型的衛(wèi)星信號利用率和體素塊利用率。文獻[8]提出了一種體素塊優(yōu)化的層析區(qū)域建模思路，在垂直方向上以1973—2012年探空水汽廓線為依據(jù)對香港地區(qū)進行4種不同的垂直分層建模，對比分析了不同高度間隔的非均勻分層模型對層析結(jié)果質(zhì)量的影響。

傳統(tǒng)的垂直分層建模方法將研究區(qū)域按等高度間距于垂直方向進行劃分[9]，該方法屬于均勻分層，優(yōu)勢在于操作簡便且分層后獲得的離散區(qū)域形狀規(guī)則[10]。由于大氣水汽在垂直方向上呈指數(shù)遞減趨勢，因此建模時采用均勻分層的不足之處在于均勻的垂直分辨率不符合大氣水汽的實際分布，會引入一定的誤差。而非均勻分層模型的垂直分辨率更加符合大氣水汽的實際分布，采用非均勻分層模型進行離散化能有效減小上述誤差。為削弱垂直分層建模過程中誤差的影響，非均勻分層模型已被廣泛運用于GNSS水汽層析技術(shù)之中[11-12]，但目前非均勻分層的建模方法仍然存在不足，非均勻分層模型的垂直分辨率設(shè)置尚未形成科學可靠且適用性強的參考標準，部分層析試驗在采用非均勻分層模型時僅根據(jù)固定的經(jīng)驗設(shè)置上疏下密的分層高度區(qū)間，并無科學可靠的分層依據(jù)[13]；另一部分試驗雖基于先驗水汽密度進行建模，但仍存在分層模型適用性不高的問題，在層析區(qū)域變更或大氣水汽垂直分布特征發(fā)生變化的情況下無法實現(xiàn)方便快捷的自適應(yīng)分層[7-8,14]。

針對上述問題，本文提出一種自適應(yīng)非均勻指數(shù)分層(adaptive non-uniform exponential vertical stratification，ANEVS)的GNSS水汽層析方法，詳細闡述了該垂直分層方法的分層步驟和具體公式。該方法基于層析區(qū)域內(nèi)大氣水汽于垂直方向上呈指數(shù)遞減的分布特征，能夠自適應(yīng)地完成任意給定層析區(qū)域的非均勻垂直分層建模，優(yōu)化各層水汽密度參數(shù)的分配，縮小底層各層水汽密度參數(shù)之間的差異，提高了垂直分層的合理性和真實性，進而改善層析結(jié)果質(zhì)量，提高層析解算精度，同時彌補了垂直分層模型適用性低的問題。

1 GNSS水汽層析基本原理

GNSS水汽層析技術(shù)將連續(xù)的研究區(qū)域離散化為數(shù)量有限的體素塊，以衛(wèi)星信號的斜路徑濕延遲(slant wet delay，SWD)為觀測值構(gòu)建大氣水汽參數(shù)、信號傳播距離與信號傳播路徑上的總水汽含量[15]之間的函數(shù)關(guān)系，反演出水汽的三維空間分布信息[16]。

GNSS水汽層析以每條信號的SWD為觀測值。忽略信號傳播過程中的彎曲量，假設(shè)傳播路徑為三維空間直線，同時假設(shè)各體素塊內(nèi)水汽密度在一個層析歷元內(nèi)為常數(shù)。在此基礎(chǔ)上建立模型為

(1)

式中，SWDi表示第i條衛(wèi)星信號的斜路徑濕延遲；lij表示第i條信號穿過第j個體素塊的直線長度；Dj表示第j個體素塊內(nèi)的水汽參數(shù)，這里為該體素塊的待求水汽密度；n表示總體素塊個數(shù)。

基于上述模型構(gòu)建的觀測方程組系數(shù)矩陣嚴重秩虧，因此在構(gòu)建層析方程組時，一般通過加入約束條件以解決層析方程組的不適定性問題[17-18]。約束條件主要包括垂直約束、水平約束和頂層約束。垂直約束采用指數(shù)函數(shù)進行約束，其模型為[19]

(2)

式中，hi和hj分別表示大氣高度；ρi和ρj分別表示高度hi和hj處的水汽密度；H表示水汽標高[20]。

水平約束采用高斯加權(quán)函數(shù)進行約束，其模型為[21]

w1D1+w2D2+…+wi-1Di-1+wi+1Di+1+

+wjDj=Di

(3)

(4)

式中，Di表示同層第i個體素塊內(nèi)的水汽密度；wi表示同層第i個水汽密度的權(quán)系數(shù)；dij表示第i個體素塊與同層第j個體素塊中點之間的距離；σ表示平滑因子[22]，其取值與體素塊水平方向的長度有關(guān)[21]，本文的平滑因子取值為15。

頂層約束將最頂層體素塊內(nèi)的水汽密度按照經(jīng)驗強制約束為某一定值，部分層析試驗將頂層水汽密度值約束為0.01 g/m3[23]，但頂層體素塊實際具有一定高度，僅頂層體素塊頂端水汽密度值接近0.01 g/m3，中間及底端水汽密度大于0.01 g/m3，因此本文以先驗探空數(shù)據(jù)為依據(jù)，采用頂層體素塊高度區(qū)間內(nèi)的探空水汽密度均值為約束值。

結(jié)合觀測方程和約束條件的層析方程組可以表示為

(5)

式中，A0、AV、AH、AT分別表示觀測值系數(shù)矩陣、垂直約束系數(shù)矩陣、水平約束系數(shù)矩陣、頂層約束系數(shù)矩陣；DW表示待求水汽密度向量；b0、bV、bH、bT分別表示對應(yīng)的觀測值向量，其中，bV、bH一般情況下為0向量；方程總個數(shù)m等于各部分方程個數(shù)o、t、v、h之和。

式(5)的解算方法有奇異值分解法、最小二乘法、代數(shù)重構(gòu)法(algebraic reconstruction technique，ART)等[24-25]，本文采用ART對層析方程組進行解算。ART避免了系數(shù)矩陣的求逆運算，具有抗噪性強的優(yōu)點，即使在觀測方程較少、信號分布不均勻和信號噪聲較大的不利條件下，也能解算出精度較高的結(jié)果[26-27]。

2 基于水汽垂直分布特征的自適應(yīng)非均勻指數(shù)分層方法

自適應(yīng)非均勻指數(shù)分層方法總體分為兩大步驟：首先，對層析區(qū)域進行大氣水汽的垂直分布建模，以該區(qū)域的先驗水汽廓線信息為依據(jù)獲得大氣水汽的垂直分布特征函數(shù)；然后，基于該水汽垂直分布特征函數(shù)，通過迭代的方式進行各層高度區(qū)間的計算，實現(xiàn)對任意層析區(qū)域進行任意層數(shù)劃分的自適應(yīng)分層。

2.1 水汽垂直分布建模

非均勻指數(shù)分層基于先驗水汽垂直分布特征函數(shù)，以水汽密度區(qū)間計算分層高度區(qū)間，在進行自適應(yīng)分層前需先對層析區(qū)域大氣水汽的垂直分布進行建模，以獲取高度與水汽密度之間的函數(shù)關(guān)系。

基于水汽密度隨高度增加呈指數(shù)遞減規(guī)律對大氣水汽的垂直分布特征進行建模[28]，高度h處的水汽密度表達為

(6)

式中，ρh表示高度h處的水汽密度，單位為g/m3；ρ0表示地表處水汽密度參數(shù)；Dv表示水汽垂直分布特征參數(shù)；H表示水汽標高，單位為m[20]。

由式(6)可知，某一高度h處的水汽密度值ρh(g/m3)可表示為高度h的函數(shù)，其中包含水汽密度參數(shù)ρ0和水汽垂直分布特征參數(shù)Dv兩個未知參數(shù)，通過對研究區(qū)域的先驗水汽密度數(shù)據(jù)進行擬合可獲得式(6)中的兩個未知參數(shù)。獲得先驗水汽密度數(shù)據(jù)的方法有多種，以探空氣象站的探空數(shù)據(jù)為例，探空數(shù)據(jù)具有高垂直分辨率的優(yōu)點，高精度的先驗探空水汽密度可通過原始探空資料中包含的溫度、露點溫度等觀測數(shù)據(jù)進行計算，水汽密度Dw計算公式為

(7)

式中，e表示水汽分壓，單位為hPa；Rv=461.5 J·kg-1·K-1；T表示溫度，單位為K。其中，e可由世界氣象組織推薦的Goff-Gratch公式計算

(8)

式中，T1=273.16 K；Td表示露點溫度，單位為K；e表示水汽分壓，單位為hPa。

上述部分詳細闡述了利用探空數(shù)據(jù)計算先驗水汽密度的公式和步驟，但由于探空數(shù)據(jù)的空間分辨率和時間分辨率較低，可能出現(xiàn)層析區(qū)域或?qū)游鰰r間段內(nèi)探空數(shù)據(jù)無法獲取的情況。針對該情況，利用ECMWF或美國國家環(huán)境預(yù)報中心(European Centre for Medium-Range Weather Forecasts，NCEP)提供的再分析資料中溫度和相對濕度等氣象數(shù)據(jù)同樣可以獲得可靠的先驗水汽密度數(shù)據(jù)[18,29]。

圖1為水汽垂直分布特征擬合圖，其中散點代表每一高度位置處所對應(yīng)的實測探空水汽密度，曲線代表最小二乘擬合后的函數(shù)曲線，設(shè)擬合后曲線的函數(shù)表達式為

圖1 水汽垂直分布特征擬合Fig.1 Vertical distribution characteristics of water vapor

(9)

式中，x表示高度，單位為m；y表示水汽密度，單位為g/m3；a和b分別表示地表處水汽密度參數(shù)ρ0和水汽垂直分布特征參數(shù)Dv擬合后的具體數(shù)值；H表示水汽標高，單位為m。式(9)反映了層析區(qū)域內(nèi)大氣水汽于垂直方向上的總體分布特征，為下一步自適應(yīng)分層提供了高度值與水汽密度值之間的具體函數(shù)關(guān)系。

2.2 自適應(yīng)非均勻指數(shù)分層

在獲得如式(9)所示的層析區(qū)域水汽垂直分布特征函數(shù)后，可通過如下方法進行自適應(yīng)垂直分層：將層析區(qū)域劃分為L層，基于文獻[4]的研究，垂直分層間隔需大于300 m[4]，首先將第1層高度區(qū)間設(shè)置為[0,300]，根據(jù)式(9)計算300 m處水汽密度為a·eb·(-300/H)(g/m3)，第1層水汽密度區(qū)間為[a,a·eb·(-300/H)]，現(xiàn)還需將剩余高度劃分為L-1個高度區(qū)間，計算剩余L-1層每層水汽密度間隔為(a·eb·(-300/H)-0.1)/(L-1)，單位為m，以該值計算各層預(yù)設(shè)水汽密度區(qū)間，得剩余L-1層的預(yù)設(shè)水汽密度區(qū)間表達式為

(10)

式中，a、b、H表示含義與式(9)相同；n表示層數(shù)(2≤n≤L)；L表示預(yù)計劃分總層次數(shù)。

將剩余L-1層對應(yīng)式(10)的預(yù)設(shè)水汽密度區(qū)間代入式(9)，反推計算出剩余L-1層的預(yù)設(shè)高度區(qū)間，其中第2層預(yù)設(shè)高度區(qū)間為

(11)

式中，a、b、H、L表示含義與式(10)相同。

若該層高度間隔大于300 m，則以此時的預(yù)設(shè)高度區(qū)間為結(jié)果，完成非均勻指數(shù)分層；若該層高度間隔小于300 m，則將該層高度區(qū)間繼續(xù)以最小高度間隔300 m設(shè)置為1層，例如第2層高度區(qū)間為[300,600]，密度區(qū)間為[a·eb·(-300/H),a·eb·(-600/H)]，剩余層數(shù)為L-2，每層密度間隔為(a·eb·(-600/H)-0.1)/(L-2)，以新的密度間隔計算新的剩余層預(yù)設(shè)密度區(qū)間，在進行第m次以300 m為高度間隔劃分最底層后，剩余層的預(yù)設(shè)密度區(qū)間可表示為

(12)

式中，a、b、H表示含義與式(11)相同；m表示以300 m為高度間隔進行劃分的層次數(shù)量；n表示剩余第m+1層至L層的層號(m+1≤n≤L)；L表示預(yù)計劃分總層次數(shù)。

每次迭代，將式(12)的剩余層預(yù)設(shè)密度區(qū)間代入式(9)反推新的預(yù)設(shè)高度區(qū)間，判斷其中最低層高度間隔是否大于300 m，最低層高度間隔表示為

(13)

式中，a、b、H、m、L含義與式(12)相同；Hi表示剩余L-m層中最底層高度間隔，單位為m，若Hi大于300 m則完成非均勻指數(shù)分層，反之則繼續(xù)以300 m為高度間隔劃分剩余L-m層中的最底層并以此類推。

上述的自適應(yīng)非均勻指數(shù)分層如圖2所示。ANEVS進行垂直分層時，先通過所需分層層數(shù)計算水汽密度區(qū)間，再通過水汽密度區(qū)間和先驗水汽垂直分布特征函數(shù)反推高度區(qū)間，通過迭代判斷并調(diào)整最底層高度間隔。該方法能根據(jù)層析區(qū)域的大氣水汽垂直分布特征和所需分層層數(shù)自動調(diào)整分層高度區(qū)間，同時自適應(yīng)地控制各層水汽密度差異在一個很小的范圍內(nèi)，做到細分底層區(qū)域密度不同的水汽、合并高層區(qū)域密度相近的水汽，保證了垂直分層的精度和真實性，且該方法適用于任意層析區(qū)域，即使在水汽垂直分布特征發(fā)生改變的情況下，也能迅速根據(jù)改變后的水汽垂直分布特征自動進行所需分層數(shù)的自適應(yīng)非均勻分層，具有極強的適用性。

圖2 自適應(yīng)非均勻指數(shù)分層 Fig.2 Adaptive non-uniform exponential vertical stratification

3 GNSS水汽層析試驗驗證

3.1 試驗區(qū)域及時段

本文以香港為試驗區(qū)域，圖3為試驗區(qū)域內(nèi)參考站及探空氣象站分布圖，其中實心圓圖形表示19個參考站所在位置，實心三角形圖形表示探空氣象站所在位置。本文利用19個連續(xù)運行參考站(continuously operating reference system，CORS)的GPS觀測數(shù)據(jù)進行水汽層析試驗，以京士柏探空氣象站(HKKP)的探空(radiosonde，RS)數(shù)據(jù)為參考，將HKKP位置天頂方向的層析水汽密度與探空水汽密度進行對比分析。

圖3 香港地區(qū)參考站及探空站地理位置Fig.3 Geographic distribution of the Hong Kong SatRef with GNSS reference stations and radiosonde

試驗時段選取2019年8月(31d)，試驗時段內(nèi)香港地區(qū)正值夏季，氣候炎熱，且由于受熱帶氣旋“韋帕”及“白鹿”的影響，共發(fā)生3次大規(guī)模降雨，因此香港2019年8月同時具備了暴雨和炎熱無雨的不同水汽情況，十分適宜進行水汽層析試驗。

3.2 試驗方案設(shè)置

觀測值SWD的獲取使用GAMIT10.7對GPS觀測數(shù)據(jù)進行處理，為消除多路徑效應(yīng)的影響，將衛(wèi)星截止高度角設(shè)置為15°；引入固體潮模型和海潮模型進行改正；映射函數(shù)模型設(shè)置為VMF1；引入BJFS、CHAN、USUD 3個外網(wǎng)輔助IGS站進行聯(lián)合解算以削弱天頂對流層延遲(zenith total delay，ZTD)的相關(guān)性。

層析方程組的解算采用ART技術(shù)。該技術(shù)中迭代初值的精度對層析結(jié)果的準確性具有重要影響，高精度的迭代初值可同時提高迭代速度和層析結(jié)果質(zhì)量。本文根據(jù)香港地區(qū)2010—2018年探空數(shù)據(jù)計算獲得的不同年份8月全部垂直水汽廓線均值作為迭代初值進行層析試驗[27]。

由于在試驗時段內(nèi)只能獲取每日UTC 00∶00時和UTC 12∶00時的探空數(shù)據(jù)作為參考，故每日只選取UTC 00∶00時和UTC 12∶00時兩個時段進行水汽層析試驗。本文以水汽密度等于0.1 g/m3的高度作為對流層層頂，通過試驗時段內(nèi)62次探空頂層高度的平均值確定試驗區(qū)域頂層高度，經(jīng)計算頂層高度平均值為10 770 m，標準差(standard deviation，STD)為935 m。水汽層析試驗區(qū)域水平方向范圍為22.204°N—22.564°N與113.844°E—114.384°E之間，水平分辨率為0.090°×0.090°[23]，垂直方向范圍為地面至10 770 m之間，垂直分辨率由分層方案決定。

方案A：該方案基于傳統(tǒng)均勻分層，均勻分層高度間隔在800 m時，可以獲得較好的層析結(jié)果[5]，因此該方案將10 770 m高的層析區(qū)域分為13層，使每層間隔為828.5 m，最接近800 m。頂層約束的先驗水汽密度以先驗探空數(shù)據(jù)為依據(jù)，由于頂層約束值作為層析方程組的觀測值直接參與三維水汽密度場的反演，不能采用試驗時段內(nèi)的探空數(shù)進行先驗水汽密度的計算，且大氣水汽存在周期性變化規(guī)律[30-33]，故采用試驗?zāi)攴萸耙荒?2018年)同時段(夏季8月)的探空數(shù)據(jù)計算先驗水汽密度，經(jīng)計算，水汽密度均值為0.23 g/m3，故將頂層約束設(shè)置為0.23 g/m3。

方案B：該方案基于ANEVS，分層數(shù)設(shè)置為與方案A相同的13層。大氣水汽的垂直分布特征函數(shù)采用先驗探空水汽密度數(shù)據(jù)進行擬合。由于探空氣象數(shù)據(jù)具有高精度、高垂直分辨率的優(yōu)點，同時，以探空站點處的水汽密度信息代表探空站周圍一定區(qū)域內(nèi)大氣水汽密度真值的可靠性已經(jīng)被證實[34-36]，故探空站的垂直水汽密度廓線常被用作層析區(qū)域內(nèi)水汽密度的參考真值，以輔助分析水汽層析的反演結(jié)果[7-8,37]。本文的層析試驗區(qū)域?qū)儆谛》秶鷮游鰠^(qū)域，基于高精度的先驗探空水汽密度建立的水汽垂直分布模型在該層析區(qū)域內(nèi)具備有效性。由于水汽垂直分布特征函數(shù)直接影響了觀測方程系數(shù)矩陣的數(shù)值，因此作為先驗信息的水汽垂直分布特征函數(shù)不能采用試驗時段內(nèi)的探空水汽密度數(shù)據(jù)進行擬合，基于大氣水汽的周期性變化特征[30-33]，以試驗?zāi)攴萸耙荒?2018年)同時段(夏季8月)的探空水汽密度作為先驗信息在時間域上具備代表性，擬合后的指數(shù)函數(shù)表達式為

(14)

式中，Dw表示水汽密度值，單位為g/m3；h表示高度，單位為m。探空數(shù)據(jù)擬合后的R2為0.977 4；均方根誤差(root mean squared error，RMSE)為0.919 g/m3?；谑?14)的水汽垂直分布特征函數(shù)對試驗區(qū)域運用ANEVS進行垂直分層后，將0～1500 m以300 m為高度間隔劃分5層，1500～10 770 m劃分為剩余8層，每層水汽密度間隔約為1.70 g/m3，第6層為1500～1838 m，高度間隔338 m，大于300 m，經(jīng)前一年同時段探空數(shù)據(jù)計算頂層水汽密度均值為0.87 g/m3，故將頂層約束設(shè)置為0.87 g/m3。

分別以方案A和方案B對試驗區(qū)域進行垂直分層建模，圖4為兩種分層方案的分層建模示意圖。

圖4 方案A及方案B垂直分層模型的垂直分辨率設(shè)置Fig.4 The vertical resolution of the vertical stratification model of project A and project B

4 GNSS水汽層析試驗結(jié)果分析

4.1 不同垂直分層方案的解算精度分析

圖5為試驗時段內(nèi)UTC 00∶00時與UTC 12∶00時層析歷元下兩種方案的單天RMSE值對比圖。其中方案A的RMSE數(shù)值最大出現(xiàn)在DOY 213 UTC 00∶00時(2.15 g/m3)，數(shù)值最小出現(xiàn)在DOY 214 UTC 12∶00時(0.70 g/m3)；方案B數(shù)值最大出現(xiàn)在DOY 233 UTC 00∶00時(1.83 g/m3)，數(shù)值最小出現(xiàn)在DOY 214 UTC 12∶00時(0.54 g/m3)。方案B有近50%的試驗RMSE小于1 g/m3，而方案A僅4組試驗RMSE小于1 g/m3，同時方案B有90%以上的RMSE小于方案A，故單天解精度對比中方案B明顯優(yōu)于方案A。

圖5 層析時段內(nèi)兩種方案層析結(jié)果的RMSE對比Fig.5 Comparison of the RMSE of the tomography results derived from 2 projects during the tomography periods

本文以高精度和高垂直分辨率的探空水汽密度廓線信息為參考值，對整個試驗數(shù)據(jù)集的層析結(jié)果進行整體精度評定，表1給出了兩種方案層析結(jié)果的RMSE和MAE(mean absolute error)。由表1可知，方案B的兩種統(tǒng)計量相比于方案A均具有更低數(shù)值，其中RMSE降低了0.401 g/m3，MAE降低了0.223 g/m3，說明ANEVS的整體解算精度高于均勻分層。

表1 兩種方案的層析結(jié)果與探空數(shù)據(jù)對比的平均RMSE和平均MAETab.1 Average RMSE and average MAE of the tomography results derived from 2 projects using radiosonde data as a reference (g/m3)

圖6為試驗時段內(nèi)，兩種方案層析水汽密度與探空數(shù)據(jù)之間的散點圖。如圖6所示，方案B相較于方案A具有更強的聚集性，同時方案B具有更高的R2、更低的RMSE，擬合線性函數(shù)的斜率也更接近1，表明方案B解算結(jié)果質(zhì)量優(yōu)于方案A。對比兩方案的散點圖還可發(fā)現(xiàn)，方案B散點在20～25 g/m3區(qū)間內(nèi)的聚集程度高于方案A，這是由于均勻分層于近地區(qū)域的水汽分層粗糙，底層間距大，水汽密度跳躍殘缺，衛(wèi)星信號利用率低，大幅降低了方案A底層的解算精度，導致底層解算精度與中間層解算精度相當，因此方案A的散點圖于20～25 g/m3區(qū)間內(nèi)的聚集程度十分接近5～20 g/m3區(qū)間，而方案B的底層水汽分層細致，信號利用率高，大幅提高了底層解算精度。下面本文將對兩種方案在不同天氣條件下和不同高度處的層析結(jié)果做進一步分析。

圖6 兩種方案層析反演水汽密度與探空水汽密度散點圖Fig.6 Scatter plots of water vapor density between radiosonde and tomography derived from 2 projects

4.2 不同天氣層析結(jié)果分析

圖7為層析水汽密度廓線與探空水汽密度對比圖，其中紅色曲線和藍色曲線分別表示方案A和方案B的垂直水汽廓線信息，綠色散點表示探空水汽密度信息。本文選取DOY 221 UTC 00∶00時和DOY 238 UTC 00∶00時的層析水汽廓線為例與探空數(shù)據(jù)進行驗證分析，其中DOY 221最高溫度達34.8℃，無降水活動發(fā)生，屬于高溫炎熱的無雨天氣；DOY 238由于受本月第二次強臺風“白鹿”的影響出現(xiàn)持續(xù)大暴雨，屬于暴雨天氣。由圖7可知，方案B與探空水汽密度具有更好的符合度，尤其在3000 m以下的高度區(qū)間內(nèi)，方案B的水汽廓線能更好地反映水汽的真實情況。此外，對比兩種方案在不同天氣條件下的垂直水汽廓線可以發(fā)現(xiàn)，方案A在近地高度區(qū)間內(nèi)僅能表現(xiàn)出大氣水汽的整體變化趨勢，其原因是在進行層析區(qū)域建模時，基于均勻分層的方案A在3000 m以下僅具有3層完整的分層區(qū)間，導致通過水汽層析技術(shù)反演獲得的水汽密度參數(shù)數(shù)量較少，無法表現(xiàn)出大氣水汽在惡劣天氣條件下的波動現(xiàn)象，而基于ANEVS的方案B在3000 m以下具有8層完整的分層區(qū)間，充足的水汽密度參數(shù)使其具備了反映近地空間大氣水汽波動現(xiàn)象的能力。

圖7 兩種方案的層析水汽廓線與探空水汽密度對比Fig.7 Comparison of tomographic water vapor profiles between 2 projects using the radiosonde as a reference

水汽密度廓線只具備表現(xiàn)天頂方向的垂直水汽密度信息的能力，并不能對三維水汽密度場進行展示，如圖8所示，本節(jié)基于方案A和方案B的層析反演結(jié)果，展示了DOY 238 UTC 00∶00時三維水汽密度場的可視化結(jié)果。方案A與方案B相比，方案A中下邊界處層與層之間的色階變化明顯，上邊界處的色階變化則較小，表明均勻分層反演的三維水汽密度信息具有下邊界處匱乏而上邊界處冗余的缺點，對比圖中方案B反演的三維水汽場可發(fā)現(xiàn)，ANEVS彌補了均勻分層的不足，其色階變化豐富，具有更精確的水汽密度信息反演能力。

圖8 兩種方案的三維水汽密度場可視化Fig.8 Visualization of three dimensional water vapor density field derived from 2 projects

表2為DOY 221和DOY 238不同方案層析水汽密度與探空水汽密度的相關(guān)系數(shù)。在無雨日(DOY 221)和降雨日(DOY 238)的對比中，兩種方案均在無雨天氣情況下具有更高的相關(guān)系數(shù)，說明惡劣天氣一定程度上影響了層析結(jié)果的質(zhì)量；在兩種方案之間的對比中，不同天氣條件下方案B的相關(guān)系數(shù)均高于方案A，表明ANEVS相比于傳統(tǒng)均勻分層具有更優(yōu)質(zhì)的大氣水汽反演效果；對比兩種天氣條件下不同方案間的相關(guān)系數(shù)之差可發(fā)現(xiàn)，在水汽不穩(wěn)定的惡劣天氣，方案B對于相關(guān)系數(shù)的提高量大于水汽穩(wěn)定的無雨天氣，說明方案B在惡劣天氣條件下能更好地優(yōu)化大氣水汽的反演效果。

4.3 不同高度層析結(jié)果分析

圖9給出了兩種方案的層析結(jié)果在每1層上與探空數(shù)據(jù)對比的RMSE。從圖中可看出，兩種分層方案的各層RMSE均呈現(xiàn)先增大再減小的變化規(guī)律，在方案B中所有分層的RMSE均小于1.5 g/m3，而方案A中1～5層的RMSE大于1.5 g/m3，其中第2層的RMSE大于2.5 g/m3，表明方案A于低海拔處的層析解算結(jié)果質(zhì)量較低。從兩種方案的各層RMSE對比中還可發(fā)現(xiàn)，方案A中9～11層的RMSE小于1 g/m3，而方案B中對應(yīng)層數(shù)的RMSE大于1 g/m3，產(chǎn)生此現(xiàn)象的原因是兩種方案在分層數(shù)相同的情況下，每1層所對應(yīng)的實際高度區(qū)間不同，分層層數(shù)越大，兩個方案之間的實際高度區(qū)間相差越大。由4.2節(jié)可知，分層數(shù)9～11層于方案A中實際對應(yīng)6628～9113.5 m，于方案B中其實際對應(yīng)2688～4982 m，而6628～9113.5 m區(qū)間內(nèi)的實際水汽含量相較于2688～4982 m區(qū)間內(nèi)的小很多，導致方案A的RMSE在9～11層小于方案B。因此，為了進一步對比兩種方案的垂直精度，圖9還給出了方案A和方案B在各高度處的RMSE對比圖。由圖9可知，方案A于4000 m以下的解算精度明顯低于方案B，這是由于方案A低海拔處分層粗糙，不符合大氣水汽的垂直分布特征，衛(wèi)星信號利用率低下，大幅降低了層析解算精度；4000 m以上兩種方案精度相當，即使方案A于高海拔處分層更精細，但由于衛(wèi)星信號分布稀疏、水汽稀薄，解算精度并未得到實質(zhì)性提高。方案B使用ANEVS由于考慮到水汽的垂直變化特征，以水汽密度區(qū)間為分層準則，該方法獲得的非均勻高度區(qū)間使各層水汽密度變化控制在一個較小的范圍內(nèi)，同時這種非均勻高度區(qū)間還適應(yīng)了衛(wèi)星信號呈倒錐形發(fā)散的特征，在保證了高海拔處解算絕對誤差較小的同時，大幅提高了低海拔處的解算精度。

圖9 兩種方案層析結(jié)果與探空數(shù)據(jù)對比的各分層RMSEFig.9 RMSE at each layer of the tomographic results derived from 2 projects using the radiosonde as a reference

大氣水汽整體表現(xiàn)出高層稀薄、底層充沛的特點，這種現(xiàn)象使高層水汽密度的絕對誤差相較于其他位置小很多，再加上傳統(tǒng)均勻分層將過多的層次劃分至絕對誤差本就較小的高層區(qū)域內(nèi)，因此使用與絕對誤差相關(guān)的統(tǒng)計量評定解算精度可能造成誤判，對此本文額外使用相對誤差(relative error，RE)對層析結(jié)果進行評定。相對誤差表示為絕對誤差與真值的比值，計算公式為

(15)

式中，RE表示相對誤差(%)；DT表示水汽層析反演水汽密度值，單位為g/m3；DRS表示利用探空數(shù)據(jù)計算得到的水汽密度值，單位為g/m3。

相對誤差能更好地考慮真值本身的數(shù)量級，用百分比的形式對誤差進行歸一化處理，因此在研究不同垂直分層方法之間解算精度時使用相對誤差具有更好的真實性和可比性。如圖10所示，兩種方案的各層平均相對誤差(mean relative error，MRE)均隨高度增加而增加，表明解算精度隨高度增加而遞減。在2500 m以下的高度區(qū)間內(nèi)，兩種方案的MRE均小于20%，解算精度較高，ANEVS的MRE在2000 m以下小于10%，在1000 m以下小于5%，精度可達傳統(tǒng)均勻分層的2倍以上，具有更高的解算精度；在2500～4500 m的高度區(qū)間內(nèi)，兩種方案的解算精度同時隨高度增加而降低，但方案B的精度仍稍優(yōu)于方案A。由圖10還可以發(fā)現(xiàn)，方案A只有1～7層MRE小于100%，8～13層MRE均大于100%，說明傳統(tǒng)均勻分層并未切實提高高層區(qū)域的解算精度；方案B僅最頂層MRE超過100%，倒數(shù)第2層MRE接近40%，其余層MRE數(shù)值更小，ANEVS細致的低層層次劃分適應(yīng)了大氣水汽的分布特征，充分高效地利用了低海拔密集的衛(wèi)星信號觀測值，水汽層析解算效率遠高于方案A。綜上，在以MRE為依據(jù)的解算精度對比上，ANEVS相比于傳統(tǒng)均勻分層有很大提高，該方法能將低海拔水汽進行細分，控制各層水汽密度差異，提高衛(wèi)星信號的利用率，將高海拔處密度相近的水汽歸并為同一層次，減少了高層水汽解算的多余工作，在自適應(yīng)地進行大氣水汽垂直分層的同時，克服了傳統(tǒng)均勻分層方案的缺點。

圖10 兩種方案層析結(jié)果與探空數(shù)據(jù)對比的各層MREFig.10 MRE at each layer of the tomographic results derived from 2 projects using the radiosonde as a reference

5 結(jié) 語

本文提出了一種基于大氣水汽垂直分布特征的自適應(yīng)非均勻指數(shù)垂直分層方法(ANEVS)。ANEVS以水汽密度區(qū)間和水汽垂直分布特征函數(shù)反推垂直分層高度區(qū)間，再通過迭代自適應(yīng)地調(diào)整高度區(qū)間，以縮小各層之間的水汽密度差異，使分層高度區(qū)間滿足實際水汽變化規(guī)律，實現(xiàn)對任意給定層析區(qū)域進行任意分層層數(shù)的自適應(yīng)非均勻分層，該非均勻分層方法在具備科學可靠分層依據(jù)的同時又具有方便快捷且自適應(yīng)的優(yōu)勢。本文利用香港CORS網(wǎng)2019年8月的GNSS觀測數(shù)據(jù)對比分析了ANEVS與傳統(tǒng)均勻分層的層析結(jié)果的質(zhì)量及精度，得到以下結(jié)論。

(1)ANEVS在層析結(jié)果的精度與質(zhì)量上均優(yōu)于傳統(tǒng)均勻分層。在散點圖的對比上，ANEVS表現(xiàn)出更強的聚集性，同時ANEVS散點圖具有更高的確定系數(shù)和更低的均方根誤差，ANEVS層析結(jié)果的均方根誤差和平均絕對誤差相較于傳統(tǒng)均勻分層分別降低0.401 g/m3和0.223 g/m3。

(2)在不同天氣條件下，ANEVS水汽廓線能更好地反映水汽的垂直擾動變化。不同天氣下ANEVS水汽密度與探空水汽密度的相關(guān)系數(shù)均高于傳統(tǒng)均勻分層，在暴雨天氣下兩種方法的相關(guān)系數(shù)之差大于無雨天氣，表明ANEVS在惡劣天氣條件下對大氣水汽反演效果的優(yōu)化更加顯著。

(3)在不同高度的層析結(jié)果對比中，ANEVS對低層高度區(qū)間內(nèi)層析結(jié)果精度的提高更加顯著。在分層數(shù)同為13層的條件下，傳統(tǒng)均勻分層的平均相對誤差有6層高于100%，僅1～7層小于50%，ANEVS的平均相對誤差僅最頂層1層高于100%，1～9層小于20%，ANEVS的均方根誤差于4000 m以下的高度區(qū)間內(nèi)顯著減小，表明ANEVS大幅提高了各層層析結(jié)果的解算質(zhì)量，其中低層層析結(jié)果質(zhì)量的提升尤為顯著。

筆者后續(xù)的研究與工作內(nèi)容主要包括兩方面：一方面是調(diào)整并擴大層析試驗區(qū)域，以便引入復(fù)數(shù)個探空氣象站的數(shù)據(jù)，對大氣水汽垂直分布特征模型的有效覆蓋范圍進行定量分析；另一方面，基于對先驗大氣水汽密度進行長時間跨度的統(tǒng)計學分析，研究并總結(jié)出一種能夠適用于不同月份的表示大氣水汽垂直分布特征的建模準則。

致謝：特別感謝香港地政總署測繪處提供的香港衛(wèi)星定位參考站網(wǎng)絡(luò)的GNSS觀測數(shù)據(jù)(https:∥www.geodetic.gov.hk/en/satref/satref.htm)，香港King 's Park天文臺提供的京士柏氣象站的無線電探空數(shù)據(jù)(http:∥weather.uwyo.edu/upperair/sounding.html)，麻省理工學院地球大氣和行星科學系提供的GAMIT/GLOBK軟件。