郭秋英,薛學(xué)龍,余 淼,張海平,李國(guó)偉,畢京學(xué)

(1.山東建筑大學(xué)測(cè)繪地理信息學(xué)院,濟(jì)南 250101;2.山東省國(guó)土測(cè)繪院,濟(jì)南 250102)

0 引言

大氣中水汽含量雖然很少,卻是大氣中最活躍的成分之一[1]。大氣水汽的時(shí)空變化是災(zāi)害性氣候形成的原因之一。因此,快速、精確地獲取大氣水汽對(duì)于氣象預(yù)報(bào)具有重要的意義[2]。

全球衛(wèi)星導(dǎo)航系統(tǒng)(global navigation satellite system,GNSS)具有高時(shí)間分辨率、高精度和全天候的探測(cè)大氣水汽能力,彌補(bǔ)了傳統(tǒng)水汽探測(cè)技術(shù)的不足,成為探測(cè)水汽的重要技術(shù)手段。地基GNSS反演大氣可降水量(precipitable water vapor,PWV)主要有相對(duì)定位和精密單點(diǎn)定位(precise point positioning,PPP)兩種模式,都能得到精度較高的對(duì)流層延遲(zenith total delay,ZTD)[3]。相對(duì)定位需要引入遠(yuǎn)距離測(cè)站進(jìn)行聯(lián)合解算,才能有效消除大部分GNSS觀測(cè)誤差,從而提高GNSS反演PWV 的精度。而PPP不需引入遠(yuǎn)距離測(cè)站進(jìn)行解算,具有成本低、效率高及處理快等優(yōu)勢(shì)[4]。此外,整周模糊度固定模式下的PPP(PPP ambiguity resolution,PPP-AR)技術(shù)的發(fā)展使得PPP 收斂時(shí)間進(jìn)一步縮短[5],從而實(shí)現(xiàn)PWV 估值精度的提升。

目前,世界上有四個(gè)全球衛(wèi)星導(dǎo)航系統(tǒng),包括中國(guó)的北斗(BeiDou navigation satellite system,BDS)、美國(guó)的GPS、俄羅斯的GLONASS和歐盟的伽利略(Galileo),這些系統(tǒng)相互兼容、相互發(fā)展。利用多模GNSS組合可以增加觀測(cè)衛(wèi)星數(shù)量,從而提高ZTD解算精度[6]。李星星等[7]使用單系統(tǒng)GNSS和多系統(tǒng)GNSS估計(jì)ZTD,實(shí)驗(yàn)表明,與單一GNSS相比,多模GNSS實(shí)時(shí)ZTD估計(jì)可以實(shí)現(xiàn)毫米級(jí)的精度和高可靠性。LU[8]基于多GNSS(GLONASS+GPS+伽利略+北斗),利用當(dāng)前GNSS的所有可用觀測(cè)結(jié)果,得到實(shí)時(shí)的ZTD/PWV 估值,可以實(shí)現(xiàn)毫米級(jí)的精度。李宏達(dá)等[9]利用GPS/BDS/GLONASS/Galileo組合獲取PWV 估值,實(shí)驗(yàn)表明,多系統(tǒng)組合PPP估計(jì)ZTD的初始化時(shí)間縮短,且能夠獲取更高精度的ZTD和PWV。胡鵬[10]利用實(shí)時(shí)PPP 技術(shù),分析了GNSS四系統(tǒng)反演水汽的精度差異,結(jié)果顯示,不同系統(tǒng)得到的ZTD/PWV 結(jié)果存在明顯的差異,多系統(tǒng)組合觀測(cè)獲取的水汽序列結(jié)果最為穩(wěn)健,但沒有對(duì)北斗三號(hào)獲取的PWV 進(jìn)行評(píng)估。夏鵬飛等[11]基于GNSS觀測(cè)數(shù)據(jù)與ERA-5資料,構(gòu)建對(duì)流層模型;新的ZTD 模型優(yōu)于后處理GNSS-ZTD 估值,研究結(jié)果表明,對(duì)浮點(diǎn)PPP求解施加ZTD 約束后,垂直方向的收斂時(shí)間也有顯著的改善。

