一種新的高精度全球對流層天頂延遲模型

黃良珂， 陳華， 劉立龍*， 姜衛平

1 桂林理工大學測繪地理信息學院, 桂林 541004 2 廣西空間信息與測繪重點實驗室, 桂林 541004 3 武漢大學衛星導航定位技術研究中心, 武漢 430079 4 武漢大學測繪學院， 武漢 430079

0 引言

對流層中蘊含著豐富的氣體，當電磁波信號穿過對流層時，會使這些電磁波信號產生延遲，稱為對流層延遲.對流層延遲是空間大地測量中最主要的誤差源之一，其也是制約高精度全球導航衛星系統(GNSS)導航定位的重要因素.在GNSS精密數據處理中，常采用對流層延遲經驗模型來估計出對流層天頂靜力學延遲(ZHD)或ZTD作為初始值，如果模型計算的ZHD或ZTD初始值不精確，最終也會影響GNSS定位的精度.由此可見，高精度的對流層延遲經驗模型顯得尤為重要，尤其在全球坐標參考框架建立等GNSS 應用中.因此，建立實時高精度全球對流層延遲模型具有重要的現實意義.

對流層延遲改正模型主要可分為兩大類，即需要實測氣象參數的模型和實時對流層延遲模型(非氣象參數模型).在依賴實測氣象參數的對流層延遲模型中，國際上比較經典的有Hopfield(1971)、Saastamoinen(1972)和Black(1978)等對流層延遲模型.由于難以獲取全球任意位置實時的實測氣象參數，從而限制了需要實測氣象參數的對流層延遲模型在GNSS導航和實時定位中的應用.為了滿足對流層延遲應用需求的實時性，利用大氣再分析資料構建區域或者全球的實時對流層延遲模型獲得廣泛關注，如UNB系列模型(Leandro et al.,2006,2008)和EGNOS模型(Penna et al.,2001)等.近年來，隨著GPS/GNSS導航定位技術的發展，諸多學者對實時對流層延遲改正模型的構建給予了極大的關注.如國外學者構建的TropGrid系列模型(Schüler,2014)和 GPT系列模型(B?hm et al.,2007,2015; Lagler et al.,2013).與此同時，中國學者在這方面也取得了豐碩的成果，Song等(2011)利用多年歐洲中尺度天氣預報中心(ECMWF)大氣再分析資料，構建了適用于中國區域的對流層天頂延遲模型(簡稱SHAO模型)，其在中國地區的改正精度相對于EGNOS模型提升顯著.黃良珂等(2014)使用數年IGS站的高精度ZTD數據對EGNOS模型在亞洲地區進行了單站精化，構建了亞洲地區的精化模型(SSIEGNOS模型)，SSIEGNOS模型能取得較好的單站ZTD改正效果，且在長期的ZTD預報中能保持穩定的性能.李薇等(2012)在分析ZTD時空特性基礎上，利用多年的美國環境預報中心(NCEP)再分析資料建立了一種全球格網ZTD新模型——IGGtrop模型，該模型在全球范圍內獲得了較好的精度.Li等(2015)對IGGtrop模型進行改進，發展了IGGtrop_ri(i=1,2,3)模型，新模型對赤道區域進行了簡化，并極大地減少了模型參數，使新模型更好地應用于北斗/GNSS的對流層延遲改正.姚宜斌等(2013)利用GGOS大氣中心提供的2°×2.5°(緯度×經度)分辨率ZTD格網數據，結合球諧函數構建了全球非氣象參數GZTD模型, GZTD模型表現出了良好的對流層修正效果.此外，姚宜斌等(姚宜斌等,2015; Yao et al.,2015)先后構建了高時間分辨率的GZTD-6h模型和ITG格網模型，新模型相比于GZTD模型其精度有一定的改善，并提高了時間分辨率.趙靜旸等(2014)利用ECMWF中心提供的ERA-Interim資料，構建了顧及垂直剖面函數的全球非氣象參數對流層延遲模型(SHAO-H模型)，該模型具有精度穩定、可適用于任意高度對流層延遲修正等優點.姚宜斌等(2016)聯合GGOS大氣格網ZTD和IGS中心提供的精密ZTD產品，以Delaunay三角網的形式存儲模型參數，進而構建了多源數據聯合的全球ZTD模型，新模型表現出全球適用和區域增強的優點.Sun等(2017)基于非線性假設，利用GGOS大氣ZTD格網產品構建了顧及ZTD年周期、半年周期、4個月變化和季度變化的全球ZTD格網模型(GZTDS模型)，并結合GZTDS和VMF1模型建立了GSTDS模型，結果表明GZTDS模型計算ZTD的性能與GPT2w模型相當.Li等(2018)分別采用指數函數和多項式函數去表征ZTD垂直剖面函數的平均值和振幅，最終構建了兩個非氣象參數的全球ZTD模型，即IGGtrop_SH和IGGtrop_rH模型，這兩個新模型均表現出良好的ZTD改正性能，IGGtrop_SH模型的性能略優于IGGtrop_rH模型.Sun等(2019)利用長期的ERA-Interim數據分析了對流層垂直遞減率的精細季節變化，在此基礎上構建了高水平分辨率的全球格網對流層模型(GTrop模型)，結果表明該模型較好地改善了在高海拔地區的ZHD估計精度.

