覃澤穎，周淼

(1.桂林好測(cè)信息科技有限公司，廣西 桂林 541000； 2.桂林市測(cè)繪研究院，廣西 桂林 541000)

1 引 言

在全球?qū)Ш叫l(wèi)星系統(tǒng)(GNSS)中，衛(wèi)星信號(hào)穿過(guò)大氣層時(shí)會(huì)受到對(duì)流層的折射而產(chǎn)生彎曲和延遲，延遲量(即對(duì)流層延遲)從天頂方向到地平方向約為 2 m～ 20 m，因此對(duì)流層延遲是GNSS導(dǎo)航定位中的重要誤差源之一[1]。對(duì)于對(duì)流層延遲值，可利用模型修正法和氣象數(shù)值預(yù)報(bào)資料計(jì)算得到，對(duì)于后者廣為使用的有歐洲中尺度預(yù)報(bào)中心(ECMWF)和美國(guó)國(guó)家環(huán)境預(yù)報(bào)中心(NCEP)提供的再分析資料。目前，國(guó)際上已有諸多學(xué)者利用ECMWF或NCEP再分析資料建立了區(qū)域或全球的對(duì)流層延遲模型[2～4]。同時(shí)，國(guó)內(nèi)學(xué)者也利用大氣再分析資料建立了豐富的區(qū)域或全球?qū)α鲗友舆t改正模型[5～8]。同時(shí)，在利用大氣再分析資料計(jì)算對(duì)流層延遲信息在區(qū)域或全球范圍的精度評(píng)估方面也開(kāi)展了相關(guān)研究。文獻(xiàn)[9]利用2013年全球365個(gè)IGS站評(píng)估GGOS ZTD在全球范圍內(nèi)的精度，其年均偏差和RMSE分別為 -0.5 cm和 1.73 cm。文獻(xiàn)[10]利用中國(guó)區(qū)域7個(gè)IGS站評(píng)估ERA-Interim資料計(jì)算ZTD的年均偏差的絕對(duì)值和RMSE分別優(yōu)于 1 cm和 2 cm，且無(wú)明顯的季節(jié)變化。文獻(xiàn)[11]利用2015年237個(gè)陸態(tài)網(wǎng)GNSS測(cè)站評(píng)估ECMWF地表資料計(jì)算ZTD的年均RMSE誤差為 3.14 cm。文獻(xiàn)[12]利用2017年26個(gè)陸態(tài)網(wǎng)GNSS測(cè)站實(shí)測(cè)ZTD對(duì)ERA5和ERA-Interim資料計(jì)算ZTD的精度進(jìn)行評(píng)估，結(jié)果表明ERA5相較于ERA-Interim精度有明顯提升。此外，文獻(xiàn)[13]則是利用219個(gè)陸態(tài)網(wǎng)GNSS測(cè)站評(píng)估了中國(guó)區(qū)域ERA5地表及分層資料計(jì)算ZTD的精度，結(jié)果表明利用分層資料計(jì)算ZTD具有更高的精度。

一般情況下，大氣再分析資料在應(yīng)用前需要采用獨(dú)立的觀測(cè)值來(lái)對(duì)其進(jìn)行精度評(píng)價(jià)。而MERRA-2和ERA5分別是由美國(guó)NASA和歐洲中尺度預(yù)報(bào)中心提供的最新一代大氣再分析資料，具有極高的時(shí)空分辨率，由于目前尚無(wú)文獻(xiàn)對(duì)MERRA-2和ERA5再分析資料在廣西區(qū)域計(jì)算ZTD和ZWD的精度進(jìn)行評(píng)估。為此，本文聯(lián)合2017年廣西區(qū)域陸態(tài)網(wǎng)6個(gè)GNSS測(cè)站ZTD產(chǎn)品和4個(gè)探空資料來(lái)驗(yàn)證MERRA-2和ERA5再分析資料積分計(jì)算ZTD和ZWD的精度，并對(duì)其誤差的時(shí)空變化特性進(jìn)行分析，可為廣西區(qū)域?qū)α鲗友舆t模型構(gòu)建和GNSS水汽探測(cè)研究對(duì)數(shù)據(jù)源的使用提供參考，因此具有重要的現(xiàn)實(shí)意義。

