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

雷海林，陳永生，張陽陽

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

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

雷海林，陳永生，張陽陽

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

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

小波分析；RTKLIB；去噪；收斂時(shí)間

0 引言

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

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

1 小波分析的基本理論

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

1.1 小波變換的定義

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

(1)

(2)

(3)

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

(4)

1.2 小波變換的塔式算法

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

(5)

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

(6)

(7)

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

(8)

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

2 小波濾波降噪方法與實(shí)例分析

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

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

(9)

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

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

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

2.2 最佳尺度的確定

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

(10)

再計(jì)算出相鄰尺度間的比值

(11)

通常總有r>1，隨著j逐漸增大，可認(rèn)為當(dāng)r≤1.1時(shí)，去噪與信號(hào)失真效果最佳，此時(shí)取最佳尺度j=k。

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

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

(12)

(13)

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

(14)

2.4 實(shí)例分析

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

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

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

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

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

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

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

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

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

3 結(jié)束語

本文對(duì)RTKLIB處理變形監(jiān)測(cè)數(shù)據(jù)解的特點(diǎn)作了分析，提出一種減小坐標(biāo)振蕩幅值，并加速坐標(biāo)收斂，提高解穩(wěn)定性的小波濾波去噪方法。通過實(shí)例分析，驗(yàn)證了當(dāng)小波基選取恰當(dāng)、分解尺度最佳、閾值選取合理時(shí)，此方法具有有效性。結(jié)合RTKLIB的逐歷元坐標(biāo)輸出的功能，可實(shí)現(xiàn)對(duì)變形監(jiān)測(cè)目標(biāo)的準(zhǔn)實(shí)時(shí)監(jiān)控，有廣闊的工程應(yīng)用前景。但此方法在各方向上加速收斂時(shí)間的效果不統(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ū)域性實(shí)時(shí)精密單點(diǎn)定位方法[J].測(cè)繪學(xué)報(bào)，2011，40(2)：226-231.

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

[4] 劉經(jīng)南，劉暉.連續(xù)運(yùn)行衛(wèi)星定位服務(wù)系統(tǒng)——城市空間數(shù)據(jù)的基礎(chǔ)設(shè)施[J].武漢大學(xué)學(xué)報(bào)·信息科學(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] 寧津生，汪海洪，羅志才.小波分析在大地測(cè)量中的應(yīng)用及其進(jìn)展[J].武漢大學(xué)學(xué)報(bào)·信息科學(xué)版，2004，29(8)：659-663.

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

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

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

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

[11]趙瑞珍，宋國鄉(xiāng).一種基于小波變換的白噪聲消噪方法的改進(jìn)[J].西安電子科技大學(xué)學(xué)報(bào)，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—)，男，湖南邵陽人，碩士生，研究方向?yàn)镚NSS數(shù)據(jù)處理及理論。

雷海林，陳永生，張陽陽.一種加快RTKLIB收斂速度的小波分析方法[J].導(dǎo)航定位學(xué)報(bào)，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