基于GPS信噪比觀測值的水位反演研究

程文

(中冶集團(tuán)武漢勘察研究院有限公司,湖北 武漢 430080)

0 引 言

GPS-R技術(shù),是利用原本干擾定位的海平面反射的導(dǎo)航衛(wèi)星偽隨機測距碼信息或者載波信號(L1或L2),提取海平面反射面的特征,實現(xiàn)對海浪、潮汐等參數(shù)的遙感與監(jiān)測．

GPS-R技術(shù)中的基于信噪比(SNR)的海面測高方法,可以直接使用已有的普通接收機和單個天線來同時接收直射信號和反射信號,不需要對原有接收機和天線進(jìn)行改裝,實際操作更為簡便．

從1993年第一次提出利用海面反射的GPS信號測高(PARIS)概念至今,反射信號的應(yīng)用研究已經(jīng)有了幾十年的時間．使用由反射表面反射的衛(wèi)星反射信號來實現(xiàn)目標(biāo)監(jiān)測和反演的應(yīng)用,在國內(nèi)外測繪、農(nóng)業(yè)、林業(yè)、海洋學(xué)、大氣科學(xué)等諸多行業(yè)內(nèi)得到了推廣[1-12]．一些早期的概念、實驗和假設(shè)逐漸開始引起人們對GNSS反射信號的研究興趣,由于利用全球衛(wèi)星導(dǎo)航系統(tǒng)(GNSS)不僅可以全天候進(jìn)行觀測而且省時省力、成本低、精度高,所以美國,意大利,日本,西班牙,英國和澳大利亞等國家紛紛展開了對GPS-R技術(shù)的研究．

我國學(xué)者對GPS-R領(lǐng)域的研究比國外的一些高校、研究機構(gòu)要晚．但我國的研究發(fā)展速度迅猛,許多高校、研究機構(gòu)對該領(lǐng)域進(jìn)行了不同程度的研究,我國學(xué)者先后利用GPS-R技術(shù)對土壤濕度、植被生物量、積雪厚度、水位變化進(jìn)行了研究,并取得了一定的研究成果．隨著中國北斗衛(wèi)星導(dǎo)航系統(tǒng)BDS的建成和完善,中國也在這一領(lǐng)域也開始發(fā)揮主導(dǎo)作用[1-6]．

利用GPS-R技術(shù)反演水位變化研究時很多學(xué)者用到了美國華盛頓Friday港灣的SC02站,并對此站數(shù)據(jù)進(jìn)行分析,得到了很好的成果[3,5,9]．為了推廣此技術(shù)的應(yīng)用,我們此次選擇了位于美國洛杉磯Cameron港灣的Calcasieu Pass GPS站．它起初安裝時采用的是雙頻載波相位的Trimble NetRS型號GPS接收器,到2014年10月時原接收器被Trimble NETR9型號的接收器取代,新接收器可以同時接收GPS和GOLNASS信號．

1 GPS SNR反演水位變化原理

圖1是利用信噪比(SNR)法探測水位變化的示意圖,從圖中可以看出GPS接收機同時接收導(dǎo)航衛(wèi)星直接信號Ad和海水面反射信號Am．h表示天線相位中心到瞬時海水面的距離,也稱垂直反射距離;θ表示直射信號與海水面之間的角度,Ad表示直射信號振幅強度,Am表示瞬時海水面反射信號振幅強度．

由于復(fù)合信號中直射信號占主導(dǎo),為突顯經(jīng)海平面發(fā)射后的反射信號的影響,可使用低階多項式去除信噪比趨勢項Ad,余下的SNR殘差系列可表示為δ(SNR)．

圖1 SNR法探測水位變化的示意圖

SNR與信號振幅有以下關(guān)系[5]:

(1)

式中,Ac為Ad和Am組成的合成信號的振幅,可以看出Ac與Ad和Am及其夾角的余弦值有關(guān)．

圖1中Am可表示為[7-9]

(2)

式中:E為衛(wèi)星高度角;γ為載波波長;φ為反射信號振幅的相位值,若計:

t=sinE,

則式(2)可簡化為標(biāo)準(zhǔn)的余弦函數(shù),表達(dá)式為

Am=Acos(2πft+φ).

(3)

從式中可以看出f包含了h,γ等信息．由于sinE為非等間隔采樣,致使SNR殘差序列無法保證整周期截斷,因此可采用Lomb-Scargle(L-S)譜分析方法解決快速傅里葉變換無法解決的非等間隔采樣問題[11]．

由此通過f可獲取天線相位中心至瞬時潮位面的垂直反射距離h．

2 測站介紹

