基于小波分析的多模多頻GNSS-MR潮位反演

陳 殊 何秀鳳 王笑蕾 宋敏峰

1 河海大學地球科學與工程學院，南京市佛城西路8號，211100

全球氣候變暖導致海平面不斷升高，這不僅會影響氣候變化，也對人類生命安全造成嚴重威脅。因此，監測海平面高度變化意義重大[1]。近年來，GNSS多路徑反射測量(GNSS multipath reflectometry, GNSS-MR)被用于海面高度探測，目前利用GNSS-MR技術進行潮位監測在國內外已有一定的研究基礎。Larson等[2]利用大地測量型接收機獲取的GPS實測數據對潮位進行反演，實驗表明GPS-MR反演結果與驗潮站實測數據具有較好的一致性。L?fgren等[3]分別使用GPS和GLONASS系統L1和L2頻段信號進行潮位反演，證實GLONASS系統也可用于潮位反演。Roussel等[4]基于最小二乘方法，對GPS和GLONASS系統反演的潮位值進行融合。張雙成等[5]利用不同時間段的GPS信噪比數據進行海面高度反演，研究表明GPS-MR反演結果與實測潮位較差約為10 cm。Jin等[6]首次使用Beidou系統SNR(signal-to-noise ratio)和相位組合方法估計海面變化。Wang等[7]利用基于穩健估計的數據融合方法進行4個系統(GPS、GLONASS、Galileo和Beidou)聯合潮位反演，結果表明多模多頻GNSS-MR可以有效提高潮位反演精度和時間分辨率。從國內外研究來看，GNSS-MR技術可以用于監測潮位變化。

目前GNSS-MR潮位反演的原理和方法已經相對成熟，即利用LSP[8]方法提取SNR中頻率信息，再將該信息轉換為垂直反射高度RH值，最后統一到潮位基準上[9]。但是，基于LSP的潮位反演方法在一段SNR序列中通常只能得到一個時刻的潮位值，這樣可能會出現反演點缺失且難以刻畫潮位起伏變化的情況。雖然多模多頻觀測可以提供更多有效衛星弧段的SNR數據，但存在個別信號質量較差，從而影響潮位反演結果的精度。針對GNSS-MR技術反演時間分辨率不足的問題，Wang等[10]提出利用小波分析方法提取SNR序列瞬時頻率，進而得到潮位值，可大幅提高反演點的時間分辨率。但是小波分析方法過于依賴SNR數據質量，使得該方法反演潮位的精度有時比LSP方法低。

本文針對小波分析精度偏低的問題，利用多模多頻GNSS數據，分析各個頻段SNR數據質量，選取質量更好的SNR數據進行小波分析潮位反演，分別以法國BRST站和香港HKQT站15 d的數據進行實驗，并與LSP反演結果進行分析比較。

1 方法與原理

1.1 GNSS-MR反演潮位原理

當GNSS測站位于海邊時，衛星直射信號與海水表面反射信號相干的合成SNR可表示為[11]：

(1)

式中，Am和Ad分別為直射信號和反射信號振幅，ψ為直射信號與反射信號之間的相位偏差。

直射信號與反射信號相位差為：

(2)

式中，D為反射信號與直射信號的路徑差，λ為信號波長，e為衛星高度角，h為垂直反射距離。

假設反射面靜止，根據式(2)可知：

(3)

(4)

1.2 小波分析提取瞬時潮位

小波分析是一種多分辨率時頻分析方法，通過伸縮和平移形成一系列靈活窗口對信號進行多尺度細化分析[14]。在小波分析中，分析函數為一個小波ψ，小波通常被定義為[15]：

(5)

式中，m、n分別為縮放因子和平移因子，二者屬于整數集合Z={0,±1,±2,…}，t為時間。對于去除趨勢項后的SNR序列在尺度a(>0)和位置b的小波變換可表示為[10]：

W(a,b;δS(sine),ψ(sine))=

(6)

式中，δS(sine)為去除趨勢項后的SNR序列，*表示復共軛，W(a,b)為小波變換系數。MATLAB中‘centfrq’和‘scal2frq’兩個函數可以在給定小波尺度情況下計算尺度與采樣頻率之間的關系。

圖1 小波分析頻譜

