王 杰，何秀鳳，王笑蕾，宋敏峰

小波分析在GNSS-IR潮位反演中的應(yīng)用

王 杰，何秀鳳，王笑蕾，宋敏峰

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

針對全球衛(wèi)星導(dǎo)航系統(tǒng)多徑反射（GNSS-IR）反演潮位值時易出現(xiàn)時間分辨率不足，而當(dāng)前常用的多模多頻法不適用于所有測站，加窗法會使反演結(jié)果精度降低的問題，提出1種小波分析反演方法：利用小波變換分析方法提取信噪比序列的瞬時頻率，然后反演瞬時潮位值；最后分別以美國SC02站17 d數(shù)據(jù)和香港HKQT站13 d數(shù)據(jù)加以驗證，并與經(jīng)典法、加窗法的反演結(jié)果進(jìn)行分析比較。實驗結(jié)果表明，用小波變換分析法反演潮位值的精度能達(dá)分米級，并且可由原先1段信噪比序列反演1個潮位值提升至單歷元信噪比反演1個潮位值，反演點個數(shù)大幅增加。

多路徑反射；信噪比；潮位監(jiān)測；小波分析

0 引言

隨著全球變暖，人們更加重視海洋氣候變化所帶來的問題。海面變化受氣候變化影響，同時也影響著沿海城市人類的社會活動。近年來，全球衛(wèi)星導(dǎo)航系統(tǒng)（global navigation satellite system, GNSS）與遙感相結(jié)合的技術(shù)正不斷發(fā)展完善，這種新技術(shù)稱為GNSS多徑反射（GNSS interferometric reflectometry, GNSS-IR）。通過提取連續(xù)運行參考站（continuously operating reference stations, CORS）的信噪比（signal to noise ratio, SNR）數(shù)據(jù)進(jìn)行頻譜分析，進(jìn)而反演出有效高度，具有部署簡單、分布廣泛等特點[1]。2009年文獻(xiàn)[2]提出用于積雪厚度測量，接著推廣到潮位測量中[3]，此后國內(nèi)外基于GNSS-IR技術(shù)的潮位監(jiān)測開始開展相關(guān)研究。國外方面：文獻(xiàn)[4]提出在北極格林蘭島附近海域使用全球定位系統(tǒng)（global positioning system, GPS）反射信號成功反演了潮位高度變化；文獻(xiàn)[5]利用GPS信噪比進(jìn)行潮汐反射研究，精度達(dá)到10 cm；文獻(xiàn)[6]利用多路徑信號對5個GPS站進(jìn)行潮位反演，潮位變化較大的區(qū)域反演精度為43 cm，潮位變化較小區(qū)域精度能達(dá)到6 cm左右。國內(nèi)方面：文獻(xiàn)[7]首次在我國進(jìn)行了岸基GNSS海洋遙感實驗；文獻(xiàn)[8]首次采用北斗衛(wèi)星導(dǎo)航系統(tǒng)（BeiDou navigation satellite system，BDS)在浙江大洋山海域進(jìn)行海面測高實驗，反演的海面高度精度達(dá)到了亞米級；文獻(xiàn)[9]用SC02站的數(shù)據(jù)進(jìn)行了不同時間跨度的潮位反演，與驗潮站實測值相關(guān)系數(shù)均優(yōu)于0.98；文獻(xiàn)[10]利用站SC02數(shù)據(jù)反演該站附近海域14 a的潮位變化情況，反演精度為8 cm。經(jīng)過學(xué)者大量研究驗證，GNSS-IR可以有效監(jiān)測潮位變化情況。

