王笑蕾，何秀鳳，宋敏峰，陳 殊，牛紫瑾

河海大學(xué)地球科學(xué)與工程學(xué)院，江蘇 南京 211100

準(zhǔn)確的水位監(jiān)測對于水資源調(diào)控、水災(zāi)監(jiān)控及氣候氣象研究十分重要。傳統(tǒng)水位計(jì)監(jiān)測技術(shù)存在造價(jià)較高或需人工維護(hù)的缺點(diǎn)，而且需要附加技術(shù)才能實(shí)現(xiàn)基準(zhǔn)統(tǒng)一[1]。同時(shí)，衛(wèi)星測高技術(shù)無法獲得近岸水面高度且時(shí)間分辨率較低[2]。隨著全球?qū)Ш叫l(wèi)星系統(tǒng)(GNSS)的不斷發(fā)展與完善，一種GNSS干涉遙感(GNSS-interferometry reflectomety,GNSS-IR)技術(shù)被發(fā)現(xiàn)可以用來進(jìn)行水位(包括潮位)[3-7]、雪深[8-11]、土壤濕度[12-14]及地表凍融[15]等環(huán)境參數(shù)的監(jiān)測。其中，GNSS-IR水位監(jiān)測技術(shù)只基于沿岸的測量型接收機(jī)，根據(jù)其信噪比(signal-to-noise ratio,SNR)或其他觀測數(shù)據(jù)中的干涉特性，便可完成水位監(jiān)測[3-4]。利用新興的GNSS-IR技術(shù)，可以實(shí)現(xiàn)全自動、小成本、高時(shí)間分辨率、長期連續(xù)的水位監(jiān)測，并且監(jiān)測結(jié)果自動固定在穩(wěn)定框架下[5]。

文獻(xiàn)[3—4]提出GNSS-IR技術(shù)，并由此發(fā)展出一套經(jīng)典的GNSS-IR水位反演理論。目前，地基GNSS-IR技術(shù)主要有兩類誤差源：反射面高度變化和大氣折射效應(yīng)。文獻(xiàn)[4]發(fā)現(xiàn)海面起伏會引起反演結(jié)果中的高度變化誤差，并提出了相應(yīng)的改正算法。文獻(xiàn)[16]利用5個(gè)GNSS站進(jìn)行了潮位反演，并且利用經(jīng)驗(yàn)潮波系數(shù)來改進(jìn)誤差改正算法。文獻(xiàn)[17]基于潮波分析模型，提出了一種不需要經(jīng)驗(yàn)潮波系數(shù)的海面高度擬合方法，進(jìn)一步完善了高度變化誤差的改正算法。大氣折射效應(yīng)則會引起兩類誤差：一類是由于大氣折射引起的信號彎曲帶來的高度角彎曲誤差；另一類是由于大氣折射引起的對流層延遲誤差。文獻(xiàn)[18]針對大氣折射引起的高度角彎曲誤差，利用折射改正公式來改正高度角偏差。文獻(xiàn)[19]利用VMF1映射函數(shù)模型和GPT2w對流層延遲模型來改正對流層延遲誤差。

2020年6月23日，北斗三號全球?qū)Ш叫l(wèi)星系統(tǒng)星座建設(shè)完成。同時(shí)，GPS和GLONASS正在推進(jìn)現(xiàn)代化，Galileo在不斷發(fā)展，各個(gè)區(qū)域系統(tǒng)也在持續(xù)改進(jìn)和建設(shè)。多系統(tǒng)發(fā)展使得GNSS-IR技術(shù)擁有了豐富的信號源，而四系統(tǒng)的多模多頻信號已經(jīng)均被證明可以監(jiān)測水位。文獻(xiàn)[3]證明了GPS信號的水位可監(jiān)測性。文獻(xiàn)[10]證實(shí)了GLONASS系統(tǒng)的水位可監(jiān)測性。文獻(xiàn)[21]證明了北斗數(shù)據(jù)的水位可監(jiān)測性。文獻(xiàn)[22]利用GPS、GLONASS和北斗信號進(jìn)行了水位監(jiān)測試驗(yàn)。文獻(xiàn)[23]利用GPS、GLONASS、Galileo和北斗系統(tǒng)進(jìn)行了水位反演，分析了四系統(tǒng)的監(jiān)測性能，并提出了一套多模多頻融合反演算法。多模多頻信號的引入，提升了時(shí)間采樣，也為提高精度提供了可能。但是多模多頻數(shù)據(jù)也帶來了更多的誤差，需要對這些誤差進(jìn)行正確認(rèn)識及合理處理，以獲得最優(yōu)反演性能。

