基于ERA5和MERRA-2再分析資料計算廣西區(qū)域ZTD和ZWD的精度分析

覃澤穎，周淼

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

1 引 言

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

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

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

2.1 MERRA-2和ERA5再分析資料

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

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

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

圖1 GNSS測站和探空站分布圖

2.2 處理方法

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

ZTD=ZTDtop+ZTDlevel

(1)

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

(2)

(3)

(4)

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

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

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

(8)

(9)

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

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

3 精度評估

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

(10)

(11)

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

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

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

由表1可知，MERRA-2再分析資料計算ZTD的偏差值的范圍為 0.31 cm～1.01 cm，平均偏差為 0.61 cm；而RMSE誤差的變化范圍在 1.58 cm～1.77 cm，平均值為 1.72 cm；ERA5再分析資料計算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資料計算ZTD的偏差均為正值，說明MERRA-2資料在廣西地區(qū)計算ZTD的值偏大。而ERA5資料計算ZTD的偏差大多為負值，說明ERA5資料在這些測站計算ZTD的值偏小。在RMSE誤差方面，每個GNSS測站ERA5資料計算ZTD的RMSE誤差均小于MERRA-2資料。此外，MERRA-2資料計算ZTD的RMSE誤差在緯度上沒有明顯的變化趨勢，而ERA5資料計算ZTD的RMSE誤差在廣西北部呈現(xiàn)相對較小的值，在廣西南部呈現(xiàn)相對較大的值，其原因主要是沿海地區(qū)水汽含量高、變化較為劇烈，導(dǎo)致計算ZTD的誤差較大。

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

4 結(jié) 論

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

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

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

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

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