Calcasieu Pass GPS站位于美國洛杉磯Cameron港灣(29°46.1′N,93°20.6′W),簡稱CALC站．本文選取該站2019年1月份及2017年全年GPS數(shù)據(jù),實際潮位數(shù)據(jù)使用離此站約200 m處的驗潮站數(shù)據(jù)．

Calcasieu Pass GPS站站址及周圍環(huán)境如圖2所示．

3 實驗分析

將CALC站獲取的(.d)格式的文件轉(zhuǎn)換為(.o)格式的文件,然后提取SNR數(shù)據(jù),獲得SNR文件．使用MATLAB編寫程序進(jìn)行運算,結(jié)合該站實際情況,選取合適的方位角和高度角,解算得到振幅值、觀測衛(wèi)星數(shù)量、反演高度等信息,根據(jù)解算的信息繪制出相應(yīng)的時間序列圖像．

然后使用MATLAB讀取驗潮站的潮位信息,繪制出潮位變化時間序列圖．將潮位變化圖與反演水位變化圖進(jìn)行對比,進(jìn)行相關(guān)性分析．

3.1 高度角選取

圖3(a)、(b)分別為接收機高度為10 m和15 m時CALC站周圍的菲涅耳區(qū)域圖,可以看出當(dāng)高度角為5°、接收機高為10 m時,所接收到的反射信號中最遠(yuǎn)距離可達(dá)200 m左右;當(dāng)高度角為5°、接收機高度為15 m時,所能接收的最遠(yuǎn)信號是來自周圍250 m左右的信號．同時當(dāng)反射器高度固定時,若高度角的范圍小于10°,反射區(qū)域會擴大反射信號,從而將會收到很多來自外界的影響;若高度角的范圍大于10°,反射范圍隨高度角變大而縮小,達(dá)到40°時,反射信號已經(jīng)相對較差,無分析的意義．

(a)10 m (b)15 m

根據(jù)以上菲涅耳反射區(qū)域圖反映的情況,在此次實驗中選取的接收機GPS衛(wèi)星信號高度角為10°～40°．

3.2 方位角選取

當(dāng)我們以40°為分界點繪制不同方位角信噪比數(shù)據(jù)的頻譜圖時,多數(shù)年積日中320°～360°和0°～40°內(nèi)接收不到任何信號,所以選取方位角為40°～320°．

以任意年積日如第4年積日為例,通過對直接信號和反射信號進(jìn)行分析,得到不同方位角SNR數(shù)據(jù)的頻譜圖,如圖4所示．

(a)120°～160° (b)160°～200°

(c)200°～240° (d)240°～280°

圖4中所采用的衛(wèi)星是1～32號所有GPS衛(wèi)星,每張小圖都由上下兩部分組成,上半部分橫坐標(biāo)是高度角的范圍,縱坐標(biāo)是SNR的線性變化值(volts/volts);下半部分中橫坐標(biāo)是接收機高度,縱坐標(biāo)是SNR的線性變化值(volts/volts)．由圖4可以看出在所選的高度角內(nèi)衛(wèi)星接收到的信號有較強的連續(xù)性,不同的衛(wèi)星在同一個方位角范圍內(nèi)的頻譜圖都有較強的一致性．

3.3 SNR變化圖

以任意一顆GPS衛(wèi)星在任意兩個年積日為例,選取同一顆衛(wèi)星在不同年積日的數(shù)據(jù)進(jìn)行分析,繪制SNR值及其殘差隨高度角的變化圖．圖5為第16號衛(wèi)星在第1和第31年積日的數(shù)據(jù)．

(a)第1年積日 (b)第一年積日去趨勢后

(c)第31年積日 (d)第31年積日去趨勢后

圖5(a)是在第1年積日時SNR隨高度角變化的示意圖,圖5(b)是對圖5(a)去趨勢后的SNR殘差序列圖．圖5(c)是在第31年積日時SNR隨高度角變化的示意圖,圖5(d)是對圖5(c)去趨勢后的殘差序列圖．圖5(a)和圖5(c)中橫坐標(biāo)為所選取衛(wèi)星高度角E,縱坐標(biāo)為SNR(dB-Hz)．圖5(b)和圖5(d)中橫坐標(biāo)為高度角的正弦值sin E,縱坐標(biāo)為SNR的線性變化值(volts/volts)．

圖5中SNR的變化主要為低高度角的多路徑影響所致,進(jìn)而為利用多路徑信號反演地表參數(shù)提供了重要的數(shù)據(jù)輸入[5]．

3.4 L-S頻譜分析