以上研究結(jié)果表明,與單一或雙系統(tǒng)相比,多系統(tǒng)組合顯著提高了PWV 解算精度。由于北斗三號(hào)的B1C/B2a相較于其他頻點(diǎn)具有更好的信號(hào)強(qiáng)度[12],因此本文在利用北斗探測(cè)水汽時(shí)選擇B1C/B2a頻點(diǎn)。此外PPP-AR 在傳統(tǒng)PPP 的基礎(chǔ)上加入了模糊度固定的方法,可以提高定位精度和縮短收斂時(shí)間。因此本文基于最終精密星歷、采用不同方案下獲取ZTD 并進(jìn)一步轉(zhuǎn)換為PWV,從PPPZTD 與IGS-ZTD 的相關(guān)性、PPP-ZTD 的收斂時(shí)間、PPP-ZTD 的精度以及PPP-PWV 的精度4個(gè)方面評(píng)價(jià)PPP-AR 探測(cè)大氣水汽的性能,對(duì)于GNSS氣象學(xué)應(yīng)用具有重要的意義。

1 PPP反演大氣水汽的原理與方法

1.1 PPP基本觀測(cè)模型

GNSS PPP中偽距和載波相位觀測(cè)方程可表示為

式中,Pj,Lj分別表示偽距和載波相位觀測(cè)量;上標(biāo)s為觀測(cè)的衛(wèi)星號(hào);G,C,R,E 分別表示GPS,BDS,GLONASS和Galileo;下標(biāo)j為載波的頻率;ρs為衛(wèi)星s至接收機(jī)的幾何距離;c為真空中光速;dtr,dts分別為接收機(jī)和衛(wèi)星鐘差;d0為衛(wèi)星軌道誤差;tSTD為對(duì)流層斜路徑延遲;為電離層延遲;分別為接收機(jī)端和衛(wèi)星端的偽距硬件延遲;分別為接收機(jī)端和衛(wèi)星端載波相位硬件延遲;為偽距觀測(cè)值殘差;為載波相位觀測(cè)噪聲和其他誤差。

BDS-3中B1I/B3I偽距偏差可以通過無電離層組合消除,而BDS-3中B1C/B2a無法消除,為了進(jìn)行改正,可以使用差分碼偏差(differential code bias,DCB)產(chǎn)品進(jìn)行改正,公式為[13]

在PPP解算中,由于未檢校的相位延遲(uncalibrated phase delay,UPD)被模糊度吸收,導(dǎo)致模糊度不具備整數(shù)特性,無法直接固定為整數(shù)。因此,精確地將相位偏差進(jìn)行分離并改正是實(shí)現(xiàn)PPP模糊度固定的關(guān)鍵,目前通常采用整數(shù)鐘差法或UPD 法來消除相位偏差帶來的影響。

1.2 模糊度固定方法

1.2.1 寬巷模糊度固定

模糊度固定一般采用寬窄巷(wide lane-narrow lane,WL-NL)的一個(gè)固定的過程[14],melbournewubbena(MW)組合觀測(cè)值可以轉(zhuǎn)換為

為了避免偽距噪聲的影響,在計(jì)算LsMW時(shí)需要對(duì)偽距載波觀測(cè)值進(jìn)行平滑處理。具體如下

對(duì)平滑后的寬巷模糊度進(jìn)行星間單差計(jì)算,首先,根據(jù)寬巷模糊度固定殘差和高度角的關(guān)系,選擇高度角最大的衛(wèi)星作為基準(zhǔn)星,通過計(jì)算得到星間單差寬巷模度。然后,對(duì)星間單差模糊度取整,在進(jìn)行取整前后,寬巷模糊度的絕對(duì)值差值不大于0.25,則進(jìn)一步進(jìn)行下列驗(yàn)證

其中

式中,P表示模糊度固定成功率;為星間單差模糊度的實(shí)數(shù)解;為單差模糊度整數(shù)解;σ為求得的模糊度中誤差。

1.2.2 窄巷模糊度固定

經(jīng)過寬巷模糊度固定后,開始進(jìn)行窄巷模糊度固定。根據(jù)估計(jì)的模糊度實(shí)數(shù)解和星間單差寬巷模糊度可以得出星間窄巷模糊度實(shí)數(shù)解,公式如下

其對(duì)應(yīng)的協(xié)方差陣為