反演值的時間分辨率對獲取潮位變化的信息有直接影響，在保證一定精度條件下，潮位反演值分辨率越高，對于潮波系數(shù)的探測、潮位分析和海平面的準(zhǔn)確監(jiān)測越有利。然而通常一段SNR序列只能反演得到1個時刻的潮位值，這樣使得反演點數(shù)量有時不足并且不符合潮位時刻變化的特性。反演潮位常用的經(jīng)典方法是L-S譜（Lomb-Scargle periodograms, LSP）分析法，接收機采樣率越高，有效衛(wèi)星弧段的SNR值就越多，反演效果也會越好。然而這種方法處理的數(shù)據(jù)具有時間跨度，更適合反演反射面較穩(wěn)定的土壤濕度[11]、雪深[12]等，用來反演潮位便會有誤差；同時由于反演值數(shù)量受空中有效衛(wèi)星弧段的影響[13]，1個有效衛(wèi)星弧段的SNR序列只能反演1個結(jié)果，有時反演值的數(shù)量不足，不利于潮位監(jiān)測。針對反演值時間分辨率不足的問題，目前國內(nèi)外解決辦法主要歸于以下2種：第1種是從GNSS多模多頻方面著手，該方法能有效提高反演值時間分辨率，如國外使用了GPS和格洛納斯衛(wèi)星導(dǎo)航系統(tǒng)（global navigation satellite system, GLONASS）雙系統(tǒng)的不同波段反演潮位[14-15]；國內(nèi)用BDS 3個波段[16]和GPS、BDS、GLONASS 3個系統(tǒng)8個波段[17]反演潮位。但即使這樣，有時反演點仍然很少，如文獻(xiàn)[16]中26 d時間跨度反演點數(shù)量卻不到200個。同時4系統(tǒng)的CORS站并不多，遇到單系統(tǒng)的CORS站這種方法也無法使用。第2種方法稱為窗口L-S譜（windows LSP，WinLSP）分析法，即利用時間窗口或高度角窗口截取SNR序列，每個窗口得到1個反演值，由一變多。例如文獻(xiàn)[18]用的窗口長度為800個SNR（約14 min）的窗口截取SNR序列，使得每天有近80個潮位反演值。文獻(xiàn)[19]用窗口長度為5°的高度角窗口分別截取SC02和PBAY 2個站點的SNR序列，2者均能在十幾分鐘內(nèi)得到50多個反演值。然而對于第2種方法，時間窗口選取有時太短，截取的SNR序列不能滿足頻譜分析的要求，或是用高度角窗口截取后反演點增多，但反演精度會降低，而且有些測站反演潮位的有效高度角區(qū)間較小，這種方法便無法使用。

本文針對反演值時間分辨率不足的問題，利用小波分析方法提取SNR序列的瞬時頻率反演潮位值，分別以美國SC02站17d數(shù)據(jù)和香港HKQT站13 d數(shù)據(jù)加以驗證，并與LSP、WinLSP 2種方法反演結(jié)果進(jìn)行分析比較。

1 GNSS-IR反演原理及小波時頻分析

1.1 GNSS-IR反演潮位原理

圖1 GNSS-IR反演潮位原理

圖2為衛(wèi)星的SNR時序圖，在低高度角情況下SNR振蕩明顯，隨后趨于穩(wěn)定。通過對低高度角的SNR序列頻譜分析便可得到反射面參數(shù)。由于直射信號振幅遠(yuǎn)大于反射信號振幅，為減少后續(xù)頻譜分析產(chǎn)生虛假高峰，2次多項式擬合作為趨勢項的直射信號，去趨勢后得到的SNR殘差序列[20]表示為

圖2 信噪比時序

對式（2）進(jìn)行LSP分析即可得到有效高度值，結(jié)果如圖3所示。

1.2 小波分析提取瞬時潮位

1.2.1 小波分析

1.2.2 小波分析提取潮位反演值

圖4 SNR殘差序列與小波時頻分析

圖5 最大振幅值及對應(yīng)有效高度

圖5（a）中，振幅的2個峰值位于0.15和0.43附近，與圖4（b）的2個能量聚集處分別對應(yīng)。根據(jù)以上步驟即可提取SNR序列的瞬時頻率，進(jìn)一步得到反演潮位值。對于該方法反演的粗差，這里利用LSP反演結(jié)果擬合潮波曲線進(jìn)行剔除。

2 實驗與結(jié)果

2.1 香港HKQT站算例

HKQT站屬于香港衛(wèi)星定位參考站網(wǎng)，裝備TRIMBLE NETR5型接收機，位于HKQT站（如圖6所示）附近的Quarry Bay驗潮站可以提供實測數(shù)據(jù)。算例用2016年年積日第45~57天共13 d GPS L1波段數(shù)據(jù)進(jìn)行驗證。由于站點環(huán)境原因，高于9°的衛(wèi)星反射信號會被遮擋，故選取有效高度角區(qū)間4~9°，有效海面相對于接收機的方位角區(qū)間-60~105°。