以任意一顆GPS衛(wèi)星在任意年積日為例,選取同一顆衛(wèi)星在不同年積日的數(shù)據(jù),根據(jù)濾除直接信號分量后所獲取的SNR重構(gòu)序列,進(jìn)行L-S頻譜分析后反演得到海平面高度．圖6為第8號GPS衛(wèi)星在第4、第16、第22和第31年積日的數(shù)據(jù)．

(a)第4個年積日 (b)第16個年積日

(c)第22個年積日 (d)第31個年積日

圖6所示的是對第8號GPS衛(wèi)星進(jìn)行L-S頻譜分析后反演得到的不同年積日接收機天線距離瞬時海平面高度h的情況,其中橫坐標(biāo)為反射高度,表示從天線到海平面的距離;縱坐標(biāo)是SNR數(shù)據(jù)(將SNR數(shù)據(jù)從dB-Hz更改為線性單位volts/volts)．

圖6(a)所示第4年積日所獲取到的接收機天線距離瞬時海平面高度h為12.73 m,圖6(b)所示第16年積日時為12.66 m,圖6(c)所示第22年積日時為12.45 m,圖6(d)所示第31年積日時為12.44 m．

3.5 潮位對比分析

為了進(jìn)一步驗證和分析GPS SNR監(jiān)測水位變化的有效性和連續(xù)性,本文選取CALC站2019年1月份31天及2017年365天的GPS數(shù)據(jù),以及離此站約200 m處的驗潮站的實際潮位數(shù)據(jù)進(jìn)行比較．

由于通過GPS衛(wèi)星反演得到的海平面高度是從天線相位中心到海平面的高度,驗潮站所測得的結(jié)果是海平面高度對于設(shè)定基準(zhǔn)點的相對變化,因此,兩者所參考的基準(zhǔn)點不同時,對GPS反演高度與驗潮站數(shù)據(jù)進(jìn)行比較,減去各自的平均值后獲得的序列作為各自的海平面高度變化進(jìn)行比較．

圖7用的是2019年1月份的數(shù)據(jù),圖中的橫坐標(biāo)是小時(年積日×24),縱坐標(biāo)為去平均后的監(jiān)測潮位變化．其中,圓圈所表示的是根據(jù)GPS衛(wèi)星的信噪比數(shù)據(jù)反演所得的水位結(jié)果,曲線為驗潮站記錄的采樣頻率為6 min的水位變化情況．

圖7 2019年1月份反演的水位與驗潮站的水位變化對比圖

圖8采用的是2017年的數(shù)據(jù)(橫、縱坐標(biāo)同圖7),圖中上半部分曲線為SNR反演水位,下半部分曲線為驗潮站在其相應(yīng)時刻的水位變化．

(a)SNR反演水位

(b)驗潮站水位

對反演的2019年1月份31天的數(shù)據(jù)進(jìn)行粗差探測,窗口平滑,也對驗潮站數(shù)據(jù)進(jìn)行相應(yīng)的窗口平滑后對其進(jìn)行相關(guān)性分析,發(fā)現(xiàn)相關(guān)性都在0.7以上,其中相關(guān)性在0.8～0.9的天數(shù)為17天,0.9以上天數(shù)為7天,在第4年積日時出現(xiàn)過較高的相關(guān)性0.964,平均相關(guān)性是0.84．與驗潮站水位變化的均方根誤差(RMSE)為7.27 cm．

2017年全年數(shù)據(jù)反演結(jié)果與水位變化相關(guān)性系數(shù)為0.82,相關(guān)性如圖9所示,與驗潮站水位變化的RMSE為8 cm左右．

圖9 反演與驗潮站相關(guān)性

從圖7～9中可以看到,利用岸基GPS站的SNR數(shù)據(jù)獲得的水位與實際驗潮站水位具有較好的一致性,兩者的相關(guān)系數(shù)優(yōu)于0.80,RMSE在8 cm左右,驗證了GPS RMS數(shù)據(jù)能夠用于水位變化監(jiān)測理論的可行性,通過不同時段和時長的試驗分析證明了岸基GPS站可連續(xù)對水位進(jìn)行監(jiān)測[5,9,11]．

4 結(jié)束語

本文進(jìn)一步驗證了GPS信噪比數(shù)據(jù)進(jìn)行連續(xù)水位監(jiān)測的可行性,實驗結(jié)果表明該技術(shù)與實際驗潮站監(jiān)測結(jié)果具有較好的一致性．但本文實驗中尚未考慮到近岸潮汐與海況偏差等諸多因素的影響,精度有待進(jìn)一步提高．目前,此項技術(shù)不僅能應(yīng)用在潮位變化的監(jiān)測分析中,還可推廣到雪深探測,土壤濕度監(jiān)測,植被生物量監(jiān)測中,發(fā)揮更加廣泛的作用．