桑文剛，劉迎春，何秀鳳，王昭然

(1.山東建筑大學(xué) 測繪地理信息學(xué)院,濟(jì)南 250101；2.河海大學(xué) 地球科學(xué)與工程學(xué)院,南京 211100)

0 引言

水位信息對湖泊、大型水庫的水資源優(yōu)化調(diào)度與智能化管理具有極其重要的作用.近年來，水位監(jiān)測方法正逐步由傳統(tǒng)的人工水尺向浮子式水位計、壓力式水位計、超聲波水位計等升級換代，在自動化程度、測量精度方面有了較大地提高，可以滿足一般水位監(jiān)測要求[1].但其監(jiān)測時空分辨率較低，尤其是在庫區(qū)高邊坡等復(fù)雜地形以及惡劣天氣條件下，其精度、可靠性及可用性難以保障.同時，利用全球衛(wèi)星導(dǎo)航反射信號(GNSS-R)開展水面高度反演方面取得了一定的研究成果[2-7].一方面，GNSS 在覆蓋范圍、獲取數(shù)據(jù)量、所需成本上具有其他監(jiān)測方法無法比擬的優(yōu)勢[8].另一方面，國內(nèi)如三峽等大型工程，在庫岸高危邊坡及壩體上建立了數(shù)量眾多的GNSS為主的自動化變形監(jiān)測系統(tǒng).若充分利用現(xiàn)有的GNSS 測站，則能夠在不增加投入成本的基礎(chǔ)上，有效彌補(bǔ)現(xiàn)有手段時空分辨率低精度低的問題，極大地提高了GNSS 數(shù)據(jù)的利用率，對于大壩、湖泊等水位監(jiān)測具有重要意義[9].

MARTIN-NEIRA 等[10]在1993 年首次提出利用GNSS-R 反演海平面的高度變化，開辟了利用GNSS-R 監(jiān)測測站周圍環(huán)境這個新興的研究領(lǐng)域.2012 年，LARSON 等[11]利用瑞典Onsala 空間天文臺的GPS 信噪比(SNR)觀測數(shù)據(jù)，成功反演了海平面變化的時間序列，且與周圍的驗(yàn)潮站資料進(jìn)行對比分析，發(fā)現(xiàn)較差均方根(RMS)為4.8 cm，精度能達(dá)到10 cm.在國內(nèi)，采用GNSS-R技術(shù)進(jìn)行水位反演還處于起步階段，2012 年，吳繼忠等[12]給出一種利用單臺測量型 GPS 接收機(jī)反射信號測量水面高度的方法，利用SNR 數(shù)據(jù)對平靜湖面進(jìn)行了測高試驗(yàn)，結(jié)果表明利用SNR 反演水面高度在較好的條件下精度能達(dá)到±3 cm.2018 年，匡翠林等[13]基于數(shù)據(jù)的分析軟件，反演了國際GNSS 服務(wù)(IGS)站SC02 長達(dá)14 年的海平面變化時間序列，通過與該站附近的驗(yàn)潮站數(shù)據(jù)進(jìn)行對比，發(fā)現(xiàn)兩者具有很好的一致性，反演精度為8 cm.2019 年，SONG 等[14]使用GNSS-R技術(shù)來監(jiān)測雙王城水庫水位的變化，通過GPS SNR 測量估算的水位變化的每日時間序列可以達(dá)到厘米級的精度.2020 年，張弛等[15]利用平靜湖面上的 GNSS 靜態(tài)觀測數(shù)據(jù)進(jìn)行水面高度的反演，確定了GNSS-R技術(shù)反演水面高度的能力，驗(yàn)證了普通大地測量型接收機(jī)應(yīng)用于監(jiān)測水面的有效性.程文[16]利用SNR 數(shù)據(jù)成功反演了Calcasieu Pass GPS 站的水位變化，與附近驗(yàn)潮站監(jiān)測結(jié)果表現(xiàn)出良好的一致性，均方根誤差(RMSE)在8 cm 左右，同時證明了岸基GPS 站可對水位連續(xù)監(jiān)測.同年，南陽等[17]在西安灞河懸挑平臺上使用3 組不同的儀器進(jìn)行測高實(shí)驗(yàn)，指出目前L1 頻段的SNR 更適合測高反演.