為了更好地挖掘GNSS-IR水位反演誤差，本文研究4個(gè)GNSS連續(xù)跟蹤站(HKQT、SW50、SW51和SW52)的多模多頻反演結(jié)果，發(fā)現(xiàn)除了上述已經(jīng)發(fā)現(xiàn)的3類誤差，還存在一種明顯的頻間偏差。而目前對GNSS-IR頻間偏差的研究極少，尚未達(dá)成將其歸為GNSS-IR誤差源的共識。為了進(jìn)一步挖掘頻間偏差的特性、加深對該誤差的認(rèn)識，本文分析此類頻間偏差的特性，并提出相關(guān)的誤差改正方法。

1 GNSS-IR原理及系統(tǒng)誤差

1.1 GNSS-IR原理

當(dāng)僅存在一次多路徑反射條件下，直射信號與反射信號之間的路徑差為D=2h(sine)[3-5]；其中，h為反射面與天線中心之間的垂直距離，稱為有效高度(reflector height,RH)；e為高度角。因此，二者之間的相位差Δφ為[3-5]

Δφ=2πD/λ=4πh(sine)/λ

(1)

式中，λ為信號波長。相位差中隱藏了一個(gè)頻率信息f

(2)

式中，ωφ為角頻率。如果忽略輸入數(shù)據(jù)對應(yīng)時(shí)間中h的變化，即在反射面在一定時(shí)間內(nèi)靜止的假設(shè)下，式(2)通過求導(dǎo)可寫為

(3)

1.2 GNSS-IR系統(tǒng)誤差

根據(jù)已有研究，本文匯總了已有共識的3類GNSS-IR系統(tǒng)誤差，并給出相關(guān)的誤差改正或削弱方法。

1.2.1 高度變化誤差

(4)

(5)

1.2.2 高度角彎曲誤差

大氣折射效應(yīng)會引起信號彎曲，從而導(dǎo)致計(jì)算高度角與實(shí)際高度角之間的偏差，而這一高度角偏差又會進(jìn)一步導(dǎo)致反演結(jié)果的偏差。文獻(xiàn)[18]發(fā)現(xiàn)該高度角彎曲誤差，并提出利用大氣折射改正公式(6)改正高度角偏差δe的方法

(6)

式中，T為溫度，單位為℃；P為氣壓，單位為mb。也可利用其他的大氣折射公式進(jìn)行高度角偏差改正[27]。利用改正后的高度角進(jìn)行反演計(jì)算，即可削弱此誤差。該誤差與高度角有關(guān)，高度角越小，誤差越大；該誤差會引起反演結(jié)果的尺度誤差及均值誤差。對于5°～20°高度角區(qū)間，反演誤差量級在亞厘米至毫米級，可視精度要求情況進(jìn)行改正。

1.2.3 對流層延遲誤差

大氣折射效應(yīng)會引起信號延遲，由于反射信號比直射信號傳播路徑更長，因而受到的延遲影響更大，從而造成反演結(jié)果的誤差。文獻(xiàn)[19,24—26]利用GPT2w模型和VMF1映射模型來修正該誤差。先計(jì)算反射信號和直射信號間的相對延遲值τT

(7)

(8)

式中，ΔhT為延遲誤差改正量。該誤差與RH大小及高度角有關(guān)，該誤差會引起反演結(jié)果的尺度誤差及均值誤差。一般情況下(有效高度10 m以內(nèi)，最低高度角選擇5°)，影響量級在厘米級至毫米級；可視精度要求情況進(jìn)行改正。對于某些RH較大的站點(diǎn)(例如安置于燈塔上部、高壩頂部、高樓頂部等的站點(diǎn))或者極低高度角的情況，該誤差必須改正。