2 數(shù)據(jù)來(lái)源與處理方法

2.1 MERRA-2和ERA5再分析資料

MERRA-2是由美國(guó)NASA提供的最新大氣再分析資料(https：//goldsmr4.Gesdisc.eosdis.nasa.gov/data/MERRA2)，其平面分辨率高達(dá)0.5°×0.625°(緯度差和經(jīng)度差)、垂直分辨率有42層(層頂高度約為50 km)、分層資料的時(shí)間分辨率不低于 6 h、地表資料的時(shí)間分辨率為 1 h。分層數(shù)據(jù)包括氣壓、溫度、比濕和位勢(shì)高；地表數(shù)據(jù)有氣壓、溫度、比濕及地表高程。

ERA5是最新一代的ECMWF再分析資料(https：//cds.climate.copernicus.eu/)，其水平分辨率為0.25°×0.25°，垂直分辨率為37層(層頂高度約為 47 km)、分層及地表資料的時(shí)間分辨率為 1 h。分層數(shù)據(jù)包括氣壓、位勢(shì)、溫度和比濕；地表數(shù)據(jù)有氣壓、2 m露點(diǎn)溫度、位勢(shì)及比濕。

本文選用2017年廣西區(qū)域的MERRA-2和ERA5再分析資料，通過(guò)積分計(jì)算ZTD和ZWD，因此，需要獲取廣西區(qū)域6個(gè)陸態(tài)網(wǎng)GNSS站和4個(gè)探空站各站的最近四個(gè)格網(wǎng)點(diǎn)對(duì)應(yīng)的再分析資料，并利用陸態(tài)網(wǎng)6個(gè)GNSS站的實(shí)測(cè)ZTD和4個(gè)探空站的資料來(lái)評(píng)價(jià)MERRA-2和ERA5資料計(jì)算ZTD和ZWD的精度。GNSS測(cè)站和探空站位置如圖1所示。

圖1 GNSS測(cè)站和探空站分布圖

2.2 處理方法

由于MERRA-2和ERA5再分析資料層頂高度分別約為 50 km和 47 km，頂層幾乎沒(méi)有濕延遲影響，頂層之上的ZTD用Saastamoinen模型求解，再分析資料高度范圍內(nèi)ZTD和ZWD用積分方法計(jì)算，計(jì)算公式如下所示。

ZTD=ZTDtop+ZTDlevel

(1)

ZTDtop=0.002276*Ptop/1-0.00266*cos(2φ)-2.8*10-7htop

(2)

(3)

(4)

式中：Ptop是頂層的氣壓值；φ是緯度；htop是頂層高度；hgiven是待定點(diǎn)的高程；N是大氣折射指數(shù)；Nw是濕折射率。計(jì)算公式如下所示。

N=k1×(P-e)/T+k2×e/T+k3×e/T2

(5)

Nw=k2×e/T+k3×e/T2

(6)

e=q×P/0.622

(7)

式中：k1=77.604K/hPa、k2=64.79K/hPa、k3=377600K2/hPa；P為氣壓，e為水汽壓，q為比濕。

GNSS測(cè)站和探空站點(diǎn)一般不與再分析資料的格網(wǎng)點(diǎn)重合，且站點(diǎn)高程基準(zhǔn)也與格網(wǎng)點(diǎn)不同，因此需要通過(guò)內(nèi)插方式獲取公式(5)～(7)中所需的氣象參數(shù)，內(nèi)插之前必須進(jìn)行高程基準(zhǔn)的統(tǒng)一。再分析資料采用的高程系統(tǒng)為位勢(shì)高，GNSS測(cè)站的高程系統(tǒng)為大地高，探空站則為海拔高，海拔高與位勢(shì)高之間的差異對(duì)ZTD的高程改正影響較小，可忽略不計(jì)，但大地高和位勢(shì)高之間的差異則不可忽視，可采用EGM2008模型實(shí)現(xiàn)GNSS測(cè)站與格網(wǎng)點(diǎn)高程基準(zhǔn)的統(tǒng)一[14，15]。由于站點(diǎn)與其附近四個(gè)格網(wǎng)點(diǎn)的高程不一致，若直接以格點(diǎn)高程為起始高程開(kāi)始積分計(jì)算ZTD/ZWD，再內(nèi)插出站點(diǎn)位置處的ZTD/ZWD，不可避免地會(huì)降低插值精度，將會(huì)影響精度評(píng)估的結(jié)果。因此，本文直接積分計(jì)算出最近四個(gè)格網(wǎng)點(diǎn)在GNSS測(cè)站和探空站點(diǎn)高度的ZTD/ZWD值，這樣保證了GNSS測(cè)站和探空站點(diǎn)與最近四個(gè)格網(wǎng)點(diǎn)高度的一致，從而消除了ZTD/ZWD在高程方向上的影響。

再分析資料按氣壓分層，每一層對(duì)應(yīng)不同的海拔高程，如果待定點(diǎn)高程位于再分析資料的高度范圍以內(nèi)，對(duì)于氣溫和比濕來(lái)講，可利用相鄰層之間的氣象數(shù)據(jù)進(jìn)行線性內(nèi)插計(jì)算得到，氣壓則利用公式(8)、式(9)計(jì)算得到[14]，反之，需要在垂直方向上進(jìn)行一定的外推。對(duì)于氣溫，取平均遞減率 -6.5 K/km計(jì)算對(duì)應(yīng)高度上的估值[16]。對(duì)于氣壓等參數(shù)，則利用最底下三層的平均參數(shù)遞減率外推對(duì)應(yīng)高度上的估值[14]。

(8)

(9)

式中：hupper和hlower分別為上下兩層的位勢(shì)高，pupper和Plower分別為上下兩層氣壓，pz為高度z處的氣壓，pi為高度h處的氣壓。

根據(jù)上述方法計(jì)算出最近4個(gè)格網(wǎng)點(diǎn)在GNSS測(cè)站和探空站高度處的ZTD/ZWD后，采用反距離加權(quán)法(IDW)[17，18]來(lái)進(jìn)行水平方向的插值，最終獲得GNSS測(cè)站和探空站處的ZTD/ZWD值。

3 精度評(píng)估

本文以廣西區(qū)域陸態(tài)網(wǎng)2017年6個(gè)GNSS測(cè)站ZTD產(chǎn)品和4個(gè)探空站數(shù)據(jù)為參考值，評(píng)價(jià)MERRA-2和ERA5再分析資料積分計(jì)算廣西區(qū)域ZTD/ZWD的精度，并使用偏差(bias)與均方根誤差(RMSE)作為精度指標(biāo)，其公式為：

(10)

(11)

3.1 利用陸態(tài)網(wǎng)ZTD產(chǎn)品驗(yàn)證MERRA-2/ERA5資料計(jì)算ZTD的精度

利用2017年的MERRA-2和ERA5分層再分析資料積分計(jì)算得6個(gè)陸態(tài)網(wǎng)GNSS測(cè)站的ZTD時(shí)間序列，并將GNSS實(shí)測(cè)的ZTD結(jié)果作為參考值，計(jì)算每個(gè)GNSS測(cè)站ZTD的日均偏差和RMSE誤差，并進(jìn)而得到每個(gè)測(cè)站ZTD的年均偏差和RMSE誤差，結(jié)果如表1和圖2所示。

陸態(tài)網(wǎng)ZTD產(chǎn)品檢驗(yàn)MERRA-2/ERA5資料計(jì)算ZTD的精度統(tǒng)計(jì) 表1

由表1可知，MERRA-2再分析資料計(jì)算ZTD的偏差值的范圍為 0.31 cm～1.01 cm，平均偏差為 0.61 cm；而RMSE誤差的變化范圍在 1.58 cm～1.77 cm，平均值為 1.72 cm；ERA5再分析資料計(jì)算ZTD的偏差值的范圍為 -0.42 cm～0.29 cm，平均偏差為 -0.16 cm；而RMSE誤差的變化范圍在 1.22 cm～ 1.59 cm，平均值為 1.34 cm。ERA5的偏差和RMSE誤差最大值都比MERRA-2要小，ERA5的RMSE平均值比MERRA-2小約22.1%。由圖2可知，廣西地區(qū)MERRA-2資料計(jì)算ZTD的偏差均為正值，說(shuō)明MERRA-2資料在廣西地區(qū)計(jì)算ZTD的值偏大。而ERA5資料計(jì)算ZTD的偏差大多為負(fù)值，說(shuō)明ERA5資料在這些測(cè)站計(jì)算ZTD的值偏小。在RMSE誤差方面，每個(gè)GNSS測(cè)站ERA5資料計(jì)算ZTD的RMSE誤差均小于MERRA-2資料。此外，MERRA-2資料計(jì)算ZTD的RMSE誤差在緯度上沒(méi)有明顯的變化趨勢(shì)，而ERA5資料計(jì)算ZTD的RMSE誤差在廣西北部呈現(xiàn)相對(duì)較小的值，在廣西南部呈現(xiàn)相對(duì)較大的值，其原因主要是沿海地區(qū)水汽含量高、變化較為劇烈，導(dǎo)致計(jì)算ZTD的誤差較大。

圖2 MERRA-2和ERA5資料計(jì)算ZTD的偏差和RMSE分布情況

對(duì)廣西區(qū)域6個(gè)GNSS站MERRA-2/ERA5資料計(jì)算ZTD的偏差和RMSE分別作日均統(tǒng)計(jì)，圖3給出了分別位于廣西東、西、南、北四個(gè)方向的GXWZ、GXBS、GXBH和GXGL測(cè)站MERRA-2/ERA5計(jì)算ZTD的日均偏差和RMSE變化情況?？梢?jiàn)，在GXBH、GXBS、GXGL和GXWZ站，MERRA-2資料計(jì)算ZTD的日均偏差在夏季變化較大，而ERA5資料計(jì)算ZTD的日均偏差在全年的大部分時(shí)間內(nèi)相對(duì)穩(wěn)定且表現(xiàn)出較小的值。在RMSE誤差方面，MERRA-2/ERA5資料計(jì)算ZTD的日均RMSE誤差在夏季變化較大，主要原因是廣西區(qū)域夏季水汽變化較為劇烈，但其RMSE誤差基本保持在 2 cm以內(nèi)[10～12]。

圖3 GNSS測(cè)站MERRA-2/ERA5資料計(jì)算ZTD的日均偏差和RMSE變化

對(duì)廣西區(qū)域6個(gè)GNSS測(cè)站MERRA-2/ERA5資料計(jì)算ZTD的偏差和RMSE誤差分別作月均和季度統(tǒng)計(jì)，結(jié)果如圖4和表2所示。其中，圖4只給出了GXBH、GXBS、GXGL和GXWZ站MERRA-2/ERA5資料計(jì)算ZTD的月均偏差和RMSE誤差變化情況?？梢钥闯?，MERRA-2資料計(jì)算ZTD的月均偏差在大部分月份中表現(xiàn)為正偏差，說(shuō)明MERRA-2資料在這些測(cè)站計(jì)算ZTD值偏大，而ERA5資料計(jì)算的月均偏差在大部分月份表現(xiàn)為負(fù)偏差，說(shuō)明ERA5資料在這些測(cè)站計(jì)算ZTD值偏小，且在5～9月MERRA-2/ERA5計(jì)算ZTD的偏差絕對(duì)值明顯大于其他月份。該時(shí)段正處于廣西區(qū)域的雨季，因此可能與大氣中水汽含量增多和大氣對(duì)流強(qiáng)烈而引起再分析資料誤差增大有關(guān)。在RMSE誤差方面，MERRA-2/ERA5計(jì)算ZTD的RMSE誤差在夏季月份表現(xiàn)出相對(duì)較大的值，在冬季月份表現(xiàn)出相對(duì)較小的值，其原因如上所述。表2說(shuō)明MERRA-2/ERA5資料計(jì)算ZTD的偏差和RMSE誤差整體上呈現(xiàn)出了一定的季節(jié)變化規(guī)律。同時(shí)，MERRA-2資料計(jì)算ZTD的最大RMSE誤差僅在 2.3 cm左右，ERA5資料計(jì)算ZTD的最大RMSE誤差不超過(guò) 1.8 cm，進(jìn)一步表明MERRA-2/ERA5再分析資料計(jì)算的ZTD信息表現(xiàn)出了良好的季節(jié)性能。

圖4 GNSS測(cè)站MERRA-2/ERA5資料計(jì)算ZTD的月均偏差和RMSE變化

GNSS測(cè)站MERRA-2/ERA5資料計(jì)算ZTD的季度偏差和RMSE統(tǒng)計(jì) 表2

3.2 利用探空站數(shù)據(jù)驗(yàn)證MERRA-2/ERA5資料計(jì)算ZWD的精度

為了驗(yàn)證MERRA-2/ERA5再分析資料計(jì)算ZWD的精度，利用2017年廣西區(qū)域4個(gè)探空站時(shí)間分辨率為 12 h的剖面數(shù)據(jù)來(lái)檢驗(yàn)MERRA-2/ERA5再分析資料計(jì)算ZWD的精度。首先計(jì)算出廣西區(qū)域每個(gè)探空站在UTC 0：00和12：00時(shí)刻的ZWD數(shù)據(jù)，進(jìn)而得到每個(gè)探空站點(diǎn)處的日均偏差和RMSE誤差，最終統(tǒng)計(jì)得到廣西區(qū)域每個(gè)探空站MERRA-2/ERA5資料計(jì)算ZWD的年均偏差和RMSE誤差，結(jié)果如表3和圖5所示。

探空站數(shù)據(jù)檢驗(yàn)MERRA-2/ERA5資料計(jì)算ZWD的精度統(tǒng)計(jì) 表3

圖5 MERRA-2和ERA5資料計(jì)算ZWD的偏差和RMSE分布情況

由表3可知，MERRA-2再分析資料在廣西區(qū)域計(jì)算ZWD的偏差值變化范圍分別為 -0.57 cm～ 0.87 cm，平均偏差分別為 0.07 cm；ERA5再分析資料在廣西區(qū)域計(jì)算ZWD的偏差值變化范圍分別為 -1.73 cm～-0.01 cm，平均偏差分別為 -0.80 cm。由此看出，MERRA-2資料計(jì)算ZWD具有較小的偏差值，而ERA5計(jì)算ZWD出現(xiàn)了較大的絕對(duì)偏差值，但是其在廣西的平均偏差值仍較小。在RMSE誤差方面，MERRA-2計(jì)算ZWD的變化范圍分別為 1.72 cm～ 2.37 cm，平均值分別為 2.15 cm；ERA5計(jì)算ZWD的變化范圍分別為 1.60 cm～ 2.57 cm，平均值分別為 2.09 cm。由圖5可知，廣西地區(qū)ERA5資料計(jì)算ZWD的偏差均為負(fù)值，說(shuō)明ERA5資料計(jì)算ZWD的值偏小。MERRA-2/ERA5資料計(jì)算ZWD的RMSE誤差在廣西的北部和東部相對(duì)較大，其原因主要是廣西東部和桂林地區(qū)降水豐沛，致使出現(xiàn)相對(duì)較大的RMSE誤差。

對(duì)廣西區(qū)域4個(gè)探空站MERRA-2/ERA5資料計(jì)算ZWD的偏差和RMSE分別作日均統(tǒng)計(jì)，如圖6所示?？梢?jiàn)，MERRA-2/ERA5資料計(jì)算ZWD的日均偏差在全年絕大部分時(shí)間內(nèi)均表現(xiàn)為相對(duì)較大的值，主要原因是受該地區(qū)復(fù)雜氣候的影響。在RMSE誤差方面，MERRA-2/ERA5資料計(jì)算ZWD的RMSE誤差在夏季變化較大，其原因如前所述，盡管如此，RMSE誤差基本保持在 3 cm以內(nèi)。