圖6 HKQT站點及其周邊環(huán)境

分別用LSP、WinLSP和小波變換反演13 d的潮位值，其中WinLSP用長度為4°，步長為1°的高度角窗口截取，即[4° 8°]、[5° 9°]。反演結(jié)果分別與驗潮站實測值進(jìn)行對比，殘差時間序列和反演結(jié)果如圖7、圖8所示。

圖7 Wavelet與LSP對比結(jié)果

圖8 Wavelet與WinLSP對比結(jié)果

圖7（a）中：灰色曲線為驗潮站實測值；黑色圓點為LSP反演結(jié)果，精度為12.39 cm；淺灰色圓點為Wavelet法反演結(jié)果，精度為16.08 cm。圖8（a）中，深灰色圓點為WinLSP反演結(jié)果，精度為14.67 cm。圖7（b）、圖8（b）為3者與驗潮站值較差序列圖，可以看出3者與實測值的較差均在0處上下波動，LSP的波動范圍最小，精度更高。Wavelet與WinLSP波動幅度相當(dāng)，2者精度均低于LSP。WinLSP方法對高度角加窗后每個窗口中高度角的區(qū)間變小，即將1個長的SNR殘差序列分割為幾個短的SNR殘差序列，因此在進(jìn)行頻譜分析后其結(jié)果質(zhì)量會降低。而Wavelet方法極其依賴于SNR的質(zhì)量，受噪聲影響大，反演結(jié)果質(zhì)量也會降低。3者與實測值對比結(jié)果如表1所示。

表1 HKQT站3種方法反演結(jié)果

從表1可以看出：LSP的反演精度最高；WinLSP的反演精度高于Wavelet，而Wavelet的反演點數(shù)量最多，反演值時間分辨率有大幅提高。

圖9為2月26日當(dāng)天3種方法反演值數(shù)量對比細(xì)節(jié)圖，當(dāng)日1時、4時及23時內(nèi)，天空中因GPS衛(wèi)星無有效弧段可用，3種方法沒有反演結(jié)果。在2時、21時內(nèi)有衛(wèi)星弧段可用，但反演值為粗差均被剔除，3種方法同樣沒有有效值。而在6時、11時、12時、15時和17時內(nèi)，LSP和Wavelet無有效值可用，WinLSP有可用反演值，得益于當(dāng)SNR質(zhì)量較差時，WinLSP能在一定程度上將SNR數(shù)據(jù)中質(zhì)量差的部分與質(zhì)量好的部分進(jìn)行分割，反演值符合要求的被保留，誤差大的被剔除。但因加窗分割使其高度角區(qū)間變小，故反演精度也會降低。在剩余時間內(nèi)，由圖可知用Wavelet法反演結(jié)果數(shù)量較LSP法和WinLSP法明顯增多，3種方法反演的當(dāng)天潮位值在數(shù)量變化趨勢上近似，這是因為3者都是基于同1個SNR弧段進(jìn)行分析，并且Wavelet法依賴于SNR質(zhì)量，所以SNR質(zhì)量高，LSP反演結(jié)果有效，Wavelet同樣有有效反演值；SNR質(zhì)量低，LSP結(jié)果誤差較大而被剔除，則Wavelet也無有效結(jié)果可用，2者呈對應(yīng)關(guān)系。

圖9 第13天反演結(jié)果數(shù)據(jù)量對比細(xì)節(jié)

2.2 美國SC02站算例

SC02站位于美國華盛頓州，站點安裝Trimble NETRS接收機，并配備Trimble TRM29659.00天線，距離該站300m遠(yuǎn)有1個潮汐測量儀，可提供實測潮位數(shù)據(jù)（如圖10所示）。

圖10 SC02站點及其周邊環(huán)境

算例用2017年年積日第128~144天共17 d GPS L1波段數(shù)據(jù)進(jìn)行驗證，高度角選取5~15°，方位角選取50~240°。3種方法反演結(jié)果如圖11所示，其中WinLSP反演用長度為5°、步長為2.5°的高度角窗口截取SNR序列，即[5° 10°]、[7.5° 12.5°]、[10° 15°]。