在LAMBDA 算法中,精密單點(diǎn)定位的窄巷模糊度搜索時(shí),其輸入量分別為浮點(diǎn)解參數(shù)個(gè)數(shù)、固定解參數(shù)個(gè)數(shù)、具有整數(shù)特性的實(shí)數(shù)解、對(duì)應(yīng)的協(xié)方差矩陣以及固定解和固定解的殘差平方和。整數(shù)最小二乘被認(rèn)為是模糊度取整方法中最嚴(yán)謹(jǐn)?shù)姆椒?獲取的解為最優(yōu)解。一般將最優(yōu)解和次優(yōu)解進(jìn)行Ratio檢驗(yàn)

1.3 PPP水汽反演的原理

斜路徑對(duì)流層延遲(slant total delay,STD)可以由映射函數(shù)、水平梯度、干延遲(zenith hydrostatic delay,ZHD)與濕延遲(zenith wet delay,ZWD)構(gòu)成

式中,Md為ZHD 映射函數(shù);MW為ZWD 映射函數(shù);ΔLg(e,a)為大氣在水平方向梯度變化引起的誤差;e為衛(wèi)星的截止高度角;a為衛(wèi)星的方位角。

一般采用Saastamoinen模型來計(jì)算ZHD

式中,p為測(cè)站下氣壓,單位為hPa;H為測(cè)站高度,單位為m;φ為以弧度為單位的緯度。

GNSS衛(wèi)星信號(hào)STD 映射到觀測(cè)站的天頂方向上得到ZTD,ZTD 減去ZHD 可以得到ZWD,然后利用ZWD 與PWV 的轉(zhuǎn)換關(guān)系獲得PWV。1994年,Bevis等[15]給出了ZWD 到PWV 的轉(zhuǎn)換公式,計(jì)算方法為

式中,Π為水汽轉(zhuǎn)換系數(shù);ρ為液態(tài)水密度;k1=77.604 K/hPa,k2=64.790 K/hPa,k3=3.776×105K2/h Pa分別為大氣折射率試驗(yàn)常數(shù);k'2=k2-16.52;Rv=0.462 J·g-1·K-1為水汽氣體常數(shù);Rd=0.287 J·g-1·K-1為干空氣氣體常數(shù);Tm為大氣加權(quán)平均溫度。其中Tm采用Bevis經(jīng)驗(yàn)公式?；谏鲜龇椒ň涂梢垣@取測(cè)站上方的PWV。

2 基于PPP-AR獲取ZTD,PWV的性能分析

2.1 數(shù)據(jù)獲取與實(shí)驗(yàn)方案

為評(píng)估多模PPP-AR在反演水汽方面的能力,選取全球范圍16個(gè)MGEX 觀測(cè)站,2022年3月11日到3月17日(年積日70—76)、5月7日到5月13日(年積日126—132)、7月8日到7月14日(年積日189—195)、10月23日到10月29日(年積日294—300)4個(gè)時(shí)間段的觀測(cè)數(shù)據(jù)。采用GPS(G)、BDS-3(C)、GPS+BDS-3(GC)、GPS+GLO+GAL+BDS-3(GREC)以及多系統(tǒng)浮點(diǎn)解(float-GREC)五種方案對(duì)ZTD與PWV進(jìn)行估值,其中GPS采用L1/L2、BDS-3 采用B1C/B2a、GLONASS 采用L1/L2、Galileo采用E1/E5a雙頻頻點(diǎn),前四種方案為PPP-AR,最后方案為多系統(tǒng)浮點(diǎn)解。實(shí)驗(yàn)選取全球范圍16個(gè)MGEX 觀測(cè)站分布圖見圖1。

圖1 全球范圍16個(gè)MGEX觀測(cè)站分布圖Fig.1 Distribution map of 16 MGEX observation stations around the world

實(shí)驗(yàn)采用武漢大學(xué)IGS數(shù)據(jù)中心提供的最終精密星歷產(chǎn)品。該產(chǎn)品具有較高的精度和時(shí)效性,可用于水汽的探測(cè)。本文使用Net_Diff 1.16軟件進(jìn)行解算,該軟件能夠支持所有GNSS系統(tǒng)和所有頻率信號(hào)、GNSS觀測(cè)數(shù)據(jù)和產(chǎn)品下載、雙頻和三頻PPP與PPP-AR 的解算等;且在數(shù)據(jù)處理過程中,該軟件可以更加直觀地選擇不同的頻點(diǎn)及其組合。具體解算策略如表1所示。