對于PBAY站，3～11 m是能夠反映潮位變化的有效高度區間[12]，區間以外的能量聚集部分則認為是噪聲。在3～11 m范圍內，有兩個明顯的能量聚集區，分別位于橫軸sine=0.15、縱軸6.1 m和橫軸sine=0.41、縱軸6.3 m處。找出每個sine對應的最大振幅值(圖2(a))，再確定其對應的有效高度，即可將sine、振幅值、有效高度的三維圖轉變成sine和垂直反射距離的二維圖(圖2(b))。

圖2 小波分析頻譜二維表示

根據以上步驟即可提取SNR序列瞬時頻率，進而得到瞬時潮位反演值，但由于SNR數據質量原因，即便選取每個歷元能量最大值對應的有效高度，依然會存在反演值出現跳變的情況。為剔除這種粗差，本文根據LSP方法獲取的RH值計算標準偏差σ和均值μ，剔除在[μ-3σ,μ+3σ]區間外的瞬時垂直反射距離，同時利用小波分析得到的瞬時潮位反演值也需要采用潮汐調和函數進行海面動態改正[10]，技術流程如圖3所示。

圖3 基于小波分析反演潮位流程圖

2 多模多頻GNSS站點和數據分析

2.1 站點介紹

BRST站(48.4°N，4.5°W)位于法國西海岸Brest海港岸邊，距離BRST站292 m的Brest驗潮站可提供實測數據[16]。BRST站配備有Trimble NetR9大地測量型接收機，可提供GPS、GLONASS、Galileo、Beidou和SBAS衛星觀測數據，采樣間隔為30 s。為獲取來自海面的反射信號，方位角區間設置為130°～270°，有效高度角區間為5°～30°。

HKQT站(22.3°N，114.1°E)位于香港北邊海岸，安裝有Trimble NetR5大地測量型接收機，可提供GPS、GLONASS、Galileo、Beidou、QZSS和SBAS衛星觀測數據，采樣間隔為1 s、5 s和30 s[17]。該站點海面反射信號所對應區域的方位角為-60°～105°，有效高度角為4°～9°。距離HKQT站2 m的Quarry Bay實測潮位數據可用于對比分析。

2.2 數據分析

以BRST站為例，本文選取BRST站2020年年積日197～211的SNR數據進行實驗，其中信噪比數據GPS有4種(S1C、S2W、S2X、S5X)，GLONASS有4種(S1C、S1P、S2C、S2P)，Galileo有5種(S1X、S5X、S7X、S8X、S6X)，Beidou有5種(S1X、S5X、S2I、S7I、S6I)。其中，Galileo中S6X和Beidou中S7I數據記錄有缺失，因此在后續實驗中未采用這兩種信號。為獲取質量更高的SNR數據，本文分別從GNSS不同SNR類型數據LSP反演結果的精度及其小波譜表現來評定SNR數據的質量等級。

表1為BRST站4個衛星系統每種SNR類型數據LSP反演結果的日均反演點數量、均方根誤差(RMSE)和相關系數(CORR)統計結果，以評估采樣率及反演精度。BRST站GPS、GLONASS、Galileo和Beidou系統中反演精度最高的 SNR類型依次為S5X、S2P、S7X和S6I。不同SNR類型數據質量等級為：對于GPS，S5X>S2X>S2W>S1C；對于GLONASS，S2P>S2C>S1P>S1C；對于Galileo，S7X>S1X>S5X>S8X；對于Beidou，S6I>S5X >S2I>S1X。

表1 BRST站每種SNR類型的日均反演點數量、RMSE和CORR

圖4為BRST站2020年年積日201的GPS PRN27、GLONASS PRN16、Galileo PRN33和Beidou PRN21衛星不同類型SNR序列小波譜圖。從圖中可以看出，小波譜圖不僅可以確定每一時刻不同垂直反射距離RH的功率，而且還可以反映SNR序列的質量。質量較好的SNR序列具有明顯的功率集中區且由噪聲引起的散亂區較少，而質量較差的SNR序列功率分布較離散且噪聲多[9]。從小波譜圖可以看出，每種SNR類型數據質量等級為：對于GPS，S5X>S2X>S2W>S1C；對于GLONASS，S2P>S2C>S1C>S1P；對于Galileo，S7X>S5X>S1X>S8X；對于Beidou，S6I>S2I >S5X>S1X。

圖4 BRST站GNSS不同類型SNR序列小波譜

