基于VMD 和WinLSP 相結合的GNSS-R海平面高度估測模型

胡媛，袁鑫泰，劉衛(wèi)，胡慶松，江志豪，鐘李程

( 1. 上海海洋大學 工程學院，上海 201306；2. 上海海事大學 商船學院，上海 201306)

1 引言

海洋對人類社會生存和發(fā)展具有重要意義。全球氣候變化造成的海平面上升嚴重威脅到人們的生命和財產安全。因此，觀測海平面高度的變化有利于監(jiān)測自然界的變化，對于海洋學的研究和人類的發(fā)展具有重要意義[1]。隨著全球導航衛(wèi)星系統(tǒng)（Global Navigation Satellite Systems，GNSS）的不斷發(fā)展，GNSS反射（GNSS-Reflectometry，GNSS-R）技術作為一種新興的衛(wèi)星遙感技術被廣泛應用于海平面測高領域。

相比于傳統(tǒng)的測量手段，基于地基GNSS-R 的海平面高度反演技術具有全天候、覆蓋面廣、時空分辨率高、信號源多和成本低等優(yōu)勢。目前基于地基GNSS-R技術的海平面高度監(jiān)測方法一般可以分為雙天線模式的偽隨機碼[2–5]、載波相位分析[6–9]方法以及單天線模式的信噪比（Signal to Noise Ratio，SNR）分析法[10–13]。雙天線模式指的是GNSS 接收機安裝有右旋圓極化（Right-Handed Circular Polarization，RHCP）天線和左旋圓極化（Left-Handed Circular Polarization，LHCP）天線。RHCP 天線和LHCP 天線分別接收衛(wèi)星直射信號和海平面的反射信號，根據兩信號中偽隨機碼或載波相位的差別進行高度反演。偽隨機碼測高方法受限于碼片寬度，反演精度最低；載波相位分析方法的精度比較好但需要專業(yè)的接收機，成本較高不適合大規(guī)模應用。SNR 分析法則僅需一根RHCP 天線同時接收直射信號和反射信號，研究兩者干涉產生的合成SNR信號進行高度反演。而且，SNR 分析法對粗糙海平面具有更好的魯棒性。

SNR 分析法的關鍵之處在于趨勢項分離和振蕩頻率提取，傳統(tǒng)模型常采用最小二乘擬合法（Least Squares Fitting，LSF）[14]和LSP（Lomb-Scargle Periodogram）頻譜分析法[15]分別進行趨勢項分離和振蕩頻率的提取。LSF 通過最小化誤差的平方和尋找數據的最佳函數匹配，一般為低階多項式。LSF 得到的趨勢項存在特定的目標函數與之匹配，這導致其與復雜的SNR 數據進行擬合時存在局限，影響反演精度。同時，由于LSP 頻譜分析法在提取振蕩頻率過程中，會先截斷輸入信號，然后再進行周期擴展處理，得到一個虛擬的無限信號。但是，SNR 數據不是周期信號，因此這種截斷過程會造成頻譜能量的泄漏，從而導致誤差[16]。

基于上述問題，本文采用瑞典翁薩拉（Onlasa）的GTGU 站2016 年DOY（Day of Year）1-80 和美國阿拉斯加州的SC02 站2021 年DOY 152-211 的GPS 觀測數據，構建一種基于變分模態(tài)分解（Variational Mode Decomposition, VMD）和凱塞窗函數改進的LSP（記為WinLSP）的GNSS-R 海平面高度估測模型，并利用翁薩拉空間天文臺和美國國家海洋和大氣管理局（National Oceanic and Atmospheric Administration，NOAA）提供的驗潮站數據進行模型精度的驗證，該模型可為區(qū)域性海平面高度測量提供參考。

2 GNSS-R 海平面高度估測原理

2.1 經典GNSS-R 海平面高度估測模型

SNR 數據分析法進行海平面高度反演只需放置一個RHCP 天線同時接收衛(wèi)星直射信號和海平面的反射信號。圖1 展示了估測模型的幾何關系圖。