2 頻間偏差分析

2.1 站點(diǎn)介紹

本文選取了4個(gè)GNSS連續(xù)跟蹤站。其中3個(gè)站點(diǎn)(SW50、SW51和SW52)來自山東雙王城水庫的GNSS大壩監(jiān)測系統(tǒng)，另一個(gè)站點(diǎn)(HKQT)來自香港衛(wèi)星定位參考站網(wǎng)(satellite positioning reference station network，SatRef)。

雙王城水庫位于山東濰坊壽光市，水庫大壩軸線總長為9.636 km，壩高約12.5 m，最大庫容量達(dá)到6150萬m3，是南水北調(diào)東線膠東干線工程的重要調(diào)蓄水庫。為了監(jiān)測水庫大壩的穩(wěn)定性，在大壩周圍建立了3個(gè)監(jiān)測站(SW50、SW51和SW52)和1個(gè)基準(zhǔn)站(SW43)，如圖1所示。該GNSS網(wǎng)配備CHC N72接收機(jī)及CHCC220GR天線，記錄每15 s一次的GPS(L1和L2)和北斗(B1、B2和B3)觀測值。相關(guān)信號對應(yīng)的頻率、波長及SNR類型見表1。監(jiān)測站SW50、SW51和SW52均可接收到來自水庫水面的反射信號。本文選取5°～15°高度角區(qū)間SNR弧段進(jìn)行反演。同時(shí)， 根據(jù)圖2中的反射區(qū)分布，選定對應(yīng)的水域方位角為：SW50的方位角區(qū)域?yàn)?0°～165°，SW51的方位角區(qū)域?yàn)?30°～230°，SW52的方位角區(qū)域?yàn)?70°～350°。算例選取2017-11-22—2018-01-10期間的數(shù)據(jù)進(jìn)行分析。在SW50附近有一個(gè)水位監(jiān)測站，提供每天一次的水庫水位測量值。

圖1 雙王城水庫GNSS監(jiān)測網(wǎng)Fig.1 GNSS net of Shuangwangcheng Reservoir

圖2 雙王城站點(diǎn)信號反射區(qū)Fig.2 FFZs of sites of Shuangwangcheng Reservoir

HKQT站位于香港鲗魚涌，屬于香港衛(wèi)星定位參考站網(wǎng)，站點(diǎn)環(huán)境如圖3所示。HKQT站配有Trimble NetR9接收機(jī)及Trimble 59800.0型天線。本文使用其接收到的四系統(tǒng)——GPS(L1C/A、L2P、L2C和L5)、GLONASS(G1C、G1P、G2C和G2P)、Galileo(E1、E5、E7和E8)和北斗(B1和B2)5 s采樣的觀測數(shù)據(jù)進(jìn)行反演。相關(guān)信號對應(yīng)的頻率、波長及SNR類型見表1。算例選取5°～15°高度角區(qū)間SNR弧段進(jìn)行反演。同時(shí)，根據(jù)圖4中的反射區(qū)分布，選定對應(yīng)的海域方位角為-60°～105°。算例選取2017年DOY 224—DOY 244期間的數(shù)據(jù)進(jìn)行分析。距離HKQT站點(diǎn)2m處有一隸屬Intergovern-mental Oceanographic Commission(IOC,http:∥www.ioc-sealevelmonitoring.org/)的驗(yàn)潮站Quarry Bay可提供實(shí)測的潮位數(shù)據(jù)。

圖3 HKQT站點(diǎn)環(huán)境Fig.3 HKQT site surroundings

2.2 雙王城站點(diǎn)反演序列

圖4 HKQT站點(diǎn)信號反射區(qū)Fig.4 FFZs of HKQT

表1 使用信號對應(yīng)的頻率、波長、碼及SNR類型