圖6 探空站MERRA-2/ERA5資料計(jì)算ZWD的日均偏差和RMSE變化

探空站MERRA-2/ERA5資料計(jì)算ZWD的季度偏差和RMSE統(tǒng)計(jì) 表4

對(duì)廣西區(qū)域4個(gè)探空站MERRA-2/ERA5資料計(jì)算ZWD的偏差和RMSE誤差分別作月均和季度統(tǒng)計(jì)，結(jié)果如圖7和表4所示。可以看出，MERRA-2計(jì)算ZWD的月均偏差，在57957和59431站的大部分月份中表現(xiàn)為正偏差，說(shuō)明MERRA-2資料計(jì)算的ZWD值偏大，而在59211和59265站的大部分月份中表現(xiàn)為負(fù)偏差，說(shuō)明其計(jì)算的ZWD值較小。ERA5計(jì)算ZWD的月均偏差，在57957、59211和59265站的大部分月份中表現(xiàn)為負(fù)偏差，尤其在59211和59265站，其在夏季的時(shí)間內(nèi)表現(xiàn)出較大的負(fù)偏差，說(shuō)明這些站在夏季受水汽劇烈變化的影響較大。在RMSE誤差方面，MERRA-2/ERA5資料計(jì)算ZWD的RMSE誤差，在5-9月份相對(duì)較大，其原因如上所述。由表4可知，MERRA-2資料計(jì)算ZWD的偏差無(wú)明顯的季節(jié)變化，其RMSE誤差表現(xiàn)出一定的季節(jié)變化；ERA5資料計(jì)算ZWD的偏差和RMSE誤差具有明顯的季節(jié)變化規(guī)律。由此表明，在廣西區(qū)域MERRA-2/ERA5資料計(jì)算ZWD時(shí)受季節(jié)變化的影響相對(duì)較大。

