崔冰波，陳熙源

（東南大學(xué) 儀器科學(xué)與工程學(xué)院，南京 210096）

基于一種改進EMD算法的GPS多徑效應(yīng)抑制方法

崔冰波，陳熙源

（東南大學(xué) 儀器科學(xué)與工程學(xué)院，南京 210096）

為了抑制短基線雙差分GPS測量中的多徑誤差，提出了一種改進的經(jīng)驗?zāi)B(tài)分解(EMD)濾波方法。首先分析了含噪信號EMD分解的噪聲傳播特性，提出應(yīng)用噪聲輔助數(shù)據(jù)分析方法對低階模態(tài)分量進行處理，實現(xiàn)將信號的噪聲壓縮至低階本征模態(tài)函數(shù)(IMF)中，進而提高含噪信號高階模態(tài)分量的信噪比。然后基于上述改進的EMD算法實現(xiàn)了模態(tài)單元閾值降噪，隨機采樣第1階IMF構(gòu)造多個具有相同信噪比的序列，最后平均處理多個序列的降噪結(jié)果以消除EMD分解的位置敏感性誤差，與小波降噪以及傳統(tǒng)EMD閾值降噪的仿真對比表明，提出的降噪算法有顯著優(yōu)勢。最后，將提出的算法用于GPS測量信號的多徑誤差提取，結(jié)果表明該方法可以用于抑制短基線GPS測量的多徑誤差。

經(jīng)驗?zāi)B(tài)分解；多徑效應(yīng)；GPS；小波降噪

隨著微電子技術(shù)和信號處理技術(shù)的發(fā)展，基于GPS接收機內(nèi)部實時信號處理的各種濾波算法正成為多徑抑制的主要方法[1]。小波變換具有多分辨率和良好的時頻特性，近年來廣泛應(yīng)用于非平穩(wěn)信號的處理，例如應(yīng)用離散小波變換[2]和自適應(yīng)小波變換[3]分析GPS 雙差分測量信號以及基于平穩(wěn)小波變換[4]抑制GPS多徑信號都獲得了較好的效果，但小波降噪需選擇合適的小波基函數(shù)與分解層次，且針對不同的信號需選擇不同的濾波參數(shù)，故而其適用性較差。經(jīng)驗?zāi)B(tài)分解(EMD)由Huang等人提出，是一種數(shù)據(jù)驅(qū)動的自適應(yīng)時域分解方法，有很強的局部自適應(yīng)分解特性，在一定程度上克服了小波分解的不足。雖然 EMD的理論基礎(chǔ)尚不充分，但是由于其容易實現(xiàn)且總是獲得良好的效果，所以其在非線性、非平穩(wěn)信號處理中獲得了廣泛的關(guān)注。Chen等人采用EMD部分重構(gòu)法對GPS的多徑效應(yīng)進行了研究[5]，由于該方法直接將含有噪聲的模態(tài)分量去除，所以有可能在高信噪比的信號濾波中造成信號變形。Kopsinis等人基于平移不變小波變換思想提出一種基于 EMD模態(tài)單元的區(qū)間閾值濾波方法(Clear Iterative EMD Interval-Thresholding,EMD-CIIT)[6]，并且通過大量仿真實驗說明該方法對低信噪比信號的除噪效果優(yōu)于小波除噪，在一定程度上緩解了濾波中信號變形的問題，但由于其實現(xiàn)過程需要借助小波分析或者其它方法對高頻分量處理，失去了EMD自適應(yīng)分解的優(yōu)勢。為了提高EMD閾值濾波方法的適用性，本文提出了一種改進的 EMD濾波方法，并對模擬數(shù)據(jù)和 GPS觀測數(shù)據(jù)進行了仿真實驗，結(jié)果表明改進后的 EMD閾值濾波方法較小波降噪和EMD-CIIT除噪性能有顯著改善。

1 多徑效應(yīng)