如圖1 所示，天線到海平面的反射器高度為h，衛(wèi)星仰角為θ。根據反射信號相對于直射信號存在的額外的路徑和仰角 θ，由幾何關系可求得反射器高度h，進而可求得相位差為

圖1 全球導航衛(wèi)星系統(tǒng)反射（GNSS-R）海平面高度估測模型幾何關系圖Fig. 1 Global navigation satellite systems-reflectometry(GNSS-R) sea level height estimation model geometric relationship diagram

式中， λ為載波波長。RHCP 天線接收直、反射信號的能量時，對于前者會最大程度的接收。對于后者則表現為抑制，但沒能完全抑制。因此，接收機接收到的SNR 信號可以表示為[17]

2.2 變分模態(tài)分解理論

VMD 是由Dragomiretskiy 和Zosso[18]提出的自適應、完全非遞歸的信號分離方法。VMD 通過迭代搜尋變分模型最優(yōu)解來確定每個分量的頻率中心及帶寬，實現原信號的頻域剖分和本征模態(tài)函數（Intrinsic Mode Function，IMF）分量的分離。此外，VMD 有效地解決了經驗模態(tài)分解（Empirical Mode Decomposition，EMD）存在的模態(tài)混疊、端點效應以及分解層數不可控的問題[19–20]。根據VMD 的分解原理可得[21]：

式中，X(t) 表 示輸入的原信號；uk(t)表示分解得到的第k個IMF 分量；K表示分解層數；rn(t)表示殘差項，一般為噪聲。VMD 分解得到的IMF 分量為一個非平穩(wěn)的調幅調頻信號，即

式（2）和式（3）說明了GNSS-R 海平面高度估測模型中需要分離出SNR 數據的趨勢項，即干擾項，才能通過頻譜分析方法提取振蕩頻率。針對傳統(tǒng)方法信號分離不佳的問題，VMD 能夠不局限于特定的表達式，而是基于原始信號本身的特征進行趨勢項分離，為估測模型提供更高質量的SNR 數據。為了說明VMD 對SNR 數據趨勢項的分離效果，選擇了GTGU站2016 年DOY 1 的PRN02 的5°～30°衛(wèi)星仰角范圍的原始SNR 數據進行VMD 分解，如圖2 所示?？梢悦黠@地看到VMD 有效地將原始SNR 數據分解成具有不同頻率的信號。由式（6）可知，VMD 分解得到的IMF 分量頻率從高到低排列分別是IMF1-IMF4，同時還有一個殘差項。通過比較可以看到IMF4 分量與原始SNR 數據的相關性最高，即為趨勢項。GNSS天線實際收集信號時，不僅僅接收衛(wèi)星的直射信號和經海平面反射的反射信號，同時還接收周圍環(huán)境產生的噪聲信號。這些噪聲是非平穩(wěn)的高斯白噪聲，可能會導致頻譜分析時提取頻率出現誤差。VMD 在重構振蕩項的同時也減弱了噪聲對反演精度的影響。

圖2 信噪比（SNR）數據的本征模態(tài)函數（IMF）分量的組成，從上到下，依次是IMF 分量（頻率從高到低）、殘差和原始SNR 數據（紅線）Fig. 2 Composition of the intrinsic mode function (IMF) component of the signal-to-noise ratio (SNR) data, from top to bottom, followed by IMF components (from high frequency to low frequency), residual and original SNR (red line)

2.3 凱塞窗函數

不同于周期圖法和Welch 算法[22]，LSP 頻譜分析法適用于對非均勻采樣的時域信號進行頻譜分析。在反演過程中， LSP 頻譜分析法在提取非均勻采樣的SNR 數據的振蕩頻率時可能因頻譜泄露導致反演誤差的增大。窗函數能使時域信號更好地滿足頻譜分析處理的周期性要求，減少因頻譜泄露導致的泄漏，使得所需的振蕩頻率被更好地提取。窗函數的主瓣寬度和旁瓣幅度分別影響到頻譜分析過程中頻率的分辨率過濾干擾信號的能力。最理想的窗函數是主瓣寬度最窄，旁瓣最小衰減度最大。但在實際的頻譜分析過程中很難找到同時滿足這兩個要求的窗函數。