盡管上述區域或全球對流層延遲模型均表現出各自的優越性，但仍存在模型方程未同時顧及高程和緯度因子的影響，且模型構建僅使用單一格網點數據等不足.相關研究表明ZTD與高程和緯度均存在顯著的相關性，為此，本文提出引入滑動窗口算法將全球剖分為大小一致的規則窗口，使用多年的全球GGOS大氣格網ZTD產品及對應的橢球高數據，構建全球每個窗口同時顧及高程、緯度和季節變化的ZTD模型，最終建立一種新的高精度全球ZTD格網模型(簡稱GGZTD模型)，并以未參與建模的全球GGOS大氣格網ZTD數據及全球IGS站精密ZTD產品為參考值，驗證GGZTD模型在全球的適用性.

1 ZTD時空特性分析與GGZTD模型構建

1.1 ZTD時空特性分析

相關文獻研究表明ZTD時間序列主要表現為年周期和半年周期變化，同時在部分地區也表現出一定的日周期變化(姚宜斌等，2013)，為了進一步驗證ZTD值的年周期、半年周期和日周期振幅在全球的分布，利用2010—2015年水平分辨率為2°×2.5°、時間分辨率為6 h的全球GGOS大氣格網ZTD數據，采用式(1)來計算出全球各格網點ZTD值不同周期的振幅值，結果如圖1所示.

(1)

式中，A0表示ZTD年均值，(A1,A2)、(A3，A4)和(A5，A6)分別表示ZTD的年周期、半年周期和日周期的系數，doy表示年積日，hod表示UTC時間.

圖1 GGOS格網ZTD數據計算全球ZTD平均值、年周期、半年周期和日周期分布Fig.1 Distributions of the mean value, annual, semi-annual and diurnal amplitudes of ZTD derived from GGOS Atmosphere gridded ZTD products over globe

由圖1可知，ZTD平均值的分布與高程和緯度存在較大關系，其在青藏高原、南極和格陵蘭等高海拔地區的值較小，在海洋地區其平均值較大；對于ZTD年周期振幅，其在北半球的年周期振幅值高于其在南半球的值，在太平洋西部、亞洲東部和南部、北美洲東南部、非洲、澳大利亞北部等地區表現出顯著的年周期振幅值，在南半球的中高緯度地區其值相對較?。粚τ赯TD的半年周期振幅，其在北半球的值仍然大于南半球，在中國東北部、亞洲東南部和低緯度的部分地區存在顯著的半年周期振幅值，在南半球的中高緯度地區仍表現出相對較小的值；對于ZTD的日周期振幅，其在低緯度的部分地區表現出相對較大的日周期振幅值，尤其在南美洲中部和非洲南部地區存在4～8 mm的日周期振幅值，但是在全球大部分地區ZTD的日周期變化較小，其值低于2 mm.由于在低緯度地區，ZTD的年周期和半年周期相對較為顯著，此外，該地區的氣候較為復雜，其引起的ZTD的變化也較為劇烈，盡管在該地區的部分區域存在4～8 mm的日周期變化，那么在ZTD建模時，這些相對較小的ZTD日周期變化容易被ZTD年周期和半年周期的模型化殘差及其他非季節性信號淹沒.因此，為了優化ZTD建模的模型參數，提高模型的計算效率，本文模型構建時對ZTD的季節變化只考慮其年周期和半年周期.