GPS接收機通過接收來自衛(wèi)星（至少4顆）的廣播星歷信號，并測量信號傳輸延時，最終與光速相乘計算出與衛(wèi)星的距離來實現(xiàn)GPS的定位功能。多徑效應(yīng)是指天線接收到的信號不僅包含衛(wèi)星直射信號還包括接收機周圍反射物反射的衛(wèi)星信號，從而影響信號傳輸時間的測量結(jié)果，并最終影響GPS定位精度。目前抑制多徑效應(yīng)的方法主要包括：選擇合適的接收機工作環(huán)境、采用先進的天線和接收機硬件技術(shù)以及應(yīng)用自適應(yīng)濾波技術(shù)等等，例如遠(yuǎn)離反射源，加裝扼流圈天線，增加接收機的相關(guān)器個數(shù)等。這些技術(shù)都有一些的局限性，比如扼流圈的增加，不僅增加了接收機的成本、尺寸和重量，其對于來自天線上方的多徑信號（比如高層建筑的反射信號）效果十分有限。雙差分模型的 GPS定位技術(shù)能有效的規(guī)避系統(tǒng)的大部分具有相關(guān)特性的誤差，比如衛(wèi)星鐘差、接收機鐘差、電離層與對流層延遲等等，但是由于不同接收機的熱噪聲和多徑誤差不具有顯著的相關(guān)性，所以無法通過雙差分模型進行消除[7]。目前基于天線和接收機硬件的多徑抑制方法雖然有一定的效果，但是增加了系統(tǒng)的成本，限制了其在民用領(lǐng)域的使用范圍，使得多徑誤差仍然是短基線高精度差分定位的主要誤差[8]。無論在靜態(tài)基站應(yīng)用還是動態(tài)接收機應(yīng)用中，多徑效應(yīng)都帶有明顯的時變特性，所以固定結(jié)構(gòu)的濾波設(shè)計無法對其有效處理。基于大多數(shù)情況下多徑信號的變化緩慢的特性，對其采用EMD分析處理具有很大優(yōu)勢，可以近似認(rèn)為多徑信號存在于 EMD分解后的高階模態(tài)函數(shù)和余項中。

2 改進的EMD

EMD將信號自適應(yīng)的分解成具有物理意義的一組IMF，基于EMD的濾波方法都是基于IMF實現(xiàn)，所以改善IMF精確度是實現(xiàn)高效濾波的重要保證。目前 EMD濾波的應(yīng)用研究主要集中在部分重構(gòu)法、結(jié)合其它濾波方法（如小波變換）構(gòu)造的組合濾波方法以及改進的模態(tài)單元閾值濾波方法，上述方法均在一定程度上提高了信號的信噪比，但是忽略了EMD 分解缺陷導(dǎo)致高頻 IMF中含有少量信號微弱特征的情況，使它們不適用于高信噪比信號的降噪。

2.1 噪聲估計

為了提高 EMD對含噪信號弱特征信息的辨識度，必須對低階IMF的噪聲水平進行分析。本文提出了一種基于奇異譜分析（Singular Spectrum Analysis，SSA）的方法對IMF分量進行噪聲水平估計，SSA將一維時間序列映射到表征信號成分能量強度的多維空間[9]?；诙鄰叫盘柧徛艺鹗幍奶卣鳂?gòu)造信號：

其中，ω(i)為噪聲，模擬周期為幾十秒到幾十分鐘的正弦振蕩。

圖1 含噪信號及其各階IMF分量的奇異譜Fig.1 The SS of noisy signal and its IMFs

2.2 噪聲壓縮

為了解決EMD分解可能出現(xiàn)的小波消失、模態(tài)混疊等問題，Wu和Huang提出了基于整體平均經(jīng)驗?zāi)B(tài)分解（Ensemble Empirical Mode Decomposition,EEMD）的噪聲輔助數(shù)據(jù)分析方法[10]，使分解得到的IMF能更準(zhǔn)確反映真實物理意義，但EEMD需要進行幾十甚至幾百次 EMD的“篩選”運算，同時由于噪聲的加入使其IMF分量的個數(shù)增加，使EEMD算法的運算量和后續(xù)處理工作量增大。