圖11 SC02站反演結(jié)果

圖11（a）中：灰色曲線表示驗潮站實測值；淺灰色圓點表示W(wǎng)avelet反演結(jié)果，反演精度為14.47 cm；黑色圓點表示LSP反演結(jié)果，反演精度為13.09 cm。對比LSP反演結(jié)果可以看出Wavelet方法反演結(jié)果的時間分辨率明顯提高。圖11（b）中深灰色圓點表示W(wǎng)inLSP反演結(jié)果，反演精度為22.85 cm。圖11（c）為5月16日3種方法的反演值數(shù)量對比細(xì)節(jié)圖，可以看出3種方法反演結(jié)果數(shù)量變化趨勢有良好的一致性。其中，在當(dāng)日第1、3、4、6、9、11時內(nèi)LSP反演值為粗差而被剔除，第22、23、24 h內(nèi)，由于潮位較低，衛(wèi)星經(jīng)海面的反射信號強度減弱，LSP無有效反演結(jié)果。WinLSP與Wavelet結(jié)果分布情況較LSP表現(xiàn)更好，且在第1、4、11、24 h內(nèi)，LSP無結(jié)果但Wavelet仍有反演值，這期間Wavelet反演結(jié)果與驗潮站實測值的標(biāo)準(zhǔn)差為9.69 cm。3種方法與驗潮站實測值總體對比情況見表2，誤差分析及相關(guān)系數(shù)見圖12和圖13。

表2 SC02站3種方法反演結(jié)果

圖12 較差分布

圖13 反演值與驗潮站相關(guān)系數(shù)分析

結(jié)合圖12與圖13分析，用 LSP、WinLSP和Wavelet 3種方法反演的結(jié)果與實測值較差均以0為中心左右分布，說明其反演結(jié)果不存在系統(tǒng)誤差。3種方法反演值與驗潮站實測值的相關(guān)系數(shù)均優(yōu)于0.96，其中LSP方法反演結(jié)果的精度最好，Wavelet方法略低，WinLSP法最差。同時，圖中在潮位較高區(qū)域，3種方法的誤差分布明顯較為集中，這是由于隨著潮位升高，海水反射面與測站的垂直距離逐漸減小，測站接收到經(jīng)海面反射的衛(wèi)星信號增強，反演結(jié)果相比低潮位情況更理想，因此得到的有效反演值會更多。

LSP和Wavelet都是基于整段SNR序列反演潮位值，而WinLSP則是用高度角窗口將SNR序列由1段截取為3段再分別用LSP進(jìn)行反演；這樣雖然有利于潮位實時監(jiān)測，但結(jié)果精度卻降低。若對分段的SNR序列分別用小波分析處理，則其結(jié)果如圖14所示。

圖14 整體與分段Wavelet較差序列

圖14中深灰色圓點表示每個SNR整體弧段的Wavelet反演結(jié)果與實測值較差序列，淺灰色圓點表示SNR加窗分段以后的Wavelet反演結(jié)果與實測值較差序列。分段后的Wavelet反演結(jié)果精度為20.34 cm，精度較整體分析的結(jié)果下降了40.57 %，與驗潮站的相關(guān)系數(shù)由0.96下降為0.87。

經(jīng)以上案例分析，驗證了Wavelet方法可用于GNSS-IR潮位反演并且提高反演值時間分辨率，能夠獲得更多潮位變化的信息量。同時基于整段SNR序列分析的Wavelet反演精度更好，而對于分段的SNR短序列，其反演結(jié)果精度相比較差，雖然在實時性上得到了一定的提高，但精度下降嚴(yán)重。這是因為隨著弧段數(shù)據(jù)的減少，頻譜分析時的準(zhǔn)確度會降低，如果逐漸減少弧段的長度以提高潮位監(jiān)測的實時性，精度無法達(dá)到要求，不利于潮位的實時監(jiān)測。

3 結(jié)束語

