CORS站GNS S觀測數(shù)據(jù)處理方法研究

邵明政

（廣饒縣水利工程公司，山東 廣饒 257300）

0.引言

20世紀末我國開始進行CORS建設，CORS的主要作用包括提供位置服務、推動氣象學與地球動力學等多領域?qū)W科發(fā)展[1，2]。近些年，我國的城市建設不斷加快，各種城市化建設所需的信息需求也在不斷加大，大多數(shù)城市已經(jīng)建設滿足當?shù)匚恢梅招枨蟮腃ORS網(wǎng)并投入使用。CORS系統(tǒng)作為一種集成互聯(lián)網(wǎng)技術(shù)、數(shù)據(jù)通信技術(shù)、衛(wèi)星導航技術(shù)的實時導航定位服務系統(tǒng)，具有自動化、網(wǎng)絡化等優(yōu)勢[3]。目前，CORS系統(tǒng)在基礎測繪、工程測量、環(huán)境監(jiān)測等領域發(fā)揮著重要作用，因此，準確掌握CORS系統(tǒng)中站點觀測數(shù)據(jù)特征，評估站點在各種復雜環(huán)境下的穩(wěn)定性狀態(tài)，對于CORS站提供更好、更穩(wěn)定、更持續(xù)的位置服務具有重要意義。

作為一種對非線性、非平穩(wěn)性信號具有良好處理效果的處理方法，小波分析的出現(xiàn)打破了GNSS數(shù)據(jù)處理方法的壁壘，具有的時頻特征使其在GNSS數(shù)據(jù)領域受到越來越多的應用[4，5]。近些年，不斷有學者對小波分析進行深入研究，使其不斷優(yōu)化。奇異譜分析是一種對信號具有穩(wěn)定識別與強化功能的信號處理方法，近些年也被用于GNSS數(shù)據(jù)處理中[6]。本文充分依托小波分析在數(shù)據(jù)處理中的優(yōu)勢，將小波分析應用于GNSS觀測數(shù)據(jù)不同成分提取中，實現(xiàn)CORS站點穩(wěn)定性分析。最后將小波分析與奇異譜分析共同應用于CORS站點GNSS觀測數(shù)據(jù)趨勢項提取中，并對二者提取效果進行了對比。

1.研究方法

1.1 小波變換分析

小波變換方法自提出起就快速成為信號領域應用最為廣泛的方法之一，相比于傳統(tǒng)的Fourier變換信號處理方法，小波變換能夠取得更加顯著的效果。

1.1.1連續(xù)小波變換

假設存在Ψ（t）∈L2（R），其Fourier變換為φ^(ω)，當φ^(ω)滿足完全重構(gòu)或者恒等分分辨率條件時[7]：

這說明式（1）是有界的，此時Ψ（t）為基本小波，伸縮平移基本小波后得到新的小波序列：

式（2）中，a、b分別為小波序列的伸縮因子與平移因子。

連續(xù)小波變換主要有以下幾個方面的性質(zhì)：

（1）線性性質(zhì)

（2）時標定理

（3）時移共變性

（4）能量守恒

小波變換過程中不增加信號能量，信號能量始終保持穩(wěn)定。

1.1.2離散小波變換

連續(xù)小波變換中，a與b不斷變化，為減少冗余，可進行離散化處理。離散方式是對a與b作冪級數(shù)，分別取a=ajo，b=kajob0，其中j∈Z，擴展步長a0>1，離散小波函數(shù)可表示為：

重構(gòu)公式為：

式（9）中的常數(shù)c與信號無關(guān)。

1.1.3小波分析分解層次的確定

通常來說，在進行信號小波分析時，信號中的噪聲會隨著分解層次的增加被更多地濾除。然而，分解層次在一定范圍內(nèi)是可行的，超過一定范圍會將信號中真實信號剔除，造成信號缺失。因此，在進行小波分析時分解層次尤為關(guān)鍵，直接影響小波分析去噪效果。信號的信噪比對于小波分析分解層次的影響較大，若信噪比較小，信號中包含的噪聲過多，此時若將分解層次設置過小會剔除過多噪聲。

在對CORS站監(jiān)測數(shù)據(jù)進行小波去噪時，信噪比無法得知，因此若要得到最優(yōu)分解層次須通過估計的方式。首先將分解層次定為1，在此基礎上依次累加，計算在不同分解層次下去噪前后信號均方根誤差得到最優(yōu)分解層次：