為了充分利用我國大中型水庫、大壩為數(shù)眾多的GNSS 變形監(jiān)測系統(tǒng)提供的豐富的數(shù)據(jù)資源，發(fā)揮GNSS-R 反演水位在該領(lǐng)域的技術(shù)優(yōu)勢，本文詳細(xì)推導(dǎo)了GNSS SNR 數(shù)據(jù)反演水面高度的原理，給出相應(yīng)的數(shù)據(jù)處理方法，并選擇濟(jì)南臥虎山水庫開展實(shí)驗(yàn)，通過在壩體及泄洪閘口不同位置設(shè)站觀測，證明該方法的精度及可靠性.

1 GNSS-R 反演水面高度原理及數(shù)據(jù)處理方法

1.1 基于SNR 水面測高原理

由圖1 可知，位于水系岸邊的GNSS 設(shè)備所接收的是衛(wèi)星的直射信號與水面反射信號疊加的混合信號.其中反射信號成分所占比重會在SNR 得到直觀反映.因此，SNR 通常被用來評判信號質(zhì)量，其值越高表示接收到的信號質(zhì)量越好.直射信號和反射信號之間的干擾會在SNR 數(shù)據(jù)中產(chǎn)生特征變化，而該特征變化取決于接收機(jī)天線相位中心與反射介質(zhì)表面之間的距離[9].

圖1 GNSS-R 反演水面的原理

由圖1 可知，接收機(jī)天線同時接收衛(wèi)星直射及水面的反射信號，h是天線相位中心到瞬時水平面的垂直反射距離，θ 是瞬時衛(wèi)星的直射信號和水平面的夾角，即衛(wèi)星高度角.紅色部分為直射信號與反射信號之間的路徑差D，可表示為

其中，D1表示為

依據(jù)圖1 所示幾何關(guān)系，由D1可以推導(dǎo)出D2

由式(2)和式(3)可得，直射信號與反射信號之間的路徑差D為

根據(jù)反射信號與直射信號之間的路徑差可以得到相位差φ為

式中，φ是信號的載波波長 λ、衛(wèi)星高度角 θ、垂直反射距離h的函數(shù)，并且與高度角的正弦值 sin θ 成正比.

根據(jù)接收機(jī)載波相位跟蹤環(huán)SNR 和振幅之間的相關(guān)關(guān)系[18]

式中：AH是合成信號的振幅；Ad是直射信號振幅；Ar是反射信號振幅；cos θ 是直射信號和反射信號夾角的余弦值.

在合成信號中，直射信號占據(jù)主要成分[19]，且直射信號的振幅Ad遠(yuǎn)大于反射信號的振幅Ar，為了有效分析反射信號成分的影響，需采用二次多項(xiàng)式去除趨勢項(xiàng)Ad得到 SNR 殘差序列，用 δ(SNR) 來表示.

結(jié)合式(5)和式(6)，SNR 殘差序列可用式(7)表示：

式中，A、φ是 δ(SNR) 中反射信號的振幅和相位.令x=sin θ，f=，式(7)可表示為

利用頻譜分析獲取f，進(jìn)而可計算垂直距離h.考慮到衛(wèi)星高度角正弦值變化是不均勻的，根據(jù)式(7)可知 δ(SNR) 相對于高度角正弦值是非等間隔采樣的，因此需要采用Lomb-Scargle 周期算法代替?zhèn)鹘y(tǒng)的快速傅里葉變換(FFT)方法，以解決非等間隔采樣問題，來獲得殘差序列頻譜能量的主頻[20].

對于殘差序列 δSNR(tj)，j=1，2，3···，N，其功率譜可表達(dá)為頻率f的函數(shù)

式中：Pδ(f)是頻率f為的周期信號的功率；δSNR(tj)是殘差序列信號；tj是殘差序列信號所對應(yīng)的時間；N為數(shù)據(jù)統(tǒng)計量；τ 為時間平移不變量，用來調(diào)節(jié)非等間隔采樣的數(shù)據(jù)引起的相移[21-22].

故可以通過Lomb-Scargle 周期算法獲取f，進(jìn)而得到有效高度值[23]，即反射面到天線相位中心的垂直反射距離h

1.2 數(shù)據(jù)處理方法