為了克服EEMD的上述缺陷，胡愛軍等通過加入高頻諧波來改變信號的極值分布，從而消除模態(tài)中的異常事件[11]，該方法在異常事件比較復(fù)雜時(比如隨機噪聲干擾下)諧波參數(shù)的選擇需要進一步研究，算法適用性有待提高。仿真分析表明低階IMF的分解結(jié)果對后續(xù)分解的影響較大，減少低階IMF的分解誤差可以顯著的提高EMD分解精度?；诎自肼暤牧憔堤匦?，可以對含噪信號加入白噪聲，使更多的高頻噪聲在低階IMF中分解出來，同時為了減少運算量，只對初步估計為噪聲占優(yōu)的 IMF分量進行噪聲輔助數(shù)據(jù)分析。在確定了含噪信號y(t)的輔助分析階數(shù)M、噪聲標(biāo)準(zhǔn)差以及加噪次數(shù)I后，改進的分解算法的描述如下：

步驟①：初始化k,j，h(t)=y(t)，其中k是IMF的階數(shù),j是加噪次數(shù)；

步驟②：對h(t)加噪處理h(t)=y(t)+wn(t)，j=j+1；

步驟③：找出序列h(t)的全部極值，用三次樣條插值構(gòu)造序列的上下包絡(luò)線u(t)、d(t)，計算包絡(luò)均值m(t)=(u(t)+d(t))/2；

步驟④：求h(t)=h(t)-m(t)，判斷h(t)是否滿足停止條件，不滿足則返回步驟③繼續(xù)篩選，直至篩選出IMFkj，如滿足結(jié)束本次篩選，返回步驟②繼續(xù)加噪篩選過程直至j=I，；

圖2為EMD、EEMD以及本文算法對式(1)模擬數(shù)據(jù)s(t)處理后的結(jié)果，可以發(fā)現(xiàn)由于噪聲的加入EEMD較EMD增加了1個IMF，而本文算法較EMD減少了1階IMF，這是由于在對低階IMF進行噪聲輔助分析后，基于白噪聲的零均值特性，噪聲被壓縮到低階IMF中，使得其對信號真實信息的影響減小，從而使分解得到的IMF物理意義更明確，這為后續(xù)提出的閾值濾波方法提供基礎(chǔ)。

圖3為三種方法分解后IMF分量的奇異譜效果，右上角為第 1個奇異值占該 IMF分量能量貢獻(xiàn)占96.5%以上的放大效果圖，可以發(fā)現(xiàn)EEMD分解得到的各IMF分量的區(qū)分度優(yōu)于傳統(tǒng)EMD分解結(jié)果，紅色橢圓標(biāo)記部分包含加入噪聲引入的多余IMF分量。本文算法由于增加了低階IMF的輔助噪聲處理，使信號中的噪聲轉(zhuǎn)移到前2階，并且由于低階（1～2階）IMF噪聲特性更加明顯，奇異值差別不大，使其奇異譜變得更加平坦。由于只對低階IMF進行噪聲輔助數(shù)據(jù)分析處理，本文算法得到的高階IMF間區(qū)分度不如EEMD，但是算法的運算量明顯少于EEMD，鑒于本文主要研究混有白噪聲的 GPS雙差分信號中多徑信號的提取（假設(shè)處理的多徑信號較噪聲而言頻率很低）的目的，不對低頻的多種多徑信號進行區(qū)分，故而此處對降噪效果并無影響。

3 閾值濾波

基于小波分析的閾值除噪是對不同尺度下小波系數(shù)（常為高頻系數(shù)）的幅值進行閾值處理，由于IMF數(shù)據(jù)在過零點附近的幅值會小于任一大于零的閾值，所以幅值閾值濾波方法不適用于EMD閾值濾波。本文采用的是一種基于IMF相鄰過零點區(qū)間極值特征的閾值濾波，即考慮IMF兩過零點間的極值大小，如果極值大于閾值則保留該區(qū)間的全部幅值，否則該區(qū)間定性為噪聲數(shù)據(jù)直接舍去，避免了基于IMF幅值進行閾值濾波造成的數(shù)據(jù)變形[6]。設(shè)含噪信號s(t)經(jīng)改進后的EMD分解得到如下的形式:

式中，L為分解得到的IMF階數(shù)，r(t)為信號的余項，由于多徑信號主要集中在高階IMF中，所以可以通過下式求解除噪后的信號：

式中，fthr(·)為閾值濾波函數(shù)，m和n為需閾值處理的IMF階次上下限，即等式右邊第1項為閾值濾波后的IMF分量，其余項為信號中未經(jīng)閾值處理的分量，其中k＜m的IMF直接舍去?？紤]到EMD分解對隨機噪聲產(chǎn)生的位置具有敏感性，本文采用隨機取樣第1階IMF數(shù)據(jù)的方式，構(gòu)造多個與原數(shù)據(jù)有相同信噪比的新序列來消除位置敏感性帶來的誤差，本文沒有借助其它濾波方法處理第1階IMF，而是直接基于改進的EMD分解結(jié)果進行閾值濾波，其實現(xiàn)過程如下：

圖2 EMD、EEMD及本文算法分解效果Fig.2 Decomposition effects of EMD、EEMD and proposed method

圖3 三種方法IMF分量的奇異譜Fig.3 The SS of IMF components decomposed by the three methods

① 對s(t)進行改進的EMD分解，得到L階IMF；

② 基于除IM1F外的其余L-1階 IMF重構(gòu)信號；

③ 通過隨機采樣(隨機打亂或者循環(huán)移位)產(chǎn)生新的序列 alter(IMF1(t))，并構(gòu)造與原序列有相同信噪比的；

⑤ 重復(fù)③和④K-1次，K為設(shè)定的迭代次數(shù)即為最終降噪后的信號；

4 仿真分析

4.1 模擬數(shù)據(jù)仿真

為了驗證本文算法的效果，對四組不同頻率特征的數(shù)據(jù)進行降噪處理，前三組為matlab中wnoise函數(shù)產(chǎn)生的doppler、blocks以及heavy sine，第四組為本文式(1)建立的仿真數(shù)據(jù),序列的長度均為 1024，仿真過程分別加入不同噪聲水平的高斯白噪聲，模擬產(chǎn)生不同信噪比的含噪信號。分別采用 EMD-CIIT、小波降噪以及本文算法對信號進行降噪，三種方法均采用硬閾值法，閾值選擇采用固定閾值(sqtwolog)（雖然這種閾值方案不是最優(yōu)的，但是作為三種方法的除噪效果對比是可行的選擇）。EMD-CIIT和本文算法的迭代數(shù)均為20，其中EMD-CIIT第1階IMF的濾波選用小波降噪，本文算法對低階IMF處理中的迭代次數(shù)為200次，輔助噪聲方差0.01，噪聲輔助處理的IMF階數(shù)為2。小波降噪選擇了3種小波基，其對應(yīng)的表1中的結(jié)果，均是該小波基下基于最優(yōu)分解層次的除噪效果。本文采用兩種方法來對降噪效果進行評判，分別是均方誤差和信噪比，其表達(dá)式分別如下所述：

① 均方誤差

式中，x(i)為不含噪聲的純凈信號，x～(i)為濾波降噪后的估計信號，均方誤差越小，表明降噪效果越好。

② 信噪比

式中，分子為信號有效能量，分母為噪聲有效能量，信噪比越高表明信號殘留的噪聲越少，其單位為dB。

仿真結(jié)果如表1所示，可以看出同一種含噪信號，在不同小波基下降噪效果是不同的，對不同頻率的含噪信號，達(dá)到最優(yōu)濾波效果的小波基也是不同的，所以小波降噪的參數(shù)選取與信號本身的特征有關(guān)。