式（10）中，f（n）與fj（n）分別表示去噪前后信號，計算分解層次為j+1與分解層次為j時的均方根誤差比：

當存在最接近1的r值時，表明噪聲已被全部濾除，此時分解層次可定為j或者j+1。

1.2 奇異譜分析

奇異譜分析是由多元統(tǒng)計方法、時間序列分析方法等多種技術(shù)方法組成的，能夠?qū)Ψ欠€(wěn)定性信號進行有效處理，提取得到信號中包含的振動信息與變形信息，實現(xiàn)變形體變形趨勢的有效判斷[8]。

實現(xiàn)奇異譜分析功能的主要途徑為重建成分RC（Recons truction Component），假設存在某序列xi，該序列的重建成分RC可由T-EOF與T-PC組成，可表示為：

原始信號可以看成是若干個重建成分RC的集合，即將所有重建成分RC相加就是原始信號序列，表示為：

奇異譜分析中，若要判斷重建成分RC是否為趨勢項，可通過Kendall非參數(shù)檢驗的方式，主要步驟為：

（1）當滿足xi，k

（2）若檢驗得到xk為非趨勢項，表明τ服從正態(tài)分布，均值為0、均方差為S：

（3）置信度α取0.05，當滿足τ<-1.96S或τ>1.96S時，表明第k個重建成分RC為趨勢項成分。

2.CORS坐標序列解算

作為一款GNSS數(shù)據(jù)解算領域較為成熟穩(wěn)定的軟件，GAMIT/GLOBK軟件自問世以來就得到了廣泛應用。本文選擇GAMIT/GLOBK軟件對CORS站點觀測數(shù)據(jù)進行解算，解算模型為雙差模型，通過雙差模型與衛(wèi)星星歷可達到毫米級坐標解算精度。

通過GAMIT和GLOBK這兩個模塊，GAMIT/GLOBK軟件可快速高精度實現(xiàn)GNSS觀測數(shù)據(jù)解算，通過GAMIT得到單天解h文件并通過GLOBK對h文件進行后處理，基于GAMIT/GLOBK軟件進行GNSS數(shù)據(jù)解算的主要流程，如圖1所示。

圖1基于GAMIT/GLOBK軟件的GNSS數(shù)據(jù)解算

實驗數(shù)據(jù)選擇BRTD站點2015年1月1日至2018年12月31日觀測GNSS數(shù)據(jù)，通過GAMIT/GLOBK軟件解算得到實驗數(shù)據(jù)單日松弛解，使用標準化均方根殘差NRMS與基線重復性對解算結(jié)果進行評價。得到NRMS最小值與最大值分別為0.138 31與0.197 22，相關(guān)基線解算規(guī)范中規(guī)定NRMS值小于0.25[9]，結(jié)果表明GAMIT基線解算結(jié)果滿足精度指標。

計算得到基線重復性常數(shù)部分與比例系數(shù)如表1所示。通過表1可知：基線解算重復性滿足精度要求。

表1基線向量重復性

GAMIT基線解算完成后，通過GLOBK軟件得到CORS站點在ITRF 2005坐標系框架系的坐標時間序列如圖2所示。

圖2 CORS站ITRF 2015坐標系下解算坐標

通過圖2可知：N、E方向上坐標時間序列呈現(xiàn)的是線性變化，并且在N方向隨時間推遲上逐漸遞減，在E方向上隨時間推遲逐漸遞增，U方向上的變化趨勢為一種周期性并伴跳躍。

3.結(jié)果與分析

3.1 趨勢項成分提取

作為表述外界環(huán)境（如，地質(zhì)變化、土地變化）的一種狀態(tài)，趨勢項變化是GNSS觀測數(shù)據(jù)中是一個重要組成部分。通過提取CORS站監(jiān)測數(shù)據(jù)中的趨勢項，可以實現(xiàn)站點的趨勢性變形特征。本文通過小波分析方法提取站點觀測數(shù)據(jù)中的趨勢項，首先對站點3個方向的坐標時間序列進行小波分析，計算不同分解層次下的均方根誤差比值，當分解層次為9與8時比值為1.06，可將分解層次定為8；以U方向時間序列為例進行小波分解，U方向小波分解逼近信號與細節(jié)信號分別如圖3、圖4所示。

圖3 U方向逼近信號

圖4 U方向細節(jié)信號