圖5中，L1、B1、B2及B3的反演與水位測量對應(yīng)極好；而L2的反演序列表現(xiàn)出雙序列現(xiàn)象，這是由于L2信號在LSP譜圖中表現(xiàn)出的雙波峰現(xiàn)象導(dǎo)致[18]。由表2可知，L1、B1、B2及B3的反演精度很好，RMSE在1.74～11.63 cm，且大部分集中在3～6 cm。這種高精度結(jié)果主要是由于水庫表面平靜、反射條件好、SNR多徑特性明顯導(dǎo)致。表2中，對于3個(gè)監(jiān)測站，L2信號的常數(shù)比L1信號大0.19 m，B1的常數(shù)與L1相似，B2的常數(shù)比L1信號大0.20～0.21 m，B3的常數(shù)比L1信號大0.13～0.15 m；這意味著，L2信號的RH比L1信號整體小0.19 m，B1的RH與L1相似，B2的RH比L1信號小0.20～0.21 m，B3的RH比L1信號小0.13～0.15 m。比較各個(gè)信號間的偏差值及表1中各個(gè)信號的波長值，可以發(fā)現(xiàn)偏差大小與波長呈顯著線性關(guān)系(CORR=98.82%)。因此，偏差可表示為a×δλ，其中a為常數(shù)。根據(jù)擬合，雙王城站點(diǎn)3個(gè)站反演結(jié)果頻間偏差的a=3.6。這種偏差量級在分米級，明顯不是相位中心偏移或是隨機(jī)誤差導(dǎo)致；而且其呈現(xiàn)了明確的系統(tǒng)特性(與波長相關(guān))。因此，可以確定有一類與波長相關(guān)的系統(tǒng)誤差確實(shí)存在。而波長λ等于光速與信號頻率的比值，即該偏差也與頻率呈現(xiàn)系統(tǒng)特性，故其為一種頻間偏差。這種頻間偏差可能是由電磁偏差、菲涅耳反射系數(shù)的組合相互作用及天線方向圖和兩個(gè)圓極化之間的功率偏移引起的[28]，這些硬件造成的誤差，理論上與頻率線性相關(guān)。

圖5 雙王城水庫站點(diǎn)反演序列 Fig.5 Retrievals of sites of Shuangwangcheng Reservoir

表2 雙王城站點(diǎn)RH反演結(jié)果統(tǒng)計(jì)表

2.3 香港鲗魚涌站點(diǎn)反演序列

經(jīng)過與2.2節(jié)所述相同的處理后，獲得了HKQT站各個(gè)信號的RH值。根據(jù)1.1節(jié)所述原理，可以通過站點(diǎn)坐標(biāo)及基準(zhǔn)轉(zhuǎn)換參數(shù)獲得特定常數(shù)，然后通過常數(shù)減去RH就可以得到潮位反演值。IOC(Intergovernmental Oceanographic Commission,http:∥www.ioc-sealevel monitoring.org/)提供了HKQT站的三維坐標(biāo)結(jié)果以及相應(yīng)框架轉(zhuǎn)換參數(shù)，根據(jù)相關(guān)數(shù)據(jù)計(jì)算獲得用于估計(jì)潮位的常數(shù)值為7.797 m。然而，根據(jù)2.2節(jié)的分析，頻間偏差是存在的。因此，通過RH與潮位觀測值之間的對應(yīng)關(guān)系，計(jì)算獲得對應(yīng)的常數(shù)。在計(jì)算對應(yīng)常數(shù)的過程中，發(fā)現(xiàn)不同的信號對應(yīng)的常數(shù)不同，且差距較大。本文統(tǒng)計(jì)了不同信號對應(yīng)的常數(shù)、反演個(gè)數(shù)、反演值與實(shí)測值間的RMSE和CORR。同時(shí)為了更好地表現(xiàn)信號間的頻間偏差，統(tǒng)計(jì)了各信號常數(shù)與L1信號常數(shù)之間的差距，相關(guān)結(jié)果記錄見表3。各個(gè)信號的反演序列及反演序列與實(shí)測序列的偏差如圖6所示。

表3 香港鲗魚涌站點(diǎn)RH反演結(jié)果統(tǒng)計(jì)表