本文針對GNSS-IR反演潮位時間分辨率不足的問題，從處理SNR時間序列方面著手，利用小波變換分析提取單歷元SNR的瞬時頻率進(jìn)而反演出潮位的變化情況。同時，分別用LSP、WinLSP 2種方法反演潮位結(jié)果作為對比。與驗潮站實測值相比，對于HKQT站：LSP反演結(jié)果標(biāo)準(zhǔn)差為12.39 cm，有214個反演點；WinLSP反演結(jié)果標(biāo)準(zhǔn)差為14.67 cm，有553個反演點，Wavelet反演結(jié)果標(biāo)準(zhǔn)差為16.08 cm，反演點個數(shù)為38 570；3者與實測值相關(guān)系數(shù)均優(yōu)于0.95。對于SC02站：LSP反演結(jié)果標(biāo)準(zhǔn)差為13.09 cm，精度最高，有432個反演點；Wavelet反演結(jié)果標(biāo)準(zhǔn)差為14.47 cm，精度稍低，反演點個數(shù)達(dá)29703個；WinLSP反演結(jié)果標(biāo)準(zhǔn)差為22.85 cm，有1770個反演點，精度最低；3者與實測值相關(guān)系數(shù)均優(yōu)于0.96。

對于潮位動態(tài)改正，常需要利用反演值擬合潮波函數(shù)求解反射面高度變化率，改正效果依賴于潮波函數(shù)的擬合程度；而用小波變換反演1段SNR序列得到連續(xù)潮位值，可直接進(jìn)行反射面高度變化率的求解。并且用小波變換分析法能大幅提高反演值時間分辨率，有利于獲取更多的潮位變化情況。然而該方法極易受噪聲影響，對質(zhì)量不好的SNR分析結(jié)果會存在大量粗差，有待進(jìn)一步對信號降噪處理。同時這種方法不能離開LSP法單獨使用，因為其粗差的剔除需要根據(jù)LSP結(jié)果的整體擬合來處理。因此遇到反演值時間分辨率不足的問題時，小波變換分析可以作為LSP反演潮位值的補充，進(jìn)而有利于更好地發(fā)揮GNSS-IR技術(shù)在潮位變化監(jiān)測中的重要作用。

[1]李惟, 朱云龍, 王峰, 等. GNSS多徑信號模型及測高方法[J]. 北京航空航天大學(xué)學(xué)報, 2018, 44(6): 126-132.

[2]LARSON K M, GUTMANN E D, ZA VOROTNY V U, et al. Can we measure snow depth with GPS receivers?[J]. Geophysical Research Letters, 2009, 36(17): L17502.

[3]ROUSSEL N, RAMILIEN G, FRAPPART F, et al. Sea level monitoring and sea state estimate using a single geodetic receiver[J]. Remote Sensing of Environment, 2015(171): 261-277.

[4]SEMMLING A M, BEYERLE G, STOSIUS R , et al. Detection of Arctic Ocean tides using interferometric GNSS-R signals[J]. Geophysical Research Letters , 2011, 38(4) : 155-170.

[5]LARSON K M , L?FGREN J S, HAAS R. Coastal sea level measurements using a single geodetic GPS receiver[J]. Advances in Space Research, 2013, 51(8): 1301-1310.

[6]LOFGREN J S, HAAS R, SCHERNECK H G. Sea level time series and ocean tide analysis from multipath signals at five GPS sites in different parts of the world[J]. Journal of Geodynamics, 2014(80): 66-80.

[7]王鑫, 孫強, 張訓(xùn)械, 等. 中國首次岸基GNSS-R海洋遙感實驗[J]. 科學(xué)通報, 2008, 53(5): 589-592.

[8]ZHANG Y, TIAN L, MENG W, et al. Feasibility of code-level altimetry using coastal BeiDou reflection (BeiDou-R) setups[J]. IEEE Journal of Selected Topics in Applied Earth Observations & Remote Sensing, 2015, 8(8): 4130-4140.

[9]張雙成, 南陽, 李振宇, 等. GNSS-MR技術(shù)用于潮位變化監(jiān)測分析[J]. 測繪學(xué)報, 2016, 45(9): 1042-1049.

[10]匡翠林, 劉凱, 周要宗. 基于GPS信噪比數(shù)據(jù)觀測海平面變化研究[J]. 海洋測繪, 2018, 38(6): 40-43.

[11]LARSON K M, SMALL E E, GUTMANN E, et al. Using GPS multipath to measure soil moisture fluctuations: initial results[J]. GPS Solutions, 2008, 12(3) : 173-177.