表1 數(shù)據(jù)解算策略Tab.1 Data solution strategy

自1998年以來,國(guó)際GNSS服務(wù)(international GNSS service,IGS)采用地面觀測(cè)站數(shù)據(jù)并結(jié)合最終精密衛(wèi)星軌道、衛(wèi)星截止高度角為7°,基于PPP技術(shù)獲取每日時(shí)間分辨率為5 min的ZTD 產(chǎn)品,精度約為4 mm[17]。這些對(duì)流層產(chǎn)品以每日文件的形式提供,覆蓋了IGS網(wǎng)絡(luò)中的350多個(gè)全球?qū)Ш叫l(wèi)星系統(tǒng)(GNSS)站點(diǎn)。因此,一般將IGS-ZTD 作為參考值。

2.2 多模GNSS探測(cè)水汽的性能分析

為分析多模PPP-AR 在探測(cè)水汽方面的性能,對(duì)不同時(shí)間段和不同衛(wèi)星組合獲取ZTD 與PWV的估值進(jìn)行了比較。從PPP-ZTD 與IGS-ZTD 的相關(guān)性、PPP-ZTD 收斂時(shí)間、PPP-ZTD 估值精度及PPP-PWV 估值精度4個(gè)方面評(píng)價(jià)多模PPP-AR 探測(cè)水汽的性能。

2.2.1 GNSS PPP-ZTD與IGS-ZTD的相關(guān)性分析

以IGS 發(fā)布的對(duì)流層產(chǎn)品(IGS-ZTD)作為參考值,評(píng)估不同方案下得到ZTD 估值與IGS-ZTD的相關(guān)性。圖2展示了2022年3月11日到3月17日(年積日70—76)ALIC 測(cè)站5 種方案下獲取ZTD 估值與IGS-ZTD 值的相關(guān)性。

圖2 ALIC測(cè)站不同星座組合的相關(guān)性Fig.2 Correlation of different constellation combinations of ALIC stations

由圖2可以看出,單系統(tǒng)獲取ZTD與IGS-ZTD的相關(guān)系數(shù)(0.902)與單系統(tǒng)C雙頻B1C/B2a的相關(guān)系數(shù)(0.895),兩者相關(guān)系數(shù)相當(dāng)。雙系統(tǒng)GC相關(guān)系數(shù)(0.936)高于單系統(tǒng)G 和單系統(tǒng)C雙頻B1C/B2a的相關(guān)系數(shù)。這表明使用兩個(gè)衛(wèi)星系統(tǒng)進(jìn)行組合可以提高相關(guān)性。此外,多系統(tǒng)的相關(guān)系數(shù)(0.938)比雙系統(tǒng)的相關(guān)系數(shù)更高,意味著使用更多衛(wèi)星系統(tǒng)進(jìn)行組合時(shí),相關(guān)性會(huì)進(jìn)一步增強(qiáng)。相較于float-GREC組合,GREC也展示了更好的相關(guān)性。

2.2.2 GNSS PPP-ZTD 收斂時(shí)間分析

以IGS 發(fā)布的對(duì)流層產(chǎn)品(IGS-ZTD)作為參考值,評(píng)估ZTD 序列的收斂時(shí)間。當(dāng)GNSS PPPZTD 與IGS-ZTD 差值小于20 mm,并且連續(xù)20個(gè)歷元差值保持在20 mm 內(nèi),認(rèn)為收斂完成,這個(gè)過程需要的時(shí)間稱為收斂時(shí)間。圖3為DGAR 測(cè)站4個(gè)不同時(shí)間段、5種方案下的收斂時(shí)間。

圖3 DGAR 測(cè)站不同時(shí)間段、不同組合下的收斂時(shí)間Fig.3 Convergence time of DGAR station under different time periods and combinations