對于ZTD在高程上的變化關系方面，常用的ZTD垂直剖面函數主要有二次多項式和負指數函數.Song等(2011)指出在10 km以下用二次多項式來表達ZTD的垂直剖面可獲得較好的ZTD高程改正效果.趙靜旸等(2014)提出了ZTD的分段垂直剖面函數，在17 km以下的大氣層范圍內利用多項式來表達ZTD的垂直剖面，17 km以上采用負指數函數來進行表達，該分段ZTD垂直剖面模型獲得了良好的ZTD高程改正效果.由于本文構建全球ZTD模型時采用的是全球GGOS大氣格網ZTD產品及對應的地表橢球高數據，為了進一步探究地表至17 km高度范圍內ZTD在高程上的變化，本文選取了全球具有代表性的MERRA-2分層大氣資料的4個格網點數據，利用積分法計算出4個格網點(地表至17 km高度范圍內)2015年1月1日0時刻的ZTD分層剖面信息，并采用二次多項式對其進行擬合，結果如圖2所示.圖2表明在地表至17 km高度范圍內，利用二次多項式可較好地反映ZTD值在高程上的變化.同時，為了進一步分析ZTD在緯度上的變化，選取經度為170°W和170°E的所有格網點數據，針對每個格網點均從500 m開始積分到MERRA-2再分析資料最頂層，進而分析ZTD隨緯度的變化關系，結果如圖3所示.由圖3可知，在一定緯度范圍內ZTD在緯度上的變化可表示為近似線性關系.因此，在本文ZTD建模中將采用二次多項式來描述ZTD的垂直剖面改正，ZTD在緯度上變化采用線性函數.

圖2 ZTD在高程上變化關系Fig.2 ZTD changes with altitude

圖3 ZTD在緯度上變化關系Fig.3 ZTD changes with latitude

1.2 GGZTD模型建立

針對當前全球對流層經驗模型的模型方程未同時顧及高程和緯度的變化及模型構建僅使用單一格網點數據等問題，結合上述對ZTD的時空特性分析，引入滑動窗口算法來解決以上問題.在此，先對滑動窗口的算法做如下簡述.

首先以平面分辨率為2°×2.5°(緯度×經度)進行全球格網剖分，得到與GGOS大氣格網產品相同大小水平分辨率的格網.滑動窗口算法的關鍵需要確定其窗口的大小及其滑動步長，滑動窗口大小的確定需顧及全球窗口剖分個數的整數性、窗口的連續性及窗口內模型參數的可求解性等原則.根據以上原則，本文以3行3列(4°×5°的區域范圍)為一個滑動窗口大小來舉例說明滑動窗口算法，其流程如圖4所示.具體過程為：首先利用格網左上角第一個窗口N1內(每個框表示一個滑動窗口大小)的數據求出其窗口內的相關模型參數，并將其作為窗口N1中心格網點(框內黑點)的結果；然后將窗口向緯度東向移動2個格網點，求解新窗口N2內的相關模型參數，將其作為窗口N2中心格網點的結果，依次類推求出這一緯度上所有窗口內的相關模型參數；然后窗口移動到下一緯度(向下移動兩個格網點)，以同樣的方法求出該緯度所有窗口內的相關模型參數，以此類推直到求出全球所有窗口內的相關模型參數.最終獲得全球所有窗口內的相關模型參數，并將其作為各自窗口中心格網點的結果，最后將全球所有窗口中心格網點組建新的全球格網，如圖4中的黑點和虛線所示.

圖4 滑動窗口算法示意圖Fig.4 The schematic of sliding window algorithm

本文僅探索使用固定的窗口大小(4°×5°)和滑動步長來構建全球ZTD模型，即每個窗口包含了9個格網點，針對每一個窗口，GGZTD模型在時間尺度上主要考慮ZTD的年周期和半年周期特性，在空間尺度上同時顧及ZTD在緯度和高程上的變化，由此可確定全球每個窗口內的ZTD模型方程如下：

ZTD(φ,h,doy)=f1(φ,h)+f2(doy),

(2)

其中φ,h分別表示緯度和高程，doy表示年積日，f1(φ,h)表示緯度和高程對ZTD的共同影響，f2(doy)表示ZTD的季節變化，其分別可表示為