圖7 探空站MERRA-2/ERA5資料計(jì)算ZWD的月均偏差和RMSE變化

4 結(jié) 論

本文分別以廣西地區(qū)4個(gè)探空站2017年探空數(shù)據(jù)計(jì)算的ZWD及陸態(tài)網(wǎng)6個(gè)GNSS測(cè)站2017年的ZTD為參考值，評(píng)估利用廣西地區(qū)MERRA-2和ERA5再分析資料計(jì)算ZTD/ZWD的精度，通過(guò)統(tǒng)計(jì)分析ZTD/ZWD的偏差和RMSE誤差的時(shí)間和空間分布，論證利用MERRA-2和ERA5再分析資料計(jì)算ZTD/ZWD的可行性，得到以下結(jié)論：

(1)以GNSS ZTD為參考值，MERRA-2和ERA5再分析資料計(jì)算ZTD的年平均偏差和RMSE誤差分別為 0.61 cm/1.72 cm和 -0.15 cm/1.34 cm。ZTD的偏差和RMSE誤差在廣西的東部和沿海地區(qū)出現(xiàn)相對(duì)較大的誤差。同時(shí)，ZTD的RMSE誤差呈現(xiàn)明顯的季節(jié)變化規(guī)律，總體上夏季大，冬季小。