3 海平面高度反演實驗

3.1 GTGU 站海平面高度反演實驗

位于海拔高度40.420 m 的GTGU 站（57.392 954 9°N，11.913 488 6°E）由瑞典查爾姆斯理工大學翁薩拉空間天文臺運營。天線的平均垂直位置約為距平均海平面上方4 m，海平面波動約為2 m。如圖3 所示，該站點配備了兩個Leica GRX1200 接收器和Leica AR 25 天線，可提供時間分辨率為1 s 的GPS 觀測數據[25]。為了最大化衛(wèi)星軌道的數量，GTGU 站天線安裝在正南方向。GTGU 站位置避開了海浪，海平面變化主要受潮汐的影響。本實驗采用GTGU 站2016 年DOY 1-80 的GPS L1 頻段觀測數據用于海平面高度反演。為確保SNR 數據來自海平面反射，方位角范圍選取70°～210°。距離GTGU 站約1 km 的翁薩拉驗潮站每分鐘測量一次海平面高度，可用于GNSS-R 海平面高度估測模型精度的驗證。

圖3 GTGU 站的環(huán)境和位置Fig. 3 Environment and location of GTGU station

分離趨勢項后，LSP 頻譜分析法可用于提取SNR振蕩項的振蕩頻率，并利用下列準則進行有效測高的質量控制：

（1）LSP 頻譜的最大幅度高于 20 V（1 s 的采樣數據）；

（2）LSP 頻譜的最大幅度高于2 倍平均背景噪聲；

（3）LSP 頻譜至少有10 個以上的峰值；

（4）LPS 頻譜的最大幅度高于第二峰值幅度的1.5 倍；

（5）有效海平面高度閾值[2, 6]。

3.1.1 基于VMD+LSP 的估測模型精度分析

實驗數據選取GTGU 站2016 年 DOY 1-80 的GPS L1 頻段信號的SNR 數據。為了驗證VMD 在處理SNR數據中的可行性，選取了一段SNR 數據進行VMD與LSF 趨勢項擬合效果對比分析，其結果如圖4 所示?？梢钥闯龌贚SF 得到的擬合結果對于復雜多變的SNR 數據分離效果不佳；而VMD 得到的SNR 趨勢項是基于SNR 數據本身進行的擬合，其分離效果較佳而且保留了SNR 數據的局部特征，為接下來LSP 頻譜分析提供了質量較高的時域信號。

圖4 趨勢項擬合效果對比Fig. 4 Comparison of the trend term fitting effect

實驗期間選取5°～15°、5°～25°和5°～30°不同仰角范圍的SNR 數據用于探究低、高仰角對反演結果的影響，實驗結果記錄在表1 中。從表1 中可以看到，基于VMD+LSP 的GNSS-R 海平面高度反演結果的均方根誤差（Root Mean Squared Error，RMSE）約為5 cm，相關系數達到0.98。從結果上看，低仰角范圍的反演精度更高，但受限于SNR 數據導致有效的反演點數較少。隨著衛(wèi)星仰角增大，多徑信號與直射信號的相干成分越少，干涉效應越不明顯，反演結果的RMSE 越來越大，但由于SNR 數據增多得到的有效反演高度也越來越多。因此，5°～30°仰角范圍的SNR數據被用于接下來的實驗。

表1 基于變分模態(tài)分解（VMD）算法的不同仰角范圍的海平面高度反演結果Table 1 Inversion results of sea level height in different elevation angle ranges based on variational mode decomposition (VMD)