f1(φ,h)=a0+a1φ+a2h+a3h2,

(3)

(4)

最終可將式(2)改寫為

ZTD(φT,λT,h,doy)=C0(φi,λi)+C1(φi,λi)·φT+C2(φi,λi)·h+C3(φi,λi)·h2

(5)

式中，(φi,λi)表示第i個窗口的幾何中心的緯度和經度，單位均為度；φT,λT和h分別表示目標點的緯度(單位是°)、經度(單位是°)和高程(單位是m).C0表示ZTD的年均值，C1表示緯度改正系數，C2和C3表示高程改正系數，(C4，C5)和(C6，C7)分別表示ZTD的年周期和半年周期系數.利用2008—2015年水平分辨率為2°×2.5°的全球GGOS大氣格網ZTD地表數據和對應的全球橢球高數據，根據式(5)采用最小二乘算法即可求解出全球每個窗口的模型系數，最終GGZTD模型的模型系數以4°×5°的格網形式存儲.GGZTD模型在全球包含了3240個格網點.

GGZTD模型的使用非常簡便，其使用過程如下：(1)用戶提供doy和目標點位置信息，根據目標點位置信息查找與目標點最近的模型參數格網點；(2)根據獲得的最近的模型參數格網點的模型參數，并利用式(5)即可計算出目標點的ZTD信息.

2 GGZTD 模型精度驗證

為了驗證 GGZTD 模型在全球的精度和適用性，選取未參與建模的全GGOS 大氣格網 ZTD 產品和全球 IGS 中心提供的 GNSS ZTD 產品為參考值來驗證 GGZTD 模型的精度，并與目前精度標稱最優的 GPT2w 模型(B?hm et al., 2015)和全球廣泛使用的 UNB3m 模型(Leandro et al., 2008)進行精度對比，由于GPT2w模型可提供水平分辨率分別為1°×1°和5°×5°的格網參數，便于后續描述，分別將其簡寫為GPT2w-1和GPT2w-5.同時，以 bias(偏差)和 RMS 誤差作為精度指標來評判模型的精度.

2.1 與全球 GGOS 大氣格網 ZTD 產品比較

以GGOS中心提供的未參與建模的2016年時間分辨率為6 h(UTC 00∶00，06∶00，12∶00和18∶00)的全球GGOS大氣格網ZTD地表產品為參考值，同時計算出GGZTD、GPT2w-1、GPT2w-5和UNB3m模型在這些格網點處對應時間的ZTD值，由此可獲得這些模型計算ZTD的bias和RMS值，結果如表1、圖5和圖6所示.

表1 2016 年全球 GGOS 大氣格網 ZTD 產品驗證 GGZTD、GPT2w-5、GPT2w-1和UNB3m模型的精度統計(單位：cm)Table 1 Error statistics of GGZTD、GPT2w-5、GPT2w-1 and UNB3m models tested using GGOS Atmosphere gridded ZTD data in 2016 (Unit: cm)

由表1可知，GGZTD模型表現出較小的bias值和最小的RMS，平均bias值和RMS分別為-0.23 cm和3.58 cm，其bias值和RMS變化范圍分別為-5.06 cm到2.88 cm和0.85 cm到8.00 cm；GPT2w-5表現出了最大的bias和RMS值，分別為24.77 cm和27.55 cm，GPT2w-1的性能相對較為穩定且精度優于GPT2w-5，由此表明GPT2w模型通過提高其模型參數的空間分辨率對其精度具有一定的改善；而UNB3m模型表現出最大的平均bias和RMS值分別為3.08 cm和6.19 cm，說明UNB3m模型在全球存在顯著的系統偏差，主要原因是UNB3m模型提供的模型參數的空間分辨率較低.在RMS方面，相對于GPT2w-1、GPT2w-5和UNB3m模型，GGZTD模型在全球計算ZTD的精度分別提升了0.8 mm、1.3 mm和26.1 mm(42%).由此可見，GGZTD模型的精度相對于UNB3m模型具有顯著的提高，同時相對于目前性能優異的GPT2w模型也具有一定的改善.