通過圖3、圖4可知：信號的波動隨著分解層次的增加逐漸緩和，本文小波分析方法中選擇db4小波基進行8層分解。3個方向經(jīng)重構(gòu)后得到的重構(gòu)信號，如圖5所示，其中藍色細線表示原始時間序列，紅色粗線表示經(jīng)小波分析方法提取趨勢項成分。方向偏移與向北偏移；在U方向上呈現(xiàn)向上偏移。相較于U方向，N、E方向的偏移量更大并且3方向的偏移速率均有所不同。

3.2 周期項成分提取

周期性變化是變形體在一定時間內(nèi)呈現(xiàn)的規(guī)律性的往返變化，CORS站點GNSS觀測數(shù)據(jù)在外界環(huán)境影響下也存在穩(wěn)定周期性變化。作為一種能夠?qū)τ^測時間序列中的周期項成分進行準確提取的方法，本文選擇使用功率譜估計方法提取CORS站點3方向周期項成分。該方法包括兩個部分，先進行傅里葉變換，再使用變換后信號除以樣本量，可表示為[10]：

上式中，N表示時間序列長度；f表示時間序列頻率。

觀測數(shù)據(jù)經(jīng)過傅里葉變換后會表現(xiàn)出周期性變化，將周期性變化制作得到周期圖（PSD），CORS站點3方向上的PSD結(jié)果如圖6所示。

圖5 3方向趨勢項成分提取結(jié)果

通過圖5可知：經(jīng)小波分析提取趨勢項與原始時間序列的變形趨勢基本保持一致，能夠?qū)υ甲冃乌厔葸M行反映。

對經(jīng)小波分析提取的變形趨勢進行量化，得到變形趨勢如表2所示。

表2 CORS站點提取變形趨勢

通過表2可知：CORS站點在N、E方向上分別呈現(xiàn)向西

圖6 CORS站點3方向周期圖

通過圖6可知：N、E方向譜峰值均為292 d、512 d、682 d；U方向譜峰值分別為186 d、341 d。作為功率譜估計中的關(guān)鍵參數(shù)，窗口長度的選取對于譜估計結(jié)果影響較大，一般取在序列長度三分之一以內(nèi)。

計算在不同窗口長度下經(jīng)小波重構(gòu)后信號與原始信號均方根誤差的變化如圖7所示。

圖7 CORS站點3方向重構(gòu)信號與原始信號均方根誤差隨窗口長度變化

通過圖7可知：當N、E、U方向的窗口長度分別為292 d、292 d、341 d時，均方根變化圖會產(chǎn)生拐點，此對應均方根誤差分別為0.095 12 mm、0.105 7 mm、0.002 147 mm。統(tǒng)計3方向上的周期性如表3所示。

表3 CORS站點3個方向的周期

通過上述圖表得知：CORS站點觀測數(shù)據(jù)在U方向的周期性要強于N、E方向。

3.3 與奇異譜分析方法的比較

將奇異譜分析應用到CORS站點觀測時間序列趨勢項提取中，并將提取結(jié)果與小波分析提取結(jié)果進行對比，結(jié)果如圖8所示。

圖8 CORS站點3方向小波分析與奇異譜分析提取結(jié)果

通過圖8可知：小波分析方法與奇異譜分析方法均不同程度將原始序列中包含噪聲剔除，提取得到時間序列主要變形趨勢。然而相較于奇異譜分析，小波分析去噪后序列更加平滑，反映出更多的時間序列特征。

4.結(jié)束語

本文通過對CORS站點觀測GNSS數(shù)據(jù)進行解算以及趨勢項、周期項提取，結(jié)果表明：使用GAMIT/GLOBK軟件可快速高精度實現(xiàn)GNSS觀測數(shù)據(jù)解算，解算結(jié)果的標準化均方根殘差NRMS與基線重復性均滿足規(guī)范精度要求。小波分析提取CORS站點GNSS觀測數(shù)據(jù)N、E方向上有-0.0326 mm/d的負向趨勢運動、0.0779 mm/d的正向趨勢運動；在U方向上有0.0026 mm/d的正向趨勢運動。通過小波分析結(jié)合周期圖提取CORS站點GNSS觀測數(shù)據(jù)周期項，表明在N、E方向上年周期較為明顯；在U方向上季度、半年、年周期較為明顯。最后將小波分析與奇異譜分析提取CORS站觀測GNSS數(shù)據(jù)趨勢項結(jié)果進行對比，表明兩種方法均能將時間序列主要特征表現(xiàn)出來，但是小波分析提取序列更加平滑，去噪效果更好。