在驗證基于VMD+LSP 的GNSS-R 海平面高度估測模型精度同時，本文使用LSF+LSP 和EMD+LSP 等模型的反演結果進行精度對比分析，如圖5 所示。從圖5 可見，不同趨勢項擬合方法的反演結果都與驗潮站數據高度一致，具體的反演結果記錄在表2 中。從表2 中，可以看到基于VMD+LSP 的新模型反演精度最好，其RMSE、相關系數和有效的反演點數分別為5.50 cm、0.98 和5 268，在RMSE 和反演點數上相比于LSF+LSP 的6.69 cm 和4 909 分別提高了約17.8%和7.3%。EMD+LSP 的反演精度次之，其RMSE、相關系數和有效的反演點數分別為5.58 cm、0.97 和4 641，相比于LSF+LSP 在精度也有16.6%的提升，但反演點數則有所減少。總體上，VMD+LSP 的反演精度優(yōu)于EMD+LSP 和LSF+LSP，這說明VMD 能夠更好地擬合SNR 數據的趨勢項，分離信號后得到了更高質量的SNR 振蕩項用于頻譜分析提取振蕩頻率。

圖5 GTGU 站基于不同趨勢項擬合方法的反演結果與驗潮站數據的對比情況（a）以及海平面高度誤差情況（b）Fig. 5 Comparison of inversion results of GTGU Station based on different trend term fitting methods with tide gauge (a)and sea level height error (b)

表2 基于最小二乘擬合+頻譜分析（LSF+LSP）、經驗模態(tài)分解+頻譜分析（EMD+LSP）和變分模態(tài)分解+頻譜分析（VMD+LSP）的GTGU 站海平面高度反演結果的精度對比Table 2 Accuracy comparison of sea level height inversion results of GTGU Station based on least squares fitting+lomb-scargle periodogram (LSF+LSP), empirical mode decomposition+lombscargle periodogram (EMD+LSP) and variational mode decomposition+lomb-scargle periodogram (VMD+LSP)

3.1.2 基于VMD+WinLSP 的估測模型精度分析

本文在構建基于VMD 的GNSS-R 海平面高度估測模型的基礎上，利用WinLSP 頻譜分析法結合提高海平面高度的反演精度。圖6 展示了LSF+WinLSP、EMD+WinLSP 和VMD+WinLSP 的海平面高度反演結果。從圖6 中可以看出各模型的反演結果與驗潮站數據高度一致，結果記錄在表3 中。對比表2 和表3，可以看出不同趨勢項擬合方法基于WinLSP 的反演結果的精度相比于LSP 有所提升。從RMSE 上看，加窗前、后基于LSF 的海平面高度反演結果精度提升得最明顯，由原本的6.69 cm 減小到5.50 cm，減少了18%。VMD 的提升效果次之，RMSE 減小了15%，由5.50 cm 減少到4.70 cm；而EMD 的RMSE 由5.58 cm減小到5.51 cm，幾乎無變化。從得到的有效反演點數上看，WinLSP 相比于LSP 提高了振蕩頻率提取的精度，使得更多的結果符合質量控制的要求。加窗前、后基于LSF、EMD 和VMD 的反演點數分別增加了11%、17%和7%。在各個組合中，基于 VMD+WinLSP的GNSS-R 海平面高度估測模型的精度最高。相比于傳統(tǒng)模型的LSF+LSP，VMD+WinLSP 的反演結果的RMSE 減小了1.99 cm，并且增加了738 點反演點數，即精度和GNSS 數據利用率分別提高了29.7%和15%。

圖6 GTGU 站基于加窗的頻譜分析（WinLSP）反演結果與驗潮站數據的對比情況Fig. 6 Comparison of inversion results based on lomb-scargle periodogram with window (WinLSP) of GTGU Station with tide gauge

表3 基于最小二乘擬合+加窗的頻譜分析（LSF+WinLSP）、經驗模態(tài)分解+加窗的頻譜分析（EMD+WinLSP）和變分模態(tài)分解+加窗的頻譜分析（VMD+WinLSP）的GTGU 站海平面高度反演結果的精度對比Table 3 Accuracy comparison of sea level height inversion results of GTGU Station based on least squares fitting+lomb-scargle periodogram with window (LSF+WinLSP), empirical mode decomposition+lomb-scargle periodogram with window(EMD+WinLSP) and variational mode decomposition+lombscargle periodogram with window (VMD+WinLSP)

3.2 SC02 站海平面高度反演實驗