圖5表明UNB3m模型在南半球高緯度地區(緯度大于60°)、印度洋、太平洋及亞洲西南部和非洲南部的部分地區出現較大的正bias值，在低緯度(30°S—30°N)的部分地區存在顯著的負bias值，其在北半球的bias值明顯小于南半球.相對于UNB3m模型，GPT2w-5和GPT2w-1均表現出相對較小的bias值，但是在太平洋東部和南半球高緯度的部分地區仍然存在相對較大的正bias值，尤其是GPT2w-5在北美洲南部和印度洋東部的部分地區存在較大的正bias值且在該地區出現最大的bias值，而在低緯度地區的大西洋、南美洲北部和非洲中部及太平洋南部的部分地區存在相對較大的負bias值，而GGZTD在全球范圍內均表現相對較小的bias值，其僅在低緯度的小部分區域存在相對較大的負bias值.

圖5 2016 年全球 GGOS 大氣格網 ZTD 產品驗證不同模型的bias值全球分布Fig.5 Distributions of the bias validated using GGOS Atmosphere gridded ZTD data for all models in 2016 over globe

由圖6可知，UNB3m模型在整個南半球區域及北半球低緯度地區存在較大的RMS，尤其在南半球高緯度地區、印度洋東部和太平洋部分地區，其在南半球的精度也明顯低于北半球，主要原因是受UNB3m模型的南北半球對稱假設的影響(Leandro et al., 2008)，此外，UNB3m模型在低緯度地區未考慮模型參數的年際變化；相對于UNB3m模型，GPT2w-5、GPT2w-1和GGZTD模型在全球范圍內表現出較小的RMS，但是在低緯度的部分地區(如大西洋、非洲、太平洋東部地區、中國東南部和澳大利亞)表現出相對較大的RMS，尤其是GPT2w-5在北美洲南部的部分地區，由于在這些地區受深對流、熱帶海洋氣候和熱帶雨林氣候等復雜氣候系統的綜合影響，進而引起ZTD的劇烈變化，因此難以對ZTD進行精確模型化；此外，相對于其他地區，低緯度地區的ZTD具有較為顯著的日周期變化(姚宜斌等, 2016)，而這三種模型均未顧及ZTD的日周期變化，從而導致所有模型在該地區均存在較大的誤差.盡管如此，相對于GPT2w模型計算的ZTD信息，GGZTD模型在太平洋東部、大西洋、北美洲南部和南極的部分地區仍具有一定的精度提升.因此，相對于UNB3m模型，GGZTD模型在全球估計ZTD的精度具有顯著的提高；而相對于目前性能優異的GPT2w模型，GGZTD模型在部分地區仍具有一定的精度改善.總之，相對于其他兩種模型，GGZTD模型在全球估計ZTD表現出穩定和優異的性能.

圖6 2016年全球GGOS大氣格網ZTD產品驗證不同模型的RMS全球分布Fig.6 Distributions of the RMS error validated using GGOS Atmosphere gridded ZTD data for all models in 2016 over globe

2.2 與全球IGS站ZTD產品比較

IGS中心能提供全球IGS站的精密ZTD產品，其時間分辨率為5 min，精度優于5 mm(Byun and Bar-Sever, 2009)，為了進一步驗證GGZTD模型在全球的精度，利用2016年全球IGS站精密ZTD產品作為參考值對其進行精度檢驗.由于部分IGS站的數據缺失較為嚴重，為了保證驗證的可靠性，選取了數據量至少含有120天的IGS站參與模型檢驗，最終在全球選取了2016年的316個IGS站，其點位分布如圖7所示.以2016年全球316個IGS站精密ZTD產品為參考值，同時計算出GGZTD模型、UNB3m模型和GPT2w模型在這些IGS站的ZTD值，最終獲得這些模型的bias和RMS值的統計結果，如表2、圖8和圖9所示.

圖7 全球316個IGS站的點位分布Fig.7 The distribution of 316 IGS stations over globe

表2 2016年全球316個IGS站ZTD產品驗證GGZTD、GPT2w-5、GPT2w-1和UNB3m模型的精度統計(單位：cm)Table 2 Error statistics of GGZTD、GPT2w-5、GPT2w-1 and UNB3m models tested using 316 IGS ZTD data in 2016 over globe (Unit: cm)

圖8 2016年全球IGS站ZTD產品驗證不同模型的bias值全球分布Fig.8 The distributions of the bias validated using IGS ZTD data in 2016 over globe

圖9 2016年全球IGS站ZTD產品驗證不同模型的RMS全球分布Fig.9 The distributions of the RMS error validated using IGS ZTD data in 2016 over globe