圖3顯示了DGAR 測(cè)站4個(gè)不同時(shí)間段的收斂過程。在PPP收斂的過程中,剛開始各方案的初始ZTD 與IGS-ZTD 之間存在較大的差值,隨著時(shí)間的推移,GNSS PPP-ZTD 與IGS-ZTD 之間的差值逐漸減小,并且變化趨勢(shì)逐漸趨于一致?？梢钥闯?雙系統(tǒng)的收斂時(shí)間小于單系統(tǒng)的收斂時(shí)間,而多系統(tǒng)的收斂時(shí)間更加迅速,固定解在浮點(diǎn)解的基礎(chǔ)之上實(shí)現(xiàn)了更加迅速的收斂。多系統(tǒng)PPP 具有更好的收斂性能,可顯著提高ZTD 估計(jì)的準(zhǔn)確性和穩(wěn)定性。通過融合多個(gè)系統(tǒng)的信息,組合PPP可以更快地獲取ZTD 的估值。為了更加直觀體現(xiàn)測(cè)站的收斂時(shí)間,統(tǒng)計(jì)了所有測(cè)站不同時(shí)間段下的收斂時(shí)間,如圖4所示。

圖4 不同時(shí)間段各測(cè)站平均收斂時(shí)間Fig.4 Average convergence time of each station in different time periods

由圖4可以看出,在4個(gè)不同時(shí)間段的測(cè)站中,基于最終精密星歷獲取的單系統(tǒng)G 和C雙頻B1C/B2a的平均收斂時(shí)間分別為32,33,34,33和31,32,34,32 min,雙系統(tǒng)GC 的平均收斂時(shí)間分別為29,30,31和30 min,多系統(tǒng)固定解GREC 平均收斂時(shí)間分別為23,23,26 和24 min,多系統(tǒng)浮點(diǎn)解GREC平均收斂時(shí)間分別為26,27,28,27 min。根據(jù)研究數(shù)據(jù)表明:與G,C,GC 相比,多系統(tǒng)固定解GREC的收斂時(shí)間分別縮短了27%,25%和20%。此外,多系統(tǒng)固定解相比于多系統(tǒng)浮點(diǎn)解的收斂時(shí)間縮短11%。

這些收斂時(shí)間表明多系統(tǒng)固定解GREC 在不同時(shí)間段的收斂速度較快,相較于單系統(tǒng)和雙系統(tǒng)組合,具有更高的效率和性能。由于多系統(tǒng)在單、雙系統(tǒng)的基礎(chǔ)之上獲得了更多的衛(wèi)星數(shù)量,觀測(cè)幾何圖形也得到了改進(jìn)。這有助于PPP的收斂,進(jìn)而改善ZTD 的收斂時(shí)間和PWV 估值精度。而多系統(tǒng)浮點(diǎn)解也表現(xiàn)出一定的優(yōu)勢(shì),雖不如多系統(tǒng)固定解GREC的收斂速度快,但相較于單系統(tǒng)G 和C固定解仍有改進(jìn)。因此,多系統(tǒng)PPP-AR 技術(shù)在大氣水汽探測(cè)中表現(xiàn)出顯著的優(yōu)勢(shì)。

2.2.3 GNSS PPP反演的ZTD 精度分析

為分析不同時(shí)間段、不同星座組合下的ZTD 估值精度,選取了ALIC,BIK0兩個(gè)測(cè)站,分別計(jì)算了不同方案下獲取的ZTD 估值與對(duì)應(yīng)參考值IGSZTD 之間的偏差。以平均偏差和均方根誤差作為精度評(píng)價(jià)指標(biāo),ALIC,BIK0測(cè)站不同時(shí)間段、不同星座組合ZTD 精度對(duì)比見表2。

表2 ALIC、BIK0測(cè)站不同時(shí)間段、不同星座組合ZTD精度對(duì)比Tab.2 Comparison of ZTD accuracy of ALIC and BIK0 stations in different time periods and constellation combinations mm

表2 統(tǒng)計(jì)了ALIC 與BIK0 測(cè)站在不同時(shí)間段、不同星座組合下的ZTD 精度,GREC 組合獲得的ZTD 估值性能比單、雙系統(tǒng)有所提高,能夠獲得更高精度的ZTD 估值?？梢?與單系統(tǒng)和雙系統(tǒng)解算相比,多系統(tǒng)PPP解算ZTD 估值精度更高。

2.2.4 GNSS PPP反演PWV 的性能分析