SC02（48.5°N，123.0°W）站位于美國西海岸的華盛頓州星期五港口，天線的平均垂直位置約在距海平面以上5.5 m，海平面波動約3 m。該測站裝置包括一個與Trimble NETRS GPS 接收器相連的Trimble TRM29659.00 天線，以15 s 為一個采樣周期記錄數據。本文使用了2021 年DOY 152-211 的GPS L1 頻段數據開展海平面高度反演實驗。為保障SNR 數據來自海平面反射，選取該站的仰角范圍為5°～15°，方位角范圍為50°～240°。Friday Harbor 驗潮站由NOAA運營，位于SC02 站以西約300 m 處并以6 min 的采樣間隔記錄海平面觀測數據。有效高度的質量控制準則與GTGU 站點類似，其中的有效海平面高度閾值改為[2, 9]。

3.2.1 基于VMD+LSP 的估測模型精度分析

與GTGU 站點一致，本文采用LSF、EMD 和VMD 3 種不同的趨勢項擬合方法結合LSP 頻譜分析法在SC02 站開展海平面高度反演實驗，如圖7 所示。從圖7 中可以觀察到SC02 站的海平面波動相比于GTGU 站更劇烈但更加符合潮水漲退的規(guī)律。基于LSF+LSP、EMD+LSP 和VMD+LSP 的海平面高度反演結果都取得了不錯的反演精度，而且似乎更多的反演結果都集中在低于實驗期間平均海平面處。此外，一些離散的反演結果仍然存在，但都在可以接受的范圍內。不同方法的反演精度記錄在表4 中。不同方法的反演結果的精度沒有顯著差別，相關系數都達到了0.99，不過可以看到VMD 相比于LSF 和EMD 在RMSE 和反演點數上具有一定的優(yōu)勢。LSF+LSP、EMD+LSP 和VMD+LSP 三 者 的RMSE 分 別 是16.36 cm、15.83 cm 和14.46 cm，反演點數分別是1 632、1 641和1 723。相比于傳統(tǒng)的LSF+LSP 估測模型，本文提出的VMD+LSP 在RMSE 上減小了1.9 cm 且增加了91的反演點數，即RMSE 和GNSS 數據的利用率分別提高了11.6%和5.6%。

圖7 SC02 站基于不同趨勢項擬合方法的反演結果與驗潮站數據的對比情況Fig. 7 Comparison of the inversion results of SC02 Station based on different trend term fitting methods with the tide gauge

表4 基于最小二乘擬合+頻譜分析（LSF+LSP）、經驗模態(tài)分解+頻譜分析（EMD+LSP）和變分模態(tài)分解+頻譜分析（VMD+LSP）的SC02 站海平面高度反演結果的精度對比Table 4 Accuracy comparison of sea level height inversion results of SC02 Station based on least squares fitting+lomb-scargle periodogram (LSF+LSP), empirical mode decomposition+lombscargle periodogram (EMD+LSP) and variational mode decomposition+lomb-scargle periodogram (VMD+LSP)

3.2.2 基于VMD+WinLSP 的估測模型精度分析

在GTGU 站的海平面高度反演實驗的結果表明，WinLSP 頻譜分析法提取SNR 振蕩項的頻率能夠減弱因頻譜泄露帶來的影響，即提高反演精度。為驗證提出的WinLSP 方法在不同觀測站的可行性，本文在SC02 站開展重復實驗。圖8 展示了實驗期間SC02站基于LSF+WinLSP、EMD+WinLSP 和VMD+WinLSP的海平面高度反演結果。與圖7 相似，從圖8 中可以看到實驗期間平均海平面上、下區(qū)域的反演結果并不均勻，大多集中在下半區(qū)域。表5 記錄了基于WinLSP的估測模型的精度分析結果。對比表5 中不同方法的反演精度，LSF+WinLSP 和EMD+WinLSP 的反演精度十分接近，VMD+WinLSP 的反演精度則較優(yōu)于前兩者。對比表4 和表5，可以看出基于WinLSP 的反演結果的精度相比于LSP 有所提升，而且適用于不同的趨勢項擬合方法。從RMSE 上看，加窗前、后LSF的海平面高度反演結果精度提升得最明顯，由原本的16.36 cm 減小到15.45 cm，減小了6%，而EMD 和VMD 的提升效果幾乎可以忽略不計。從得到的有效反演點數上看，不同方法加窗后反演點數都得到了不同程度的提升。與傳統(tǒng)的LSF+LSP 對比發(fā)現，最優(yōu)的VMD+WinLSP 估測模型的反演結果的RMSE 小了2.02 cm，并且增加了153 點反演點數，即精度和GNSS 數據利用率分別提高了約12.3%和9.4%。