由表2可知，UNB3m模型表現出最大的平均bias和RMS值，其分別為0.66 cm和5.15 cm，與GGOS大氣ZTD格網產品驗證結果一致.GPT2w-5、 GPT2w-1和GGZTD模型均表現出負bias值，說明這些模型估計值偏??；相對于UNB3m模型，GPT2w-5和GPT2w-1表現出相對較小的平均RMS值分別為3.75 cm和3.69 cm，GPT2w-1相對于GPT2w-5其精度提升了0.6 mm，進一步說明GPT2w模型通過改善模型參數的空間分辨率能在一定程度上提高其ZTD估計的精度.而GGZTD模型在所有模型中表現出最優的性能，其RMS變化范圍為1.68～7.44 cm，平均值為 3.62 cm，相對于UNB3m模型、GPT2w-5和GPT2w-1，其在全球估計ZTD的精度分別提升了15.3 mm(30%)、1.3 mm和0.7 mm.進一步表明GGZTD模型在全球范圍內的精度相對于UNB3m模型具有顯著的提升，相對于性能優異的GPT2w模型仍具有一定的改善.

由圖8可以看出，UNB3m模型在南半球中高緯度地區、澳大利亞、非洲北部、北美洲東部和西南部地區存在顯著的正bias值，在低緯度、亞洲和歐洲的部分地區存在顯著的負bias值，在南半球的bias值要高于北半球.相對于UNB3m模型，GPT2w-5和GPT2w-1在全球范圍內表現出相對較小的bias值，只有在臨近赤道的PNGM站出現了較大的負bias值，可能受赤道附近復雜的氣候影響.而GGZTD在全球范圍內均表現出穩定和相對較小的bias值，且在南北半球的bias值分布較為均勻.圖9表明，在RMS方面，UNB3m模型在南半球和低緯度地區仍表現出較大的RMS值，尤其在南極地區，且在南半球的RMS值明顯大于北半球，這與GGOS大氣格網ZTD產品驗證結果一致，進一步說明UNB3m模型的南北半球假設不合理、模型參數空間分辨率偏低和低緯度地區未顧及年際參數變化.相對于UNB3m模型，GPT2w模型和GGZTD模型在全球均展現出相對較小且穩定的RMS，但是其在低緯度地區的太平洋和印度洋的GNSS測站中存在相對較大的RMS，由于在這些地區受熱帶海洋氣候和熱帶雨林氣候等復雜氣候系統的綜合影響，從而在ZTD信息中引入了非季節性變化信號，導致了ZTD的異常擾動，因此難以對ZTD進行精確模型化.此外，相對于中高緯度地區，低緯度地區ZTD也存在較為顯著的日周期變化，然而GPT2w模型和GGZTD模型未顧及日周期變化，從而在模型中引入了未模型化的誤差.盡管如此，相對于GPT2w模型，GGZTD模型在非洲東部、中國東部、歐洲和南極地區的部分GNSS測站中仍表現出相對優越的性能.

此外，為了進一步分析這些模型的精度分布情況，對其bias和RMS結果進行了精度統計，結果如圖10和圖11所示.圖10表明，從誤差分布直方圖可看出UNB3m模型具有較小的平均正bias值，但是仍然表現出較大數量的正負bias值的個數，其較大正bias值的個數大于較大負bias值的個數，說明UNB3m模型在全球大部分區域估計的ZTD值大于IGS-ZTD值.GPT2w-5、GPT2w-1和GGZTD模型均表現出較小的平均負bias值，相對于GPT2w-5和GPT2w-1，GGZTD模型的bias值分布更為集中，而GPT2w-5相對于GPT2w-1的bias分布則較為分散，由此說明GGZTD模型比GPT2w模型和UNB3m模型更具穩定性.由圖11可知，對于UNB3m模型，其RMS在3～4 cm的范圍占據了最大的比例為34%，而小于3 cm的比例僅為10%，大于6 cm的比例高達31%.相對于UNB3m模型，GPT2w-5和GPT2w-1的RMS分布具有明顯的改善，GPT2w-5的RMS在低于3 cm的比例中占據了29%(GPT2w-1為32%)，而大于6 cm的RMS比例低至6%(GPT2w-1為5%)，其在低于4 cm的比例分別為68%和71%.在所有模型中，GGZTD模型在小于3 cm和小于4 cm的RMS分布中均占據了最大比例，其比例分別高達為35%和73%，此外，其在大于6 cm的RMS比例中具有最小值，低至4%.在低于3 cm的RMS比例中，相對于UNB3m模型、GPT2w-5和GPT2w-1，GGZTD模型分別提高了25%、6%和3%，在低于4 cm的RMS比例中，其分別提高了29%、5%和2%.由此進一步表明相對于UNB3m模型和GPT2w模型，GGZTD模型在全球范圍內具有穩定的ZTD估計性能及最優的精度.