(2)以探空數(shù)據(jù)計(jì)算的ZWD為參考值，MERRA-2和ERA5再分析資料計(jì)算ZWD的年平均偏差和RMSE誤差分別為 0.07 cm/2.15 cm和 -0.80 cm/2.09 cm。ZWD的RMSE誤差在廣西的東部和桂林地區(qū)出現(xiàn)了較大的誤差，同時(shí)，RMSE誤差也呈現(xiàn)出了類似ZTD的RMSE誤差的季節(jié)變化規(guī)律。

由此表明，MERRA-2/ERA5再分析資料在廣西區(qū)域計(jì)算ZWD/ZTD具有極高的精度和良好的穩(wěn)定性，可為廣西地區(qū)對(duì)流層延遲模型的構(gòu)建和GNSS水汽探測(cè)提供參考，也可用作精密定位中對(duì)流層的先驗(yàn)估計(jì)值。本次使用的GNSS測(cè)站和探空站數(shù)量偏少且分布不夠均勻，不能對(duì)MERRA-2和ERA5 ZTD/ZWD在空間分布的變化上做出更準(zhǔn)確的分析，接下來(lái)將利用廣西CORS站數(shù)據(jù)進(jìn)行類似研究，從而進(jìn)一步探討MERRA-2和ERA5計(jì)算ZTD/ZWD的空間分布變化。

致謝：感謝ECMWF中心提供的ERA5資料、NASA中心提供的MERRA-2資料、CMONOC資料和美國(guó)懷俄明大學(xué)提供的探空資料。