雷海林，陳永生，張陽陽

(東北大學(xué) 資源與土木工程學(xué)院，遼寧 沈陽 110819)

一種加快RTKLIB收斂速度的小波分析方法

雷海林，陳永生，張陽陽

(東北大學(xué) 資源與土木工程學(xué)院，遼寧 沈陽 110819)

針對RTKLIB在對變形體實時監(jiān)控的數(shù)據(jù)進(jìn)行處理時，存在收斂時間較長、收斂前坐標(biāo)振蕩值較大、得到的坐標(biāo)解不穩(wěn)定等問題，提出一種應(yīng)用多尺度分解的小波分析方法對解算獲得的坐標(biāo)時間序列進(jìn)行濾波去噪處理，著重研究小波基的選取原則、最佳分解尺度的確定方法，以及各分解層閾值參數(shù)的確定方法。實驗結(jié)果可以證明小波濾波降噪方法在縮短RTKLIB在PPP-Static模式下收斂時間的有效性，并能明顯降低收斂時間之前坐標(biāo)振蕩的幅值。

小波分析；RTKLIB；去噪；收斂時間

0 引言

隨著全球衛(wèi)星導(dǎo)航系統(tǒng)(global navigation satellite system，GNSS)的迅速發(fā)展，GNSS在變形監(jiān)測領(lǐng)域得到廣泛應(yīng)用。與變形監(jiān)測中數(shù)據(jù)的高精度、高可靠性要求相一致，GNSS數(shù)據(jù)處理方法一直是該領(lǐng)域的研究熱點。精密單點定位技術(shù)(precise point positioning，PPP)只需用單臺GNSS接收機(jī)結(jié)合國際GNSS服務(wù)(International GNSS Service，IGS)提供的精密星歷和鐘差文件即可實現(xiàn)與載波相位差分技術(shù)(real-time kinematic，RTK)相當(dāng)?shù)慕^對定位精度[1-4]。然而經(jīng)典的PPP處理軟件如GIPSY、Trip、IPAS、POSPacAIR、GrafNav、TerraPOS等，要么數(shù)據(jù)處理流程復(fù)雜，對從業(yè)者素質(zhì)要求高，要么是商業(yè)非開源軟件，獲取成本高。RTKLIB是一款由東京海洋大學(xué)的日本學(xué)者T.Takasu開發(fā)的兼容全球定位系統(tǒng)(global positioning system，GPS)、歐盟伽利略衛(wèi)星導(dǎo)航系統(tǒng)(Galileo navigation satellite system，Galileo)、中國北斗衛(wèi)星導(dǎo)航系統(tǒng)(BeiDou navigation satellite system，BDS)、格洛納斯衛(wèi)星導(dǎo)航定位系統(tǒng)(global navigation satellite system，GLONASS)、準(zhǔn)天頂衛(wèi)星系統(tǒng)(quasi-zenith satellite system，QZSS)及廣域差分增強系統(tǒng)(satellite-based augmentation system，SBAS)的開源程序包，其操作簡潔，圖形輸出功能強大，自帶的PPP-Kinematic、PPP-Static、PPP-Fixed等模塊可實現(xiàn)高精度的實時及事后單點定位[5]。當(dāng)選用PPP-Static模式時，它能夠逐歷元地輸出定位結(jié)果，可以實現(xiàn)近似的動態(tài)解；但其收斂時間較長，平均為10～15 min：這為工程應(yīng)用的變形監(jiān)測帶來了時間成本，并為某些變形體的準(zhǔn)實時監(jiān)測制造了障礙。

小波分析是一種多尺度的信號分析工具，具有帶通濾波的功能，能對信號有效分頻，從而在不同尺度上將粗差和噪聲分離，達(dá)到降噪的目的[6]。因此可將小波分析引入到對RTKLIB中PPP-Static模式下解算的坐標(biāo)時間序列中，通過濾波降噪的方法縮短PPP-Static模式下的收斂時間。本文將基于多尺度分解的思想結(jié)合工程實例，從小波分解最大尺度的確定、小波基函數(shù)的選取、閾值的確定等方面對坐標(biāo)時間序列進(jìn)行濾波去噪[7-10]。

1 小波分析的基本理論