圖10 2016年全球IGS站ZTD產品驗證不同模型的bias值直方分布圖Fig.10 Histogram of bias for different models tested with IGS ZTD data in 2016 over globe

圖11 2016年全球IGS站ZTD產品驗證不同模型的RMS比例分布Fig.11 Proportional distribution of RMS values for different models tested with IGS ZTD data in 2016 over globe

為了分析GGZTD模型、GPT2w模型和UNB3m模型的季節變化性能，在全球南北半球的高、中和低緯度地區分別挑選了數據較為齊全的6個具有代表性的IGS站，即KELY站、BJFS站、BJCO站、PNGM站、DUND站和SYOG站，用于參與模型的季節變化性能分析，圖12給出了所有模型估計的對應測站的ZTD日均bias時序圖.

圖12 不同模型估計的2016年ZTD值在6個IGS站的bias時間序列(括號內表示每個測站的緯度、經度和高程)Fig.12 Time series of ZTD bias for different models compared to the IGS data in 2016(the brackets indicate the latitude, longitude and altitude of each station)

由圖12可知，位于北半球高緯度地區的KELY站，所有模型在春季表現為顯著的負bias值，在其他季節的變化不明顯；對于北半球中緯度地區的BJFS站，所有模型在夏季均表現出較大的bias值(以負bias值為主)，說明北半球的夏季水汽變化較為劇烈，導致所有模型的ZTD估計值低于其實測值，在其他季節里UNB3m模型存在顯著的正bias值，而GGZTD模型和GPT2w模型的變化不明顯，且在春季和冬季表現出較小的bias值；對于低緯度地區的BJCO站和PNGM站，UNB3m模型在全年均表現為顯著的負bias值，除了在BJCO站的冬季表現為顯著的正bias值外.對于南半球中緯度地區的DUND站，所有模型在春季和冬季均表現出較大的bias值，尤其是UNB3m模型，在夏季和秋季則表現為相對平穩的bias值，可能是受南半球中緯度地區冬季和春季的水汽變化劇烈的影響；對于南半球高緯度地區的SYOG站，UNB3m模型在全年均表現出較大的正bias值，尤其在春季和冬季，GPT2w-5在春季和夏季表現出相對明顯的正bias值，而GGZTD模型和GPT2w-1在全年均體現出相對較為穩定的bias值，相對于GPT2w-1, GGZTD模型在秋季表現出更為穩定和較小的bias值.由此表明，UNB3m模型在全球范圍內估計ZTD的bias值存在明顯季節變化，進一步說明UNB3m模型的模型參數空間分辨率低下以及南北半球對稱的假設不合理，導致了UNB3m模型在南半球使用存在較大的偏差，尤其在南半球低緯度地區，而GGZTD模型和GPT2w模型估計ZTD的bias值存在相對較小的季節變化，相對于GPT2w模型，GGZTD模型在低緯度地區估計ZTD值具有顯著的精度改善，其能更好地捕捉ZTD的季節變化.

為了進一步分析所有模型估計ZTD的bias值和RMS在高程上的變化關系，對全球316個IGS站按照500 m高程間隔進行分類，最終可分為低于500 m，500～1000 m，1000～1500 m、1500～2000 m和大于2500 m等5個高程范圍，統計不同模型在每個高程范圍內的bias值和RMS，結果如圖13所示.

圖13 2016年全球IGS站ZTD產品驗證各模型的bias和RMS在不同高程范圍內的變化Fig.13 Bias and RMS error statistics for different models tested with ZTD from IGS products in different altitude range in 2016