基于以上SNR 反演水位的原理，在GNSS-R技術(shù)反演水面高度的數(shù)據(jù)處理過程中，精細(xì)化獲取包含有水面反射信號的SNR 數(shù)據(jù)，需要對采集的原始數(shù)據(jù)進(jìn)行預(yù)處理.在測站站心坐標(biāo)系下，計算觀測時段各衛(wèi)星的高度角和方位角.根據(jù)接收機(jī)與水域相對位置，初步判定反射點(diǎn)落入該水域衛(wèi)星的偽隨機(jī)噪聲碼(PRN)及其對應(yīng)時段.從觀測文件中篩選出所需的SNR 數(shù)據(jù)后，進(jìn)一步對其精細(xì)化提取.在衛(wèi)星信號頻率選擇上，SNR 數(shù)值較小會導(dǎo)致反演精度較低，經(jīng)對比，S1 均值大于S2 均值，因此選用L1 頻率的SNR數(shù)據(jù)[9].在衛(wèi)星高度角選擇上，在相對平靜的水庫湖泊，選取5°～20°高度角區(qū)間衛(wèi)星，SNR 數(shù)據(jù)因更易受到多路徑效應(yīng)的影響，周期性震蕩更明顯，有利于水面高度的反演.通過繪制高度角與SNR 關(guān)系曲線圖，通?，F(xiàn)形式有上升弧段、下降弧段以及包含完整的上升下降弧段.由于頻譜分析時，高度角的正弦值應(yīng)單調(diào)遞增.倘若兩個弧段都有記錄，則需將SNR 數(shù)據(jù)劃分為兩段，如若只含有上升或者下降單個弧段，則直接處理.完成可用于反演的SNR 信號精細(xì)化提取后，進(jìn)而用二次多項(xiàng)式去除SNR 趨勢項(xiàng)，獲得SNR 殘差序列，對該殘差序列進(jìn)行Lomb-Scargle 算法分析獲得主頻率，最后根據(jù)式(10)實(shí)現(xiàn)水面高度的反演.數(shù)據(jù)處理技術(shù)流程如圖2 所示.

圖2 GNSS-R 數(shù)據(jù)處理技術(shù)流程圖

2 GNSS-R 水位監(jiān)測實(shí)驗(yàn)及結(jié)果分析

2.1 數(shù)據(jù)來源

為了驗(yàn)證上述GNSS-R技術(shù)監(jiān)測水位變化理論方法的有效性以及探究GNSS-R 測高的可靠性和適用性，本文利用布設(shè)在水庫岸邊的兩個測站的觀測數(shù)據(jù)進(jìn)行實(shí)驗(yàn)分析.實(shí)驗(yàn)地點(diǎn)為濟(jì)南市臥虎山水庫，如圖3 所示，實(shí)驗(yàn)所用的接收機(jī)為Novatel OEM615 接收機(jī)，可以同時接收GPS 信號和GLONASS 信號.本文只使用GPS L1 波段的SNR 數(shù)據(jù)，接收機(jī)的采樣間隔設(shè)置為1 s，實(shí)驗(yàn)觀測時間為北京時間2019 年12 年1 日.實(shí)驗(yàn)期間同時實(shí)測一個水位數(shù)據(jù)，用于檢驗(yàn)GNSS-R 的水位反演結(jié)果.

圖3 測站位置、觀測環(huán)境及儀器布設(shè)

因本實(shí)驗(yàn)對反射信號的要求最大程度上來自水面，在測站的地理位置選擇上應(yīng)避免其他反射信號的干擾.其中測站1 位于臥虎山水庫大壩上，該位置附近的斜坡可以在一定程度上阻擋來自岸邊的反射信號干擾；測站2 選取的位置是在水庫閘門附近一個豎直的平臺上，最大程度保證了接收機(jī)記錄的反射信號全部來自水面的反射信號.圖3 中(a)、(b)、(c)是測站1 和測站2 的具體位置和周圍環(huán)境以及儀器布設(shè).

2.2 實(shí)驗(yàn)分析