[12]LARSON K M, SMALL E E. Estimation of snow depth using L1 GPS signal-to noise ratio data[J]. IEEE Journal of Selected Topics in Applied Earth Observations & Remote Sensing, 2016, 9(10): 4802-4808.

[13]LARSON K M , RAY R D , WILLIAMS S D P . A 10-year comparison of water levels measured with a geodetic GPS receiver versus a conventional tide gauge[J]. Journal of Atmospheric and Oceanic Technology, 2017, 34(2): 295-307.

[14]L?FGREN J S, HAAS R. Sea level measurements using multi-frequency GPS and GLONASS observations[J]. EURASIP Journal on Advances in Signal Processing, 2014, 2014(1): 50.

[15]STRANDBERG J , HOBIGER T , HAAS R. Improving GNSS-R sea level determination through inverse modeling of SNR data[J]. Radio Science, 2016, 51(8): 1286-1296.

[16]李彬彬. 北斗反射信號海面測高技術(shù)的研究[D]. 上海: 上海海洋大學(xué), 2017.

[17]陳發(fā)德, 劉立龍, 黃良珂, 等. 基于多模GNSS-MR海平面測高研究[J]. 地球物理學(xué)進(jìn)展, 2018, 33(5): 1767-1772.

[18]LARSON K M, RAY R D, NIEVINSKI F G, et al. The accidental tide gauge: a GPS reflection case study from Kachemak Bay, Alaska[J]. IEEE Geoscience and Remote Sensing Letters, 2013, 10(5): 1200-1204.

[19]WANG X L, ZHANG Q, ZHANG S C. Water levels measured with SNR using wavelet decomposition and Lomb-Scargle periodogram[J]. GPS Solutions, 2018, 22(1): 22-28.

[20]LARSON K M, SMALL E E, GUTMANN E D, et al. Use of GPS receivers as a soil moisture network for water cycle studies[J]. Geophysical Research Letters, 2008, 35(24) : L24405.

[21]馬秀紅, 曹繼平, 董晟飛. 小波分析及其應(yīng)用[J]. 計算機技術(shù)與發(fā)展, 2003, 13(8): 93-94.

[22]TORRENCE C, COMPO G P. A practical guide to wavelet analysis[J]. Bull Am Meteorol Soc, 1998, 79(1): 61-78.

Application of wavelet analysis in tidal by GNSS-IR

WANG Jie,HE Xiufeng,WANG Xiaolei,SONG Minfeng

（School of Earth Sciences and Engineering, Hohai University, Nanjing 211100, China)

Aiming at the problems that it is apt to be insufficient for the time resolution, and two usual methods, the multi-system method is not applicable to all stations, and the windowing method reduces the accuracy of the inversion results in the tide value inversion of GNSS-IR, the paper proposed an inversion method using wavelet analysis: the instantaneous frequency of signal-to-noise ratio (SNR) sequence was extracted with the wavelet transform analysis method, and the instantaneous tide level values were inverted; then the 17-day data of the US SC02 station and the 13-day data of the Hong Kong HKQT station were used to implement the verification, and the inversion results of the classical method and windowing method were used to do the comparative analysis finally. Experimental result showed that: the method could have an accuracy of dm level, and a tide value could be inverted with a single epoch SNR, improved from inverted with a sequence of SNR, which means an increased number of inverted points.

multipath reflectometry; signal-to-noise ratio (SNR); tidal level monitoring; wavelet analysis

P228.4

A

2095-4999(2020)02-0082-08

王杰，何秀鳳，王笑蕾，等. 小波分析在GNSS-IR潮位反演中的應(yīng)用[J]. 導(dǎo)航定位學(xué)報, 2020, 8(2): 82-89.（WANG Jie，HE Xiufeng, WANG Xiaolei, et al. Application of wavelet analysis in tidal by GNSS-IR[J]. Journal of Navigation and Positioning, 2020, 8(2): 82-89.）

10.16547/j.cnki.10-1096.20200214.

2019-08-15

國家自然基金項目（41830110）。

王杰（1996—），男，安徽天長人，碩士研究生，研究方向為GNSS空間環(huán)境分析。

何秀鳳（1962—），女，江蘇泰州人，博士，教授，研究方向為衛(wèi)星導(dǎo)航定位。