表5 基于最小二乘擬合+加窗的頻譜分析（LSF+WinLSP）、經驗模態(tài)分解+加窗的頻譜分析（EMD+WinLSP）和變分模態(tài)分解+加窗的頻譜分析（VMD+WinLSP）的SC02 站海平面高度反演結果的精度對比Table 5 Accuracy comparison of sea level height inversion results of SC02 Station based on least squares fitting+lomb-scargle periodogram with window (LSF+WinLSP), empirical mode decomposition+lomb-scargle periodogram with window(EMD+WinLSP) and variational mode decomposition+lombscargle periodogram with window (VMD+WinLSP)

圖8 SC02 站基于加窗的頻譜分析的反演結果與驗潮站數據的對比情況Fig. 8 Comparison of inversion results based on Lomb-Scargle Periodogram with window (WinLSP) of SC02 Station with tide gauge

3.3 討論與分析

通過在GTGU 站和SC02 站的海平面高度反演實驗的結果可見，不同觀測站的反演精度以及本文提出的估測模型的提升效果皆存在差異。GTGU 站所處海平面波動較平穩(wěn)，因此所得反演精度較高可達厘米級。相比而言，SC02 站所觀測海平面波動較大，且由于觀測位置處于港口處，信號質量受到港口處頻繁來往船只的影響。因此反演精度較低，僅為分米級。此外，GTGU 站和SC02 站的采樣頻率也存在明顯的差異，分別為1 s 和15 s。提出的估測模型的精度提升效果是多因素共同影響的結果。

4 結論與展望

高精度、低成本、長時間監(jiān)測海平面高度變化對于氣候學和水資源管理等研究具有重要的現實意義。根據傳統(tǒng)的海平面高度估測模型存在的不足，本文構建了VMD 和WinLSP 結合的GNSS-R 海平面高度估測模型，在提高反演精度的同時提高了GNSS 數據的利用率。本文的主要結論如下：

（1）傳統(tǒng)的趨勢項擬合方法多依賴于最小二乘低階多項式擬合法，受限于多項式的限制對SNR 數據的擬合效果有效。VMD 對于SNR 數據趨勢項的擬合不拘于具體形式，能夠根據中心頻率的不同對原始數據進行分解，為頻譜分析提供更高質量的SNR 振蕩項。

（2）傳統(tǒng)模型采用LSP 頻譜分析法提取SNR 振蕩項中的振蕩頻率時有可能出現頻譜泄露的問題，本文引入的基于凱塞窗函數改進的LSP 頻譜分析法減小了因頻譜泄露導致的誤差，提高了反演精度和GNSS數據利用率。

（3）通過在瑞典翁薩拉的GTGU 站和美國阿拉斯加州的SC02 站開展的海平面高度反演實驗的結果表明，與傳統(tǒng)模型精度相比，本文提出的VMD+WinLSP估測模型在GTGU 站的反演精度和GNSS 數據利用率分別提高了約29.7%和15.0%；在SC02 站的反演精度和GNSS 數據利用率分別提高了約12.3%和9.4%。

綜上所述，本文提出的基于VMD+WinLSP 的海平面高度估測模型是可行的，而且提高了反演精度和GNSS 數據的利用率。此外，本文提出的海平面估測模型仍存在一些需要解決的問題，如VMD 最佳分解層數的確定以及估測模型對于其他衛(wèi)星系統(tǒng)的適用性等。在將來，也許能夠結合GNSS 接收機與提出的估測模型構建一套高精度的海平面高度數據產品并推廣到全球范圍。