為了驗(yàn)證GNSS-R技術(shù)應(yīng)用于水庫水面的可靠性，將數(shù)據(jù)預(yù)處理軟件對數(shù)據(jù)進(jìn)行分析，基于該軟件可得到觀測時段落入反射區(qū)的衛(wèi)星的情況，對測站1 落入反射區(qū)衛(wèi)星G5 的L1 頻率所有歷元時間的有效SNR 數(shù)據(jù)可視化分析.經(jīng)對比L1 頻率的噪聲小且精度高，因此采用GPS L1 波段的SNR 數(shù)據(jù).由雙縱軸圖4 可知，衛(wèi)星高度角呈單調(diào)下降趨勢，因此不需再進(jìn)行弧段劃分.經(jīng)過去除趨勢項(xiàng)，得到SNR 殘差序列，即PRN 5 號衛(wèi)星在下降階段去除趨勢項(xiàng)后的SNR 殘差序列.最后對得到的殘差序列進(jìn)行頻譜分析，通過LSP 頻譜分析得到主頻率，然后根據(jù)式(10)得到測站1 對應(yīng)的垂直反射距離h，具體如圖5 所示.

圖4 G5 衛(wèi)星SNR 和高度角時序圖

圖5 測站1 垂直反射距離h

為了驗(yàn)證在不同測站GNSS-R 測高方法的適用性，在測站2 進(jìn)行了同樣的測高實(shí)驗(yàn).實(shí)驗(yàn)數(shù)據(jù)選用測站2 落入反射區(qū)衛(wèi)星G15 的L1 頻率所有歷元時間的有效SNR 數(shù)據(jù).同時由雙縱軸圖6 可知，衛(wèi)星高度角是下降的，所以只含有單個下降弧段，因此不需要再進(jìn)行劃分.圖7 給出了PRN 15 號衛(wèi)星在下降階段去除趨勢項(xiàng)、重采樣后的SNR 殘差序列，截至高度角5°～20°，橫坐標(biāo)為隨高度角變化的重采樣時間點(diǎn)，縱坐標(biāo)為 SNR 值的線性變化值伏特.最后對殘差序列進(jìn)行頻譜分析，根據(jù)式(10)得到對應(yīng)的垂直反射距離h，如圖8 所示.

圖6 G15 衛(wèi)星L1 SNR 和高度角時序圖

圖7 去除趨勢項(xiàng)后的SNR 殘差序列圖

圖8 測站2 垂直反射距離h

2.3 結(jié)果綜合分析

由表1 可知，將通過GNSS-R技術(shù)反演得到天線相位中心到水平面的距離與實(shí)測距離對比分析可得，垂直反射距離h的反演精度可以達(dá)到厘米級，測站1 較差為3 cm，測站2 較差為5 cm.同時由實(shí)驗(yàn)結(jié)果可知，按照最大程度使得反射信號來自于水面的需求選擇的測站1 和測站2 均可獲得有效水面高度.即無論是在臥虎山水庫大壩上的測站1，還是在水庫閘門附近豎直平臺上的測站2，均能取得較好的水位反演結(jié)果.

表1 反演水面高與實(shí)測高度對比表 m

3 結(jié)束語

為了充分發(fā)揮GNSS-R 在水位監(jiān)測上的優(yōu)勢，有效利用大壩已建成的GNSS 變形監(jiān)測系統(tǒng)測站的數(shù)據(jù)，本文在詳細(xì)推導(dǎo)GNSS-R技術(shù)監(jiān)測水面高度的原理基礎(chǔ)上，以SNR 為主要觀測量，分析精細(xì)化提取用于反演衛(wèi)星及其有用時段方法，進(jìn)而進(jìn)行去除趨勢項(xiàng)及LSP 頻譜分析，并詳細(xì)介紹了數(shù)據(jù)處理方法與技術(shù)流程.通過在濟(jì)南臥虎山水庫壩體及泄洪閘不同位置進(jìn)行實(shí)測數(shù)據(jù)反演、分析，并與現(xiàn)場實(shí)測水位高度記錄數(shù)據(jù)進(jìn)行對比，結(jié)果表明通過精細(xì)化選擇衛(wèi)星及其有效觀測時段SNR，按照本文理論方法均得到有效反演高度，通過與實(shí)測水面高度比較，兩測站較差分別為3 cm、5 cm.反演精度達(dá)到厘米級，滿足一般水位監(jiān)測的需要.充分驗(yàn)證了利用現(xiàn)有大壩和庫區(qū)已建立的GNSS 變形監(jiān)測系統(tǒng)，進(jìn)行GNSS-R 水位監(jiān)測可行性.隨著我國北斗三號(BDS-3)投入使用，可用導(dǎo)航定位衛(wèi)星增加，大大提高水位監(jiān)測時空分辨率，使得水位監(jiān)測的智慧化運(yùn)營管理更加靈活、準(zhǔn)確、多樣.