邵明政
(廣饒縣水利工程公司,山東 廣饒 257300)
20世紀(jì)末我國開始進(jìn)行CORS建設(shè),CORS的主要作用包括提供位置服務(wù)、推動氣象學(xué)與地球動力學(xué)等多領(lǐng)域?qū)W科發(fā)展[1,2]。近些年,我國的城市建設(shè)不斷加快,各種城市化建設(shè)所需的信息需求也在不斷加大,大多數(shù)城市已經(jīng)建設(shè)滿足當(dāng)?shù)匚恢梅?wù)需求的CORS網(wǎng)并投入使用。CORS系統(tǒng)作為一種集成互聯(lián)網(wǎng)技術(shù)、數(shù)據(jù)通信技術(shù)、衛(wèi)星導(dǎo)航技術(shù)的實(shí)時(shí)導(dǎo)航定位服務(wù)系統(tǒng),具有自動化、網(wǎng)絡(luò)化等優(yōu)勢[3]。目前,CORS系統(tǒng)在基礎(chǔ)測繪、工程測量、環(huán)境監(jiān)測等領(lǐng)域發(fā)揮著重要作用,因此,準(zhǔn)確掌握CORS系統(tǒng)中站點(diǎn)觀測數(shù)據(jù)特征,評估站點(diǎn)在各種復(fù)雜環(huán)境下的穩(wěn)定性狀態(tài),對于CORS站提供更好、更穩(wěn)定、更持續(xù)的位置服務(wù)具有重要意義。
作為一種對非線性、非平穩(wěn)性信號具有良好處理效果的處理方法,小波分析的出現(xiàn)打破了GNSS數(shù)據(jù)處理方法的壁壘,具有的時(shí)頻特征使其在GNSS數(shù)據(jù)領(lǐng)域受到越來越多的應(yīng)用[4,5]。近些年,不斷有學(xué)者對小波分析進(jìn)行深入研究,使其不斷優(yōu)化。奇異譜分析是一種對信號具有穩(wěn)定識別與強(qiáng)化功能的信號處理方法,近些年也被用于GNSS數(shù)據(jù)處理中[6]。本文充分依托小波分析在數(shù)據(jù)處理中的優(yōu)勢,將小波分析應(yīng)用于GNSS觀測數(shù)據(jù)不同成分提取中,實(shí)現(xiàn)CORS站點(diǎn)穩(wěn)定性分析。最后將小波分析與奇異譜分析共同應(yīng)用于CORS站點(diǎn)GNSS觀測數(shù)據(jù)趨勢項(xiàng)提取中,并對二者提取效果進(jìn)行了對比。
小波變換方法自提出起就快速成為信號領(lǐng)域應(yīng)用最為廣泛的方法之一,相比于傳統(tǒng)的Fourier變換信號處理方法,小波變換能夠取得更加顯著的效果。
1.1.1連續(xù)小波變換
假設(shè)存在Ψ(t)∈L2(R),其Fourier變換為φ^(ω),當(dāng)φ^(ω)滿足完全重構(gòu)或者恒等分分辨率條件時(shí)[7]:
這說明式(1)是有界的,此時(shí)Ψ(t)為基本小波,伸縮平移基本小波后得到新的小波序列:
式(2)中,a、b分別為小波序列的伸縮因子與平移因子。
連續(xù)小波變換主要有以下幾個(gè)方面的性質(zhì):
(1)線性性質(zhì)
(2)時(shí)標(biāo)定理
(3)時(shí)移共變性
(4)能量守恒
小波變換過程中不增加信號能量,信號能量始終保持穩(wěn)定。
1.1.2離散小波變換
連續(xù)小波變換中,a與b不斷變化,為減少冗余,可進(jìn)行離散化處理。離散方式是對a與b作冪級數(shù),分別取a=ajo,b=kajob0,其中j∈Z,擴(kuò)展步長a0>1,離散小波函數(shù)可表示為:
重構(gòu)公式為:
式(9)中的常數(shù)c與信號無關(guān)。
1.1.3小波分析分解層次的確定
通常來說,在進(jìn)行信號小波分析時(shí),信號中的噪聲會隨著分解層次的增加被更多地濾除。然而,分解層次在一定范圍內(nèi)是可行的,超過一定范圍會將信號中真實(shí)信號剔除,造成信號缺失。因此,在進(jìn)行小波分析時(shí)分解層次尤為關(guān)鍵,直接影響小波分析去噪效果。信號的信噪比對于小波分析分解層次的影響較大,若信噪比較小,信號中包含的噪聲過多,此時(shí)若將分解層次設(shè)置過小會剔除過多噪聲。
在對CORS站監(jiān)測數(shù)據(jù)進(jìn)行小波去噪時(shí),信噪比無法得知,因此若要得到最優(yōu)分解層次須通過估計(jì)的方式。首先將分解層次定為1,在此基礎(chǔ)上依次累加,計(jì)算在不同分解層次下去噪前后信號均方根誤差得到最優(yōu)分解層次:
式(10)中,f(n)與fj(n)分別表示去噪前后信號,計(jì)算分解層次為j+1與分解層次為j時(shí)的均方根誤差比:
當(dāng)存在最接近1的r值時(shí),表明噪聲已被全部濾除,此時(shí)分解層次可定為j或者j+1。
奇異譜分析是由多元統(tǒng)計(jì)方法、時(shí)間序列分析方法等多種技術(shù)方法組成的,能夠?qū)Ψ欠€(wěn)定性信號進(jìn)行有效處理,提取得到信號中包含的振動信息與變形信息,實(shí)現(xiàn)變形體變形趨勢的有效判斷[8]。
實(shí)現(xiàn)奇異譜分析功能的主要途徑為重建成分RC(Recons truction Component),假設(shè)存在某序列xi,該序列的重建成分RC可由T-EOF與T-PC組成,可表示為:
原始信號可以看成是若干個(gè)重建成分RC的集合,即將所有重建成分RC相加就是原始信號序列,表示為:
奇異譜分析中,若要判斷重建成分RC是否為趨勢項(xiàng),可通過Kendall非參數(shù)檢驗(yàn)的方式,主要步驟為:
(1)當(dāng)滿足xi,k (2)若檢驗(yàn)得到xk為非趨勢項(xiàng),表明τ服從正態(tài)分布,均值為0、均方差為S: (3)置信度α取0.05,當(dāng)滿足τ<-1.96S或τ>1.96S時(shí),表明第k個(gè)重建成分RC為趨勢項(xiàng)成分。 作為一款GNSS數(shù)據(jù)解算領(lǐng)域較為成熟穩(wěn)定的軟件,GAMIT/GLOBK軟件自問世以來就得到了廣泛應(yīng)用。本文選擇GAMIT/GLOBK軟件對CORS站點(diǎn)觀測數(shù)據(jù)進(jìn)行解算,解算模型為雙差模型,通過雙差模型與衛(wèi)星星歷可達(dá)到毫米級坐標(biāo)解算精度。 通過GAMIT和GLOBK這兩個(gè)模塊,GAMIT/GLOBK軟件可快速高精度實(shí)現(xiàn)GNSS觀測數(shù)據(jù)解算,通過GAMIT得到單天解h文件并通過GLOBK對h文件進(jìn)行后處理,基于GAMIT/GLOBK軟件進(jìn)行GNSS數(shù)據(jù)解算的主要流程,如圖1所示。 圖1基于GAMIT/GLOBK軟件的GNSS數(shù)據(jù)解算 實(shí)驗(yàn)數(shù)據(jù)選擇BRTD站點(diǎn)2015年1月1日至2018年12月31日觀測GNSS數(shù)據(jù),通過GAMIT/GLOBK軟件解算得到實(shí)驗(yàn)數(shù)據(jù)單日松弛解,使用標(biāo)準(zhǔn)化均方根殘差NRMS與基線重復(fù)性對解算結(jié)果進(jìn)行評價(jià)。得到NRMS最小值與最大值分別為0.138 31與0.197 22,相關(guān)基線解算規(guī)范中規(guī)定NRMS值小于0.25[9],結(jié)果表明GAMIT基線解算結(jié)果滿足精度指標(biāo)。 計(jì)算得到基線重復(fù)性常數(shù)部分與比例系數(shù)如表1所示。通過表1可知:基線解算重復(fù)性滿足精度要求。 表1基線向量重復(fù)性 GAMIT基線解算完成后,通過GLOBK軟件得到CORS站點(diǎn)在ITRF 2005坐標(biāo)系框架系的坐標(biāo)時(shí)間序列如圖2所示。 圖2 CORS站ITRF 2015坐標(biāo)系下解算坐標(biāo) 通過圖2可知:N、E方向上坐標(biāo)時(shí)間序列呈現(xiàn)的是線性變化,并且在N方向隨時(shí)間推遲上逐漸遞減,在E方向上隨時(shí)間推遲逐漸遞增,U方向上的變化趨勢為一種周期性并伴跳躍。 作為表述外界環(huán)境(如,地質(zhì)變化、土地變化)的一種狀態(tài),趨勢項(xiàng)變化是GNSS觀測數(shù)據(jù)中是一個(gè)重要組成部分。通過提取CORS站監(jiān)測數(shù)據(jù)中的趨勢項(xiàng),可以實(shí)現(xiàn)站點(diǎn)的趨勢性變形特征。本文通過小波分析方法提取站點(diǎn)觀測數(shù)據(jù)中的趨勢項(xiàng),首先對站點(diǎn)3個(gè)方向的坐標(biāo)時(shí)間序列進(jìn)行小波分析,計(jì)算不同分解層次下的均方根誤差比值,當(dāng)分解層次為9與8時(shí)比值為1.06,可將分解層次定為8;以U方向時(shí)間序列為例進(jìn)行小波分解,U方向小波分解逼近信號與細(xì)節(jié)信號分別如圖3、圖4所示。 圖3 U方向逼近信號 圖4 U方向細(xì)節(jié)信號 通過圖3、圖4可知:信號的波動隨著分解層次的增加逐漸緩和,本文小波分析方法中選擇db4小波基進(jìn)行8層分解。3個(gè)方向經(jīng)重構(gòu)后得到的重構(gòu)信號,如圖5所示,其中藍(lán)色細(xì)線表示原始時(shí)間序列,紅色粗線表示經(jīng)小波分析方法提取趨勢項(xiàng)成分。方向偏移與向北偏移;在U方向上呈現(xiàn)向上偏移。相較于U方向,N、E方向的偏移量更大并且3方向的偏移速率均有所不同。 周期性變化是變形體在一定時(shí)間內(nèi)呈現(xiàn)的規(guī)律性的往返變化,CORS站點(diǎn)GNSS觀測數(shù)據(jù)在外界環(huán)境影響下也存在穩(wěn)定周期性變化。作為一種能夠?qū)τ^測時(shí)間序列中的周期項(xiàng)成分進(jìn)行準(zhǔn)確提取的方法,本文選擇使用功率譜估計(jì)方法提取CORS站點(diǎn)3方向周期項(xiàng)成分。該方法包括兩個(gè)部分,先進(jìn)行傅里葉變換,再使用變換后信號除以樣本量,可表示為[10]: 上式中,N表示時(shí)間序列長度;f表示時(shí)間序列頻率。 觀測數(shù)據(jù)經(jīng)過傅里葉變換后會表現(xiàn)出周期性變化,將周期性變化制作得到周期圖(PSD),CORS站點(diǎn)3方向上的PSD結(jié)果如圖6所示。 圖5 3方向趨勢項(xiàng)成分提取結(jié)果 通過圖5可知:經(jīng)小波分析提取趨勢項(xiàng)與原始時(shí)間序列的變形趨勢基本保持一致,能夠?qū)υ甲冃乌厔葸M(jìn)行反映。 對經(jīng)小波分析提取的變形趨勢進(jìn)行量化,得到變形趨勢如表2所示。 表2 CORS站點(diǎn)提取變形趨勢 通過表2可知:CORS站點(diǎn)在N、E方向上分別呈現(xiàn)向西 圖6 CORS站點(diǎn)3方向周期圖 通過圖6可知:N、E方向譜峰值均為292 d、512 d、682 d;U方向譜峰值分別為186 d、341 d。作為功率譜估計(jì)中的關(guān)鍵參數(shù),窗口長度的選取對于譜估計(jì)結(jié)果影響較大,一般取在序列長度三分之一以內(nèi)。 計(jì)算在不同窗口長度下經(jīng)小波重構(gòu)后信號與原始信號均方根誤差的變化如圖7所示。 圖7 CORS站點(diǎn)3方向重構(gòu)信號與原始信號均方根誤差隨窗口長度變化 通過圖7可知:當(dāng)N、E、U方向的窗口長度分別為292 d、292 d、341 d時(shí),均方根變化圖會產(chǎn)生拐點(diǎn),此對應(yīng)均方根誤差分別為0.095 12 mm、0.105 7 mm、0.002 147 mm。統(tǒng)計(jì)3方向上的周期性如表3所示。 表3 CORS站點(diǎn)3個(gè)方向的周期 通過上述圖表得知:CORS站點(diǎn)觀測數(shù)據(jù)在U方向的周期性要強(qiáng)于N、E方向。 將奇異譜分析應(yīng)用到CORS站點(diǎn)觀測時(shí)間序列趨勢項(xiàng)提取中,并將提取結(jié)果與小波分析提取結(jié)果進(jìn)行對比,結(jié)果如圖8所示。 圖8 CORS站點(diǎn)3方向小波分析與奇異譜分析提取結(jié)果 通過圖8可知:小波分析方法與奇異譜分析方法均不同程度將原始序列中包含噪聲剔除,提取得到時(shí)間序列主要變形趨勢。然而相較于奇異譜分析,小波分析去噪后序列更加平滑,反映出更多的時(shí)間序列特征。 本文通過對CORS站點(diǎn)觀測GNSS數(shù)據(jù)進(jìn)行解算以及趨勢項(xiàng)、周期項(xiàng)提取,結(jié)果表明:使用GAMIT/GLOBK軟件可快速高精度實(shí)現(xiàn)GNSS觀測數(shù)據(jù)解算,解算結(jié)果的標(biāo)準(zhǔn)化均方根殘差NRMS與基線重復(fù)性均滿足規(guī)范精度要求。小波分析提取CORS站點(diǎn)GNSS觀測數(shù)據(jù)N、E方向上有-0.0326 mm/d的負(fù)向趨勢運(yùn)動、0.0779 mm/d的正向趨勢運(yùn)動;在U方向上有0.0026 mm/d的正向趨勢運(yùn)動。通過小波分析結(jié)合周期圖提取CORS站點(diǎn)GNSS觀測數(shù)據(jù)周期項(xiàng),表明在N、E方向上年周期較為明顯;在U方向上季度、半年、年周期較為明顯。最后將小波分析與奇異譜分析提取CORS站觀測GNSS數(shù)據(jù)趨勢項(xiàng)結(jié)果進(jìn)行對比,表明兩種方法均能將時(shí)間序列主要特征表現(xiàn)出來,但是小波分析提取序列更加平滑,去噪效果更好。2.CORS坐標(biāo)序列解算
3.結(jié)果與分析
3.1 趨勢項(xiàng)成分提取
3.2 周期項(xiàng)成分提取
3.3 與奇異譜分析方法的比較
4.結(jié)束語