為了實現(xiàn) EMD-CIIT，本文采用小波除噪的方法對IMF1降噪，可以看到其對高頻信號中的除噪效果并不理想，這是由于在構(gòu)造同原信號有相同信噪比的序列時，如果IMF1中含有有用信息，經(jīng)過多次迭代，將含有小幅度有用信息的值疊加到高階的IMF上，造成信號的變形，所以其濾波效果同樣依賴于第1階IMF除噪時小波參數(shù)的選取，有很大的局限性。本文算法完全依賴于信號本身，不需要其他濾波方法的輔助，所以其自適應(yīng)較好，可以看到其對不同頻率的信號均有良好的處理效果。

表1 仿真數(shù)據(jù)的降噪效果Tab.1 Results of de-noising effect of simulated data

4.2 GPS信號仿真

本文使用的GPS數(shù)據(jù)采樣周期為1 s，數(shù)據(jù)長度600個歷元，整個采集過程9顆衛(wèi)星可見(PRNs 1,6,10,16,17,21,22,26,30)，選擇仰角為的PRN-1作為參考衛(wèi)星。圖4為短基線模式下多接收機系統(tǒng)差分模型測量的基線誤差，三組數(shù)據(jù)分別為三個從接收機與主接收機的基線誤差，由于是短基線測量，可以發(fā)現(xiàn)三組基線誤差數(shù)據(jù)的趨勢項基本一致，接收機的多徑誤差非常接近，可以認(rèn)為該誤差主要為多徑誤差和熱噪聲。分別采用小波降噪、EMD-CIIT以及本文算法對基線誤差進行處理，其中小波降噪選擇 sym4小波，分解層次為3，其余兩種方法的參數(shù)和上文一致，圖5顯示了應(yīng)用三種除噪方法提取信號中多徑誤差的結(jié)果,藍(lán)實線為原始誤差數(shù)據(jù)，點劃線表示 EMD-CIIT，虛線為DWT，實黑線為本文算法。圖中方框覆蓋的部分被放大至圖形下方顯示，可以發(fā)現(xiàn)EMD-CIIT能濾除一些高頻噪聲，但是其效果不如小波變換和本文算法，本文算法在信號的高頻處不僅能保留小幅度的信號波動信息，還能盡可能的濾除高頻噪聲，較小波變換能更好的保留信號的小幅波動信息，降噪后的基線誤差即為多徑誤差，可以用來補償GPS測量數(shù)據(jù)。

圖4 多接收機測量系統(tǒng)基線誤差Fig.4 Baseline errors of multi-receiver system

圖5 基線誤差的除噪效果Fig.5 De-noising effects of baseline errors

5 總 結(jié)

本文通過研究噪聲在EMD分解中的傳播特性，提出利用噪聲輔助數(shù)據(jù)分析方法將信號中的噪聲壓縮至低階IMF中，提高了噪聲與有用信息的辨識度。由于短基線 GPS雙差分測量結(jié)果中多徑信號主要是低頻信號，低頻信號降噪時，EMD部分重構(gòu)濾波方法也能獲得很好的效果，所以在改善了EMD分解的噪聲分布后，結(jié)合區(qū)間閾值濾波和部分重構(gòu)法實現(xiàn)了一種自適應(yīng)EMD閾值降噪方法。為了消除噪聲位置對濾波結(jié)果的影響，構(gòu)造了多個與原信號有相同信噪比的序列，然后采用多次濾波運算并求濾波結(jié)果均值的方法消除EMD分解的位置敏感性誤差。模擬數(shù)據(jù)和GPS數(shù)據(jù)仿真驗證了本文算法的有效性，本文算法也可用于其它非平穩(wěn)信號的除噪。

[1]Weill L R.Conquer multipath:The GPS accuracy battle,innovation[J].GPS World,1997,8(4):59-66.

[2]Liu Lilong,Wen Hongyan,Liu Bin.Mitigation of systematic errors of GPS positioning based on wavelet denoise[C]//2010 Second WRI Global Congress on Intelligent Systems.2010:253-255.

[3]Zhong P,Ding X L,Zheng D W.Adaptive wavelet transform based on cross-validation method and its application to GPS multipath mitigation[J].GPS Solut,2008,12:109-117.

[4]Mosavi M R,Azarbad M R.Multipath error mitigation based on wavelet transform in L1 GPS receivers for kinematic applications[J].Int.J.Electron.Commun.(AEU).2013,67:875- 884.