圖13表明，UNB3m模型在海拔大于2000 m時存在較大的負bias值，GPT2w-5、GPT2w-1和GGZTD模型在低于2000 m的高程范圍內均體現出相對明顯的負bias值，所有模型的bias值在高程上無明顯的變化關系，相對于GPT2w模型，GGZTD模型在1500～2000 m高程范圍內具有較小的絕對bias值.UNB3m模型的RMS在高程上未發現明顯的變化關系，盡管如此，其在低于500 m和大于2000 m的高程范圍內存在較大的RMS，但是總體上，GPT2w-5、GPT2w-1和GGZTD模型的RMS在高程上均存在隨著高程的增大而減小.此外，位于低緯度地區的IGS站其高程大部分低于500 m，導致低于500 m高程范圍的RMS較大，且隨著高程的增加水汽的含量及其變化也相應遞減，從而更易于對ZTD模型化.

如前所述，緯度也是影響ZTD的主要變化因子.為了進一步分析不同模型估計ZTD的bias值和RMS在緯度上的變化關系，對全球316個IGS站按照15°的緯度間隔進行分類，統計不同模型在每個緯度范圍內的bias值和RMS，結果如圖14所示.

由圖14可知，UNB3m模型在南半球15°以上范圍存在較大的正bias值，尤其在南半球高緯度地區，而GPT2w-5、GPT2w-1和GGZTD模型在全球各緯度范圍內均體現出相對較小的負bias值，相對于GPT2w模型，GGZTD模型在15°S—15°N緯度范圍內表現出相對較小的絕對bias值.在RMS方面，UNB3m模型在南半球的緯度范圍內的值明顯大于北半球，尤其在南半球的高緯度地區，此外，在低緯度地區也存在相對較大的RMS.對于GPT2w-5、GPT2w-1和GGZTD模型，總體上其RMS均存在以赤道向兩極減小的趨勢.相對于UNB3m模型，GGZTD模型在整個南半球和北半球的低緯度地區具有顯著的精度改善，尤其在南半球；相對于GPT2w模型，GGZTD模型在南半球緯度大于75°區域和北半球低緯度地區仍具有一定的精度提升，尤其在南半球緯度大于75°區域.由此表明GGZTD模型在所有對比的模型中具有穩定和優越的性能，尤其在南半球高緯度地區和全球低緯度地區.

圖14 2016年全球IGS站ZTD產品驗證各模型的bias和RMS在不同緯度范圍內的變化Fig.14 Bias and RMS error statistics for different models tested with ZTD from IGS products in different latitude range in 2016

3 結論

針對已有對流層延遲模型的模型方程未同時顧及高程、緯度和季節變化以及模型構建使用單一格網數據等不足，提出了基于滑動窗口算法的全球ZTD模型構建新方法，最終建立了一種新的高精度全球ZTD模型(GGZTD).同時，聯合未參與建模的全球GGOS大氣格網ZTD產品和全球IGS站ZTD產品評估了GGZTD模型在全球的精度和適用性.結果表明：(1)以全球GGOS大氣格網ZTD產品和全球IGS站精密ZTD產品為參考值，相對于UNB3m模型和GPT2w模型計算的ZTD信息，GGZTD模型在全球范圍內表現出了最優的精度和穩定的性能，尤其是相對于UNB3m模型，GGZTD模型體現出了顯著的精度提升，即便是相對于精度優異的GPT2w模型，其在全球的部分地區也表現出了一定的精度改善.此外，在季節、高程和緯度變化方面，GGZTD模型仍體現出了最優的性能.(2)GGZTD模型具有3240個格網點，每個格網點包含了8個模型參數.相對于GPT2w-1計算ZTD時所需的模型參數，GGZTD模型的模型參數減少了99%，而相對于GPT2w-5, 其模型參數減少了69%.因此，GGZTD模型在全球計算ZTD時只需依賴相對較少的模型參數，極大地提升了模型的計算效率，進而增強了模型在全球的適用性.總之，GGZTD模型僅需要輸入位置和時間參數則能在全球獲得較高精度和穩定的ZTD信息，且只依賴于相對較少的模型參數，因此其在全球坐標參考框架建立等GNSS應用和導航中具有重要的應用.由于本文在構建GGZTD模型時只考慮了固定的窗口大小和滑動步長，在下一步工作中，將采用高分辨率的ERA5資料研究不同窗口大小和步長的選取對全球ZTD建模的影響.

致謝感謝GGOS Atmosphere機構和IGS中心提供實驗數據.