比較根據LSP反演結果對SNR數據質量進行評定可以發現，質量最好的SNR類型兩種方法評定結果一致，其他類型存在個別差異。小波譜的表現可直接決定瞬時反演結果的準確性，因此需要篩選質量較好的SNR類型進行反演以保證精度。

3 實驗與結果

3.1 BRST站GNSS-MR瞬時潮位結果

根據上節分析，本文最終選擇GPS系統S5X、S2X，GLONASS系統S2P，Galileo系統S7X，Beidou系統S6I數據進行GNSS-MR瞬時潮位反演，同時與LSP方法得到的潮位反演結果進行對比。BRST站LSP和小波分析反演結果如圖5所示。

圖5 BRST站LSP和小波分析方法潮位反演值

由圖5可知，BRST站利用LSP方法和小波分析方法獲得的潮位結果均與實際潮位對應良好。其中，LSP方法得到的結果較少，而小波分析得到的瞬時反演結果基本能描繪潮位的變化趨勢，滿足潮位監測需求。為對兩種方法進行評定并驗證小波分析方法獲取潮位的精度是否可靠，本文將小波分析結果以及LSP結果與驗潮站實測數據進行比較，表2為RMSE值、相關系數值、日均反演點數量和時間間隔統計。

從表2可以看出，小波分析方法和LSP方法相關系數均優于98%，小波分析方法的反演精度比LSP方法低3%～10%，但日均反演點數量卻比LSP方法提高約12倍，且平均時間間隔約為0.2 h，表明小波分析方法可大大提高反演結果的時間分辨率。

表2 BRST站小波分析和LSP方法反演結果統計

3.2 HKQT站GNSS-MR瞬時潮位結果

為進一步驗證基于小波分析的多模多頻瞬時潮位反演的有效性，本文又選取HKQT站2020年年積日211～225的SNR數據進行實驗。其中信噪比數據GPS有S1C、S2W、S2X、S5X，GLONASS有S1C、S1P、S2C、S2P，Galileo有S1X、S5X、S7X、S8X，Beidou有S2I、S7I。同樣進行LSP反演結果精度分析和小波譜圖質量分析，最終選擇GPS系統S5X，GLONASS系統S2P，Galileo系統S5X，Beidou系統S7I數據。LSP和小波分析反演結果如圖6所示。

由圖6可知，HKQT站利用LSP方法和小波分析方法獲得的潮位結果均與實測數據具有良好的對應關系。但LSP方法得到的結果較少，而小波分析的反演結果可以清晰反映潮位起伏變化。表3為小波分析結果和LSP結果RMSE值、相關系數值、日均反演點數量和時間間隔統計。

圖6 HKQT站LSP和小波分析方法潮位反演值

表3 HKQT站小波分析和LSP方法反演結果統計

從表3可以看出，HKQT站LSP方法反演結果的相關性優于小波分析方法，小波分析方法的反演精度比LSP方法低8%～17%，但日均反演點數量卻比LSP方法提高約13倍，除Beidou外小波分析方法的平均時間間隔也為0.2 h左右。

4 結 語

利用小波分析方法從SNR序列中獲取潮位值可以大大提高GNSS-MR技術反演潮位的時間分辨率，但小波分析方法的反演精度易受SNR數據質量的影響。本文針對該問題，利用多模多頻數據從LSP反演精度和小波譜圖質量兩個方面分析評定4種系統不同SNR類型數據的質量，進而選擇最優的SNR數據類型進行瞬時頻率潮位反演。從BRST站和HKQT站的反演結果來看，基于小波分析的多模多頻GNSS-MR潮位反演方法能夠在精度損失較小的情況下大幅提高反演點數量和數據采樣率。這種基于多模多頻 SNR數據的小波分析方法，可充分利用海面多路徑反射信號，實現對SNR數據的最大利用。

但是小波分析方法無法脫離LSP方法單獨使用，目前其粗差還需要根據LSP反演結果的標準偏差σ和均值μ計算的淘汰區間進行剔除。同時利用小波分析得到的瞬時潮位反演值也需要采用LSP反演的時間序列進行潮波函數擬合來確定海面動態變化情況。未來需要更加合理的粗差剔除方法，從小波譜圖的表現來對粗差進行探測和去除，進一步提高小波分析方法用于潮位反演的精度。