以探空站PWV(radiosondes precipitable water vapor,RS-PWV)為參考值,對(duì)上述5種方案反演的PWV 精度進(jìn)行分析。由于探空站每天只有兩個(gè)時(shí)段(UTC時(shí)間00.00和12.00)的數(shù)據(jù),因此在分析時(shí)只保留與探空站時(shí)段一致的PWV。另外,HOB2站與探空站(編號(hào)94975)相距25.1 km,WUH2與探空站(編號(hào)57494)相距5.2 km,且對(duì)流層延遲具有較強(qiáng)的時(shí)空差異性,故分析時(shí)只采用這兩個(gè)站的PWV。5種方案獲得HOB2 測(cè)站的PWV 與探空站94975獲取的PWV 對(duì)比見圖5(a);5種方案獲得WUH2 測(cè)站的PWV 與探空站57494 獲取的PWV 對(duì)比見圖5(b)。

圖5 HOB2和WUH2站與探空站(編號(hào)94975,57494)對(duì)比Fig.5 Comparison of HOB2 and WUH2 stations with sounding stations(No.94975,57494)

由圖5可知,G,C,GC,GREC,float-GREC這5種方案得到PWV 與RS-PWV 之間的差值大部分在10 mm 之內(nèi),GREC組合獲取的PWV 整體的變化更趨向于RS-PWV。以探空站RS-PWV 為參考值,5種方案計(jì)算的WUH2和HOB2站的PWV 與其相應(yīng)的探空站PWV 的均方根誤差見表3。

表3 WUH2與HOB2站PWV的精度對(duì)比Tab.3 Comparison of PWV accuracy between WUH2 and HOB2 stations

由表3可知,相較于單系統(tǒng)(G)固定解、單系統(tǒng)(C)固定解、雙系統(tǒng)固定解及多系統(tǒng)浮點(diǎn)解,多系統(tǒng)固定解獲取PWV 估值展示了更好的性能。WUH2和HOB2測(cè)站雙系統(tǒng)GC解算均方根誤差較單系統(tǒng)G,C雙頻B1C/B2a減少0.21 mm,0.29 mm 和0.10 mm,0.05 mm;相較于GC 組合,GREC 組合的PWV 精度進(jìn)一步提升,WUH2 和HOB2 測(cè)站均方根誤差減少0.02 mm 和0.10 mm;相較于float-GREC 組合,GREC組合精度也有所提升,WUH2 測(cè)站均方根誤差減少0.02 mm,HOB2測(cè)站均方根誤差減少0.08 mm。

3 結(jié)論

本文通過研究多模PPP-AR 技術(shù)反演大氣水汽的性能,得到如下結(jié)論。

1)在相關(guān)性方面:多系統(tǒng)組合ZTD 估值與IGS-ZTD 值相關(guān)性最高,時(shí)間序列更加穩(wěn)健、平滑。此外,GREC組合估計(jì)ZTD 值更加連續(xù)、可靠。

2)在收斂時(shí)間方面:在不同的時(shí)間段基于最終精密星歷獲取ZTD,多系統(tǒng)展示了更快的收斂速度,相較于單、雙系統(tǒng)收斂時(shí)間分別縮短了27%,25%和20%。多系統(tǒng)固定解相較于多系統(tǒng)浮點(diǎn)解,收斂時(shí)間縮短了11%,表現(xiàn)出良好的性能。

3)在ZTD 精度方面:以平均偏差和均方根誤差作為精度評(píng)價(jià)指標(biāo),發(fā)現(xiàn)多系統(tǒng)狀態(tài)下精度最高,由于GREC組合可觀測(cè)到更多的衛(wèi)星和更好的空間結(jié)構(gòu),由此可改善ZTD 的精度。

4)在PWV 精度方面:基于北斗三號(hào)雙頻B1C/B2a反演的PWV 與GPS的精度相當(dāng),基于多系統(tǒng)組合PPP-AR 反演的PWV 與探空站PWV 最為接近,效果最好。

以上結(jié)果表明,在水汽探測(cè)方面,北斗三號(hào)雙頻B1C/B2a獲取PWV 估值與GPS效果相當(dāng),在多系統(tǒng)的狀態(tài)下性能更佳,可以獲得更高精度的PWV 估值,為高精度的天氣預(yù)報(bào)提供有利的參考數(shù)據(jù)。