圖6中，GPS、GLONASS、Galileo和北斗各個(gè)信號的反演序列與潮位實(shí)測序列對應(yīng)良好。根據(jù)表3，HKQT站各個(gè)信號的反演精度在低分米級，RMSE在15～23 cm。表3中，波長與L1近似的信號，其對應(yīng)的常數(shù)與L1信號的常數(shù)相同或近似；而與L2波長近似的信號，其對應(yīng)的常數(shù)較L1的常數(shù)低0.12 m。這意味著，與L1波長近似的信號的RH與L1信號的RH近似，L2信號及與其波長近似的信號的RH比L1信號整體小0.12 m。比較表3中各個(gè)信號間的偏差值及表1中各個(gè)信號的波長值，可以發(fā)現(xiàn)：偏差大小與波長呈顯著線性關(guān)系(CORR=96.97%)。偏差可表示為a×δλ，其中a為常數(shù)。根據(jù)擬合，HKQT站點(diǎn)對應(yīng)的a=2.3。

圖6 HKQT站點(diǎn)反演序列及偏差Fig.6 Retrievals and bias of HKQT

3 頻間偏差改正

3.1 誤差改正方法

(9)

(10)

(11)

式(10)或式(11)可簡寫為

(12)

3.2 改正結(jié)果

利用3.1節(jié)所述方法，對第2節(jié)獲得的各個(gè)信號的RH結(jié)果進(jìn)行處理，獲得各個(gè)窗口的hi輸出結(jié)果，并將其轉(zhuǎn)換至水位基準(zhǔn)，獲得對應(yīng)的水位序列及其不確定度。值得注意的是，由于雙王城站點(diǎn)的GPS L2P信號有雙波峰現(xiàn)象，因此相關(guān)RH反演值不參與改正方程組的建立。本文相關(guān)站點(diǎn)誤差改正后的多模多頻融合水位反演結(jié)果如圖7所示。

圖7 各站誤差改正后的多模多頻融合水位反演序列Fig.7 Multi-GNSS combined retrievals with error correction of sites

經(jīng)過誤差改正處理后，各個(gè)測站的水位反演結(jié)果與實(shí)測水位結(jié)果的對應(yīng)關(guān)系較好(CORR>99%)。各個(gè)測站的誤差改正后得出多模多頻融合水位反演序列的RMSE，SW50為1.56 cm，SW51為0.67 cm，SW52為0.99 cm，HKQT為6.98 cm。與第2節(jié)所述的各個(gè)信號的反演結(jié)果的RMSE相比，改正后的融合反演值的精度提高了30%～80%。同時(shí)，本文計(jì)算了未顧及頻間偏差， 只顧及其他誤差的改正結(jié)果，得出相關(guān)的多模多頻融合水位反演序列的RMSE，SW50為3.15 cm，SW51為5.90 cm，SW52為2.24 cm，HKQT為18.97 cm。與顧及頻間偏差的改正結(jié)果相比，未顧及頻間偏差的反演結(jié)果的RMSE高1.5～12 cm。具體而言，SW50和SW52站點(diǎn)改正頻間偏差后提高了約1.5 cm，SW51站提高了約5 cm，而HKQT站提高了約12 cm。這是由于SW51站的北斗反演值個(gè)數(shù)顯著大于SW50和SW52站，因此，在反演結(jié)果改正中，頻間偏差對SW51的誤差影響更大；而HKQT站，由于其多系統(tǒng)多頻信號均可用，因此在反演結(jié)果改正中，頻間偏差對融合反演結(jié)果的誤差影響遠(yuǎn)大于雙王城各站點(diǎn)。圖8顯示了HKQT站未顧及頻間偏差的改正結(jié)果，可以看出，未改正頻間偏差時(shí)，融合反演結(jié)果的精度變差、粗差增多。

圖8 HKQT站未顧及頻間偏差的多模多頻融合水位反演結(jié)果及其不確定度Fig.8 Multi-GNSS combined retrievals without consideration of inter-frequency bias and corresponding uncer-tainties at HKQT

4 結(jié)論與討論

本文分析了頻間偏差的相關(guān)特性，發(fā)現(xiàn)了頻間偏差與波長之間存在線性關(guān)系這一規(guī)律；同時(shí)，本文顧及高度變化誤差、對流層延遲誤差及頻間偏差，提出了相關(guān)的誤差改正方法。本文結(jié)果表明：顧及頻間偏差的改正結(jié)果比未顧及頻間偏差的反演結(jié)果的RMSE高1.5～12 cm；同時(shí)，改正后的融合反演值的精度較未改正的獨(dú)立信號的反演值提高了30%～80%；精度的提高得益于多模多頻信號提供的大量冗余數(shù)據(jù)及對各類誤差(包括系統(tǒng)誤差、粗差和隨機(jī)誤差)的正確處理。值得注意的是，改正方法中的相關(guān)方程組是經(jīng)典平差函數(shù)模型的表達(dá)方式，可借鑒測量數(shù)據(jù)處理方法進(jìn)行模型擴(kuò)充和解算方法擴(kuò)充；在有新的誤差被發(fā)現(xiàn)后，也可以根據(jù)誤差特性進(jìn)行擴(kuò)展。第3.1節(jié)所述誤差改正模型具有很好的擴(kuò)展性。

除了系統(tǒng)誤差，本文還表現(xiàn)出了一類L2P信號獨(dú)有的反演結(jié)果特性。雙王城站點(diǎn)中，L2的反演序列表現(xiàn)出雙序列現(xiàn)象，這是由于L2信號在LSP譜圖中表現(xiàn)出的雙波峰現(xiàn)象導(dǎo)致[18]。而HKQT站則并未表現(xiàn)出L2P雙波峰現(xiàn)象導(dǎo)致的反演結(jié)果很差的情況。對于此類表現(xiàn)在L2P上的獨(dú)特情況——有時(shí)出現(xiàn)、有時(shí)不出現(xiàn)“雙波峰”的情況，其根本原因還沒有定論，猜測與P碼調(diào)制、L1干擾或?qū)2的半無碼跟蹤有關(guān)。這類獨(dú)特現(xiàn)象還需要更深一步的研究。

此外，本文只研究了GNSS-IR水位反演的情況，后續(xù)可以對雪深、冰厚等參數(shù)進(jìn)行GNSS-IR反演，以探討在其他參數(shù)反演中是否存在頻間偏差。同時(shí)，對于形成頻間偏差的具體原因還未被發(fā)現(xiàn)；而該誤差成因?qū)τ诟玫卣J(rèn)識及改正頻間偏差非常重要。第3.1節(jié)中，利用a×δλ來表示頻間偏差大小，而表2和表3中，某些信號并不完全符合a×δλ的數(shù)值，有約1～3 cm的差距。文獻(xiàn)[28]指出，頻間偏差的產(chǎn)生可能與電磁偏差、菲涅耳反射系數(shù)的組合相互作用及功率偏移等接收機(jī)硬件有關(guān)，而這些硬件造成的誤差理論上與頻率線性相關(guān)。因此筆者認(rèn)為，對于特定站點(diǎn)，a取決于接收機(jī)硬件，它對于不同信號是唯一的，相關(guān)差距應(yīng)該是隨機(jī)誤差導(dǎo)致的。但是，該結(jié)論僅基于現(xiàn)有研究推斷而來，a與信號特征是否有關(guān)依然待研究。

目前，GNSS-IR技術(shù)仍處在快速發(fā)展階段，還有許多的誤差特性還未理清，例如地形引起的誤差[30]，不同高度角SNR質(zhì)量引起的誤差[31]等，也有許多的誤差現(xiàn)象還沒有被發(fā)現(xiàn)成因。需要后續(xù)研究者對此進(jìn)行深入研究，以實(shí)現(xiàn)更好的GNSS-IR監(jiān)測性能，推進(jìn)GNSS-IR技術(shù)的實(shí)際應(yīng)用進(jìn)程。

致謝:感謝香港衛(wèi)星定位參考站網(wǎng)提供的GNSS數(shù)據(jù)，感謝IOC(Intergovernmental Oceano-graphic Commission, http:∥www.ioc-sealevel monitoring.org/)提供的實(shí)測潮位數(shù)據(jù)。