[5]Chen Jun,Shang Xiaodong.GPS multipath effect mitigation algorithm based on empirical mode decomposition[J].Earth and Space,2010:2395- 2404.

[6]Kopsinis Y,McLaughlin S.Development of EMD- based denoising methods inspired by wavelet thresholding [J].Transactions on Signal Processing,2009,57(4):1351-1362.

[7]李豹，許江寧，曹可勁，等.北斗二代導(dǎo)航信號抗多徑性能分析與仿真[J].中國慣性技術(shù)學(xué)報，2012，20(3)：339-342.LI Bao,XU Jiang-ning,CAO Ke-jin,et al.Analysis and simulation on anti-multipath performance of Beidou2 navigation[J].Journal of Chinese Inertial Technology,2012,20(3):339-342.

[8]宋超，李景森，皮壽熹.GPS局域增強系統(tǒng)基準(zhǔn)站多路徑抑制消除技術(shù)[J].中國慣性技術(shù)學(xué)報，2010，18(2)：204-208.SONG Chao,LI Jing-sen,PI Shou-xi.Technique of multipath- free to GPS local area augmentation system [J].Journal of Chinese Inertial Technology,2010,18(2):204-208.

[9]Hassani H.Singular spectrum analysis based on the minimum variance estimator[J].Nonlinear Analysis:Real World Application,2010(11):2065- 2077.

[10]Wu Zhaohua,Huang N E.Ensemble empirical mode decomposition:a noise assisted data analysis method[J].Advances in Adaptive Data Analysis,2009,1(1):1-41.

[11]胡愛軍，孫敬敬，向玲.經(jīng)驗?zāi)B(tài)分解中的模態(tài)混疊問題[J].振動、測試與診斷，2011，31(4)：429-434.HU Ai-jun,SUN Jing-jing,XIANG Ling.Mode Mixing in empirical mode decomposition[J].Journal of Vibration,Measurement &Diagnosis,2011,31(4):429-434.

[12]Flandrin P,Gon?alvès P,Rilling G.Empirical mode decomposition as a filter bank[J].IEEE Signal Process.Lett.,2004,11(2):112-114.

Multipath mitigation of GPS signal based on improved EMD algorithm

CUI Bing-bo,CHEN Xi-yuan
(School of Instrument Science and Engineering,Southeast University,Nanjing 210096,China)

An improved filter method based on EMD is proposed to mitigate the multipath effects in short baseline double-difference GPS observation.Noise-assist data analysis is adopted to compress the noise for low-order intrinsic mode function (IMF) after analyzing the propagation of noise in EMD.Signal-to-noise ratio (SNR) of the higher order IMF is increased which improves the accuracy of the EMD decomposition.A new threshold de-noising method is developed based on the proposed EMD algorithm.By altering in a random way the position of the samples of the first IMF,different noise versions of the original signal are obtained,and the position sensitive error of EMD is mitigated by averaging the multiple de-noising results.A comparative evaluation is made on simulation signals among wavelet de-noising method,conventional EMD threshold de-noising and our proposed method,and the results show that our proposed method outperforms the other two.At last,the method is applied in short baseline double difference GPS observation to separate the multipath errors and the results show that the new method is effective.

EMD;multipath effect;GPS;wavelet de-noising

U666.1

A

1005-6734(2014)03-0346-06

10.13695/j.cnki.12-1222/o3.2014.03.013

2014-01-20；

2014-04-21

國家自然科學(xué)基金（51375087，50975049）；高等學(xué)校博士學(xué)科點專項科研基金資助課題（20110092110039）,海洋公益性行業(yè)科研專項經(jīng)費項目（201205035-09）

崔冰波（1986—），男，博士研究生，從事衛(wèi)星導(dǎo)航、組合導(dǎo)航研究。E-mail：cuibingbo@163.com

聯(lián) 系 人：陳熙源（1962—），男，教授，博士生導(dǎo)師。E-mail：chxiyuan@seu.edu.cn