小波分析是一種窗口面積固定但其形狀可改變的時頻局部化分析方法，通過在時域上對不同頻率的取樣步長進(jìn)行調(diào)節(jié)，可以對信號進(jìn)行由粗到細(xì)的逐步觀察，實現(xiàn)多分辨分析，被譽為“數(shù)學(xué)顯微鏡”。

1.1 小波變換的定義

小波變換的基本思想是用一族函數(shù)不斷逼近信號本身，由小波母函數(shù)經(jīng)平移因子a與伸縮因子b后生成依賴于參數(shù)(a，b)的小波序列，即

(1)

(2)

(3)

相應(yīng)離散信號的離散小波變換為

(4)

1.2 小波變換的塔式算法

尺度函數(shù)φ(t)∈V0是標(biāo)準(zhǔn)正交函數(shù)，由多分辨分析確定的正交小波子空間分解為

(5)

對應(yīng)V1?V0，W1?W0，φ(t)∈V0時，有雙尺度方程

(6)

(7)

系數(shù)分解結(jié)構(gòu)如圖1所示。

(8)

系數(shù)合成結(jié)構(gòu)如圖2所示。

2 小波濾波降噪方法與實例分析

變形監(jiān)測的噪聲信號通?？梢杂眉臃ㄐ栽肼暷Ｐ兔枋?，用函數(shù)式描述為

f(t)=s(t)+n(t)+φ(t)t=0，…，N-1。

(9)

采用塔式算法選取合適小波基對隨時間變化的含噪信號f(t)進(jìn)行j級分解與重構(gòu)時，由于工程實踐中噪聲主要集中在高頻信號上；因此可通過對尺度j的選取，分層計算不同分解層上的高頻系數(shù)閾值后重構(gòu)可實現(xiàn)信號的去噪。

2.1 小波基函數(shù)的選取

合理選用小波基是信號去噪的關(guān)鍵；但小波函數(shù)具有不唯一性，對不同的數(shù)據(jù)來源和處理目的，使用不同小波基會得到不同的分析結(jié)果。通常只能根據(jù)工程實際選用幾個小波基，比較它們在某一應(yīng)用方面的效果差異，此方法有很大的片面性；但如對這些經(jīng)典小波基在對稱性、正則性、緊支性、衰減性、消失距等數(shù)學(xué)特點方面進(jìn)行分析比較，可為工程實踐中選用小波基提供理論參考。變形監(jiān)測時間序列中的隨機(jī)項是一種無規(guī)律的變化，選取小波基應(yīng)具有正則性和高度緊支性，以滿足信號刻畫的光滑性、穩(wěn)定性及較好的局部性，dbN小波顯然滿足以上條件；但dbN小波不具有對稱性，在信號分解和重構(gòu)時產(chǎn)生一定的相位失真：因此本文選用由dbN小波構(gòu)造的近似對稱正交緊支集小波，即symN小波。

2.2 最佳尺度的確定

由式(7)及式(8)可知，最大尺度j越大，被濾掉的噪聲越大，重構(gòu)時需要置零的Dj信號也越多，信號失真也會越大。一般信噪比如果較大時選用較小的尺度j即可有效濾掉噪聲；但實測數(shù)據(jù)的信噪比一般不能預(yù)先知道：所以只能逐漸增大尺度j，根據(jù)均方根誤差(RMSE)值的變化是否趨于穩(wěn)定來確定最大尺度j。當(dāng)j分別取1、2、3,…時可得

(10)

再計算出相鄰尺度間的比值

(11)

通?？傆衦>1，隨著j逐漸增大，可認(rèn)為當(dāng)r≤1.1時，去噪與信號失真效果最佳，此時取最佳尺度j=k。

2.3 小波高頻系數(shù)的閾值選取

變形監(jiān)測信號受到階段性的強噪聲影響，分解后會在各小波層呈現(xiàn)不同強度的反映；因此可在各小波分解系數(shù)中分層確定降噪閾值。由于隨著尺度的增大，小波變換系數(shù)是逐漸變化的，而傳統(tǒng)的軟閾值法中在不同的尺度下取一個常數(shù)的做法顯然不盡合理。文獻(xiàn)[11]參照軟閾值法提出了改進(jìn)公式

(12)

(13)

即尺度為1時與軟閾值的δ一致；而當(dāng)j>1時，δ隨著尺度增大逐漸減小。式中σj可通過小波變換自適應(yīng)法求得，即在分解尺度為j時，高頻信號dj(t)主要是噪聲，根據(jù)分解的高頻信號dj(t)本身來確定σj，其中σj的計算公式為

(14)

2.4 實例分析

在撫順西露天礦北部邊坡上建設(shè)有7個監(jiān)測站、2個基準(zhǔn)站，用GPS接收機(jī)在30s采樣間隔下進(jìn)行24h連續(xù)不間斷的變形數(shù)據(jù)采集。本著從準(zhǔn)實時監(jiān)測的角度出發(fā)，本文取8007號監(jiān)測點當(dāng)前時刻前30min的觀測數(shù)據(jù)，用RTKLIB下的PPP-Static模式處理，得到變形監(jiān)測坐標(biāo)時間序列。

圖3是3個方向下坐標(biāo)時間序列減去收斂近似值后的時間收斂曲線。由圖可知X、Y方向上約在30個歷元(即15 min)后，Z方向上則在20個歷元處趨于收斂；在收斂之前，坐標(biāo)值振蕩較大，表明噪聲較明顯。

表1是取不同尺度時，由式(10)、式(11)計算的各坐標(biāo)方向上原始信號與去噪后估計信號之間的均方根誤差結(jié)果。其中RMSE(k)為均方誤差值，rk為相鄰尺度間比值。按前述最大尺度確定原則，認(rèn)為X方向上j取4，Y、Z方向上j取5最佳。

坐標(biāo)方向計算參數(shù)層數(shù)j=1j=2j=3j=4j=5j=6XRMSE(k)010160135702009022130254903383rk1335614804110101151813270YRMSE(k)0122802004036110532305558rk16319180191474010440ZRMSE(k)0116301411032160547905828rk1213022791703610640

圖4、圖5、圖6依次是X、Y、Z方向上的坐標(biāo)時間序列用sym2小波按最佳尺度j分解后的細(xì)節(jié)信號與近似信號。其中X=a4+d4+d3+d2+d1；Y=a5+d5+d4+d3+d2+d1；Z=a5+d5+d4+d3+d2+d1。

表2是顧及各方向上坐標(biāo)時間序列各分解層的小波系數(shù)分布特征，按式(12)、式(13)計算求得的X、Y、Z方向上坐標(biāo)時間序列的不同分解層上的閾值。

圖7、圖8、圖9分別是X、Y、Z方向上的坐標(biāo)時間序列經(jīng)sym2小波在最佳尺度下分解后，按小波變換系數(shù)的自適應(yīng)法分層計算閾值，重構(gòu)獲得的去噪后的時間收斂對照圖。圖中data2曲線為經(jīng)小波去噪后的新坐標(biāo)時間序列。從圖中易看出：

層數(shù)坐標(biāo)方向XYZj=1026603400305j=2024103980214j=3044408910835j=4033313401535j=503430465

X、Y、Z方向上在收斂時間之前的振蕩期，去噪后有明顯的向收斂方向逼近的趨勢，且振幅明顯減小，收斂時間均有不同程度的提前；在收斂時間后，去噪后與原始的坐標(biāo)序列基本平行，在局部有小振蕩的歷元處更光滑，不再有偏離收斂方向的振動趨勢。X方向上收斂時間有明顯提前，在第12個歷元處即開始趨于收斂，有效縮短的收斂時間約為7 min。Y方向從第4個歷元處即開始向收斂方向逼近，且振幅極小，到20個歷元處開始趨于收斂，有效縮短的收斂時間約為4 min。Z方向上因原收斂時間較X、Y方向上較早，雖收斂趨勢較為明顯，但有效提前時間并不明顯，約為2 min。

3 結(jié)束語

本文對RTKLIB處理變形監(jiān)測數(shù)據(jù)解的特點作了分析，提出一種減小坐標(biāo)振蕩幅值，并加速坐標(biāo)收斂，提高解穩(wěn)定性的小波濾波去噪方法。通過實例分析，驗證了當(dāng)小波基選取恰當(dāng)、分解尺度最佳、閾值選取合理時，此方法具有有效性。結(jié)合RTKLIB的逐歷元坐標(biāo)輸出的功能，可實現(xiàn)對變形監(jiān)測目標(biāo)的準(zhǔn)實時監(jiān)控，有廣闊的工程應(yīng)用前景。但此方法在各方向上加速收斂時間的效果不統(tǒng)一，仍有待于進(jìn)一步研究。

[1] KOUBA J，HéROUX P.Precise point positioning using IGS orbit and clock products[J].Gps Solutions，2001，5(2)：12-28.

[2] 易重海，陳永奇，朱建軍，等.一種基于IGS超快星歷的區(qū)域性實時精密單點定位方法[J].測繪學(xué)報，2011，40(2)：226-231.

[3] 王鍇華，劉鴻飛.單站單頻偽距實時動態(tài)差分定位精度及其可用性分析[J].導(dǎo)航定位學(xué)報，2014，2(2)：89-92.

[4] 劉經(jīng)南，劉暉.連續(xù)運行衛(wèi)星定位服務(wù)系統(tǒng)——城市空間數(shù)據(jù)的基礎(chǔ)設(shè)施[J].武漢大學(xué)學(xué)報·信息科學(xué)版，2003，28(3)；259-264.

[5] Tokyo University of Marine Science and Technology.RTKLIB ver.2.4.2 Manual[EB/OL].(2013-4-29)[2016-03-21].http：//www.rtklib.com/prog/manual_2.4.2.pdf.

[6] 寧津生，汪海洪，羅志才.小波分析在大地測量中的應(yīng)用及其進(jìn)展[J].武漢大學(xué)學(xué)報·信息科學(xué)版，2004，29(8)：659-663.

[7] 李宗春，鄧勇，張冠勇，等.變形測量異常數(shù)據(jù)處理中的小波變換最佳級數(shù)的確定[J].武漢大學(xué)學(xué)報·信息科學(xué)版，2011，36(3)：285-288.

[8] 尹暉，朱鋒.時序數(shù)據(jù)去噪中的小波策略及評價指標(biāo)[J].武漢大學(xué)學(xué)報·信息科學(xué)版，2012，37(11)：1374-1377.

[9] 朱建軍，章浙濤，匡翠林，等.一種可靠的小波去噪質(zhì)量評價指標(biāo)[J].武漢大學(xué)學(xué)報·信息科學(xué)版，2015，40(5)；688-694.

[10]文莉，劉正士，葛運建.小波去噪的幾種方法[J].合肥工業(yè)大學(xué)學(xué)報·自然科學(xué)版，2002，25(2)：167-172.

[11]趙瑞珍，宋國鄉(xiāng).一種基于小波變換的白噪聲消噪方法的改進(jìn)[J].西安電子科技大學(xué)學(xué)報，2000，27(5)：619-622.

A method of wavelet analysis for accelerating convergence rate of RTKLIB

LEIHailin，CHENYongsheng，ZHANGYangyang

(College of Resources and Civil Engineering，Northeastern University，Shenyang，Liaoning 110819，China)

Aiming at the problem of the long convergence time，big coordinates oscillation before convergence，and unstable coordinate solution in the data processing of real-time deformation monitoring by RTKLIB，the paper proposed a method using wavelet analysis of multi-scale decomposition to filter and denoise the resolved coordinate time series，focusing on the selection principle of the wavelet base，the determination methods of optimal decomposition scale and threshold parameters of each decomposition level.Experimental result showed the feasibility of the proposed method in shortening the convergence time with RTKLIB PPP-Static mode，and significantly reducing the amplitude of coordinates oscillation before the convergence time.

wavelet analysis；RTKLIB；denoising；convergence time

2016-07-21

雷海林(1989—)，男，湖南邵陽人，碩士生，研究方向為GNSS數(shù)據(jù)處理及理論。

雷海林，陳永生，張陽陽.一種加快RTKLIB收斂速度的小波分析方法[J].導(dǎo)航定位學(xué)報，2017，5(1)：115-120.(LEI Hailin，CHEN Yongsheng，ZHANG Yangyang.A method of wavelet analysis for accelerating convergence rate of RTKLIB[J].Journal of Navigation and Positioning，2017，5(1)：115-120.)

10.16547/j.cnki.10-1096.20170124.

P228

A

2095-4999(2017)01-0115-06