龔旭崢，汪香梅，王凱時(shí)

（1. 浙江省測(cè)繪科學(xué)技術(shù)研究院，浙江 杭州 310023；2. 開化縣資源規(guī)劃數(shù)據(jù)中心，浙江 衢州 324300）

隨著全球衛(wèi)星導(dǎo)航系統(tǒng)（GNSS）的不斷建設(shè)與完善，目前全球多地已建立了GNSS 連續(xù)跟蹤觀測(cè)站。GNSS 測(cè)站的建立與衛(wèi)星數(shù)據(jù)觀測(cè)不僅為全球板塊運(yùn)動(dòng)、地殼變形監(jiān)測(cè)等領(lǐng)域研究提供了重要的基礎(chǔ)數(shù)據(jù)，也為國(guó)際地球參考框架的建設(shè)提供了有力的技術(shù)支撐。受多種外界復(fù)雜環(huán)境影響，GNSS 測(cè)站觀測(cè)數(shù)據(jù)通常包含大量噪聲，嚴(yán)重阻礙了GNSS 數(shù)據(jù)的進(jìn)一步分析與應(yīng)用，因此研究如何獲取更“干凈”的GNSS時(shí)間序列、構(gòu)建統(tǒng)一嚴(yán)密的GNSS時(shí)間序列噪聲模型具有重要意義。王健[1]等研究GNSS 測(cè)站數(shù)據(jù)發(fā)現(xiàn)，GNSS 數(shù)據(jù)中的噪聲主要包括有色噪聲與白噪聲，且在N、E、U 方向上的噪聲特性存在差異。GNSS 測(cè)站坐標(biāo)時(shí)間序列包含噪聲的多樣性和復(fù)雜性給降噪方法研究提出了更高的要求，目前常用的降噪方法包括奇異值分解（SVD）、經(jīng)驗(yàn)?zāi)B(tài)分解法（EMD）、小波變換法等[2-4]。SVD方法降噪的關(guān)鍵在于重構(gòu)奇異值個(gè)數(shù)的選取，但對(duì)包含趨勢(shì)項(xiàng)的信號(hào)進(jìn)行降噪時(shí)，奇異值個(gè)數(shù)選取會(huì)受趨勢(shì)分量的影響；EMD方法降噪時(shí)受端點(diǎn)效應(yīng)、模態(tài)混疊等因素影響，降噪效果受限；小波變換降噪時(shí)受分解層次、小波基、閾值的影響較大，降噪效果在參數(shù)選取不合適時(shí)較差。作為一種能對(duì)非線性、非平穩(wěn)性信號(hào)進(jìn)行自適應(yīng)分解的降噪方法，局部均值分解（LMD）不會(huì)產(chǎn)生如EMD方法一樣的負(fù)頻率。LMD方法降噪是通過(guò)重構(gòu)分解得到的低頻乘積函數(shù)（PF）分量和余量，剔除高頻分量[5]，由于沒有進(jìn)一步處理高頻分量，導(dǎo)致高頻分量中包含的有用信息丟失。因此，本文在LMD方法的基礎(chǔ)上，引入SVD方法提取高頻分量中的有用信息，建立了LMD-SVD 組合降噪模型，從而提高有用信息的利用率；并通過(guò)5 個(gè)IGS 測(cè)站坐標(biāo)時(shí)間序列進(jìn)行實(shí)驗(yàn)分析，以驗(yàn)證本文方法的有效性與優(yōu)越性。

1 基本原理與方法

1.1 LMD基本原理

假設(shè)存在一原始信號(hào)x(t)，其LMD過(guò)程為[6-7]：

1）計(jì)算原始信號(hào)中局部極值點(diǎn)ni和對(duì)應(yīng)極值x(ni)，則局部幅值ai和局部均值mi表示為：

對(duì)ai與mi進(jìn)行處理，得到包絡(luò)估計(jì)函數(shù)a11(t)與局部均值函數(shù)m11(t)。

2）將m11(t)從x(t)中剔除，則有：

通過(guò)a11(t)解調(diào)h11(t)，即

計(jì)算得到s11(t)包絡(luò)估計(jì)函數(shù)a12(t)，判斷a12(t)是否為1；若a12(t)不為1，則將s11(t)視為輸入信號(hào)，計(jì)算s11(t)局部均值函數(shù)m12(t)，從s11(t)中剔除m12(t)得到h12(t) ，再利用a12(t) 解調(diào)h12(t) ；重復(fù)上述步驟，直到得到純調(diào)頻信號(hào)s1n(t)。

最終包絡(luò)信號(hào)為：

式中，q為第q次迭代過(guò)程。

3）將純調(diào)頻信號(hào)與包絡(luò)信號(hào)的乘積表示第一個(gè)乘積函數(shù)：

4）將x(t) 中剔除PF1(t) ，得到新的信號(hào)u1(t) ，由于PF1(t)為原始信號(hào)中頻率最高分量，因此剔除PF1(t)分量后得到的u1(t)可表示為原始信號(hào)x(t)的平滑版本。重復(fù)上述步驟，直至得到一個(gè)單調(diào)函數(shù)uk(t)，停止迭代。

5）x(t)可表示為k個(gè)PF與一個(gè)殘余分量uk(t)的疊加，即

式中，p為第p個(gè)PF 分量，由基于LMD 的信號(hào)分解過(guò)程可知，LMD 分解的基礎(chǔ)為信號(hào)自身極值，不斷迭代得到相應(yīng)PF，信號(hào)最終被分解為殘余信號(hào)和部分PF。

1.2 端點(diǎn)效應(yīng)消除

信號(hào)在使用LMD方法進(jìn)行分解時(shí)，信號(hào)的分布特征會(huì)被信號(hào)的突變點(diǎn)、間斷點(diǎn)等端點(diǎn)破壞，從而使分解結(jié)果存在誤差，稱之為端點(diǎn)效應(yīng)[8]。本文使用噪聲輔助法（向原始信號(hào)中添加受控高斯白噪聲）降低端點(diǎn)效應(yīng)，則有：

式中，xi(t)為第i次加噪后信號(hào)；di(t)為第i次加入高斯白噪聲；ai為噪聲幅度；q為1、2。

對(duì)加噪后信號(hào)進(jìn)行LMD分解得到：

取N次結(jié)果的平均值，得到基于噪聲輔助法的LMD分解結(jié)果為：

式中， PFk為每個(gè)PF 分量平均結(jié)果；為每個(gè)余項(xiàng)的平均結(jié)果。

1.3 LMD-SVD信號(hào)降噪方法

SVD方法降噪的關(guān)鍵是奇異值個(gè)數(shù)的確定，本文采用奇異值差分譜準(zhǔn)則確定奇異值個(gè)數(shù)[9]，但該方法受趨勢(shì)分量影響較大，因此SVD方法適合對(duì)信號(hào)中高頻分量進(jìn)行降噪。高頻分量與低頻分量分界點(diǎn)的確定是實(shí)現(xiàn)LMD方法降噪的關(guān)鍵。本文通過(guò)計(jì)算連續(xù)均方根誤差（CMSE）方法確定分界點(diǎn)[10]，即

式中，n為PF分量個(gè)數(shù)；N為信號(hào)長(zhǎng)度。

第p階PF 分量能量密度同樣可通過(guò)式（11）表示，計(jì)算得到相鄰PF 分量CMSE 后，定義CMSE 全局極小值為：

將CMSE 全局極小值對(duì)應(yīng)的p值作為高頻分量與低頻分量的分界點(diǎn)，直接剔除高頻分量，保留低頻分量與余量，但會(huì)導(dǎo)致高頻分量中有用信息的丟失，因此本文在LMD方法的基礎(chǔ)上引入SVD，構(gòu)建LMD-SVD方法。該方法實(shí)現(xiàn)信號(hào)降噪的具體步驟為：①對(duì)某一原始信號(hào)，使用LMD方法進(jìn)行分解，得到若干個(gè)PF分量與余量；②計(jì)算相鄰PF分量之間的CMSE值，CMSE全局極小值即為高頻分量與低頻分量的分界點(diǎn)；③針對(duì)高頻分量，使用SVD方法進(jìn)行分解，得到相應(yīng)奇異值并根據(jù)奇異值差分譜準(zhǔn)則確定信號(hào)重構(gòu)奇異值個(gè)數(shù)；④重構(gòu)經(jīng)SVD方法去噪后的高頻分量、低頻分量與余量得到最終降噪信號(hào)，并對(duì)降噪效果進(jìn)行定量評(píng)價(jià)。

2 實(shí)例分析

本文選擇IGS 網(wǎng)站中LHAZ、BJFS、SHAO 等5 個(gè)IGS 站點(diǎn)2016—2019 年U 方向坐標(biāo)時(shí)間序列進(jìn)行處理與分析，坐標(biāo)時(shí)間序列采樣間隔為1 d。為降低原始序列中數(shù)據(jù)缺失以及粗差對(duì)后續(xù)數(shù)據(jù)處理與分析的影響，首先對(duì)原始坐標(biāo)時(shí)間序列進(jìn)行插補(bǔ)、粗差剔除等預(yù)處理[11]。受篇幅限制，本文僅將SHAO 站U 方向時(shí)間序列作為實(shí)驗(yàn)數(shù)據(jù)。

2.1 LMD方法坐標(biāo)時(shí)間序列降噪

受篇幅限制，本文僅以LHAZ 測(cè)站U 方向坐標(biāo)時(shí)間序列為例進(jìn)行說(shuō)明，首先通過(guò)LMD方法分解時(shí)間序列，得到8個(gè)PF分量與1個(gè)余量；然后計(jì)算相鄰PF分量之間的CMSE 值，結(jié)果見表1，可以看出，PF 分量在第3 階時(shí)出現(xiàn)CMSE 全局極小值，因此可將第3 階PF分量認(rèn)定為高頻分量與低頻分量的分界點(diǎn)；剔除前3 個(gè)高頻PF 分量，對(duì)剩余PF 分量和余量進(jìn)行重構(gòu)，得到降噪后時(shí)間序列（圖1）。

表1 相鄰PF分量之間的CMSE值/10-8

由圖1可知，LHAZ站U方向坐標(biāo)時(shí)間序列降噪后時(shí)間序列基本反映了原始時(shí)間序列的走勢(shì)，但二者之間差異較大。對(duì)原始時(shí)間序列與降噪后時(shí)間序列殘差進(jìn)行統(tǒng)計(jì)，結(jié)果見圖2。

圖2 LMD方法降噪后殘差統(tǒng)計(jì)

2.2 LMD-SVD坐標(biāo)時(shí)間序列降噪

為進(jìn)一步提取高頻噪聲中的有用信息，提高原始信號(hào)中有用信息的利用率，本文利用SVD方法對(duì)前3個(gè)PF分量進(jìn)一步降噪；再將經(jīng)SVD降噪后得到的有用信息與剩余PF分量、余量進(jìn)行重構(gòu)，獲取最終降噪時(shí)間序列。構(gòu)建高頻分量Hankel矩陣并進(jìn)行SVD，計(jì)算前100個(gè)奇異值差分譜（圖3），可以看出，在第5個(gè)位置出現(xiàn)譜峰值，而后奇異值差分譜值降低，因此選擇前5個(gè)奇異值進(jìn)行重構(gòu)得到高頻分量降噪后序列。

圖3 LMD-SVD奇異值差分譜

重構(gòu)降噪后序列、剩余PF 分量和余量得到降噪后U 方向時(shí)間序列（圖4），可以看出，LMD-SVD 方法降噪得到的時(shí)間序列不僅可以準(zhǔn)確反映原始時(shí)間序列走勢(shì)，更能縮小與原始時(shí)間序列差異。統(tǒng)計(jì)降噪后序列的殘差（圖5），計(jì)算得到殘差振幅平均絕對(duì)值為1.7 mm。相較于單一LMD 方法，LMD-SVD 方法降噪結(jié)果的殘差平均值明顯降低，能夠進(jìn)一步提取高頻PF分量中有用信息，得到更好的降噪結(jié)果。

2.3 定量分析對(duì)比

通過(guò)上述方法對(duì)LHAZ、BJFS、SHAO等5個(gè)IGS站U方向坐標(biāo)時(shí)間序列進(jìn)行降噪處理，將信噪比、相關(guān)系數(shù)和均方根誤差作為時(shí)間序列降噪效果的評(píng)價(jià)指標(biāo)。信噪比反映有用信號(hào)與噪聲的比值，值越大，表示降噪效果越好；相關(guān)系數(shù)反映降噪后信號(hào)與原始信號(hào)波形的相似度，值越大，表示降噪效果越好；均方根誤差反映降噪前后信號(hào)差異，值越小，表示降噪效果越好[12-13]。

信噪比的計(jì)算公式為[14]：

相關(guān)系數(shù)的計(jì)算公式為：

均方根誤差的計(jì)算公式為：

式中，n為信號(hào)長(zhǎng)度；σx和σx?分別為x和?的標(biāo)準(zhǔn)差；cov(x'?)為x和x?的協(xié)方差；x(i)為真實(shí)信號(hào)；為降噪后信號(hào)。

對(duì)各站點(diǎn)U 方向坐標(biāo)時(shí)間序列降噪結(jié)果進(jìn)行統(tǒng)計(jì)，結(jié)果見表2。由圖1、圖4和表2可知，LMD-SVD方法對(duì)各站點(diǎn)坐標(biāo)時(shí)間序列降噪結(jié)果的信噪比、相關(guān)系數(shù)明顯提高，均方根誤差明顯降低。統(tǒng)計(jì)5 個(gè)測(cè)站各評(píng)價(jià)指標(biāo)的平均值發(fā)現(xiàn)，相較于LMD 方法降噪結(jié)果，LMD-SVD 方法的信噪比與相關(guān)系數(shù)分別提高了34.28%與17.11%；均方根誤差降低了51.31%。綜上所述，本文提出的LMD-SVD 方法的降噪結(jié)果優(yōu)于單一LMD方法。

表2 兩種方法降噪結(jié)果統(tǒng)計(jì)

3 結(jié) 語(yǔ)

本文在LMD方法的基礎(chǔ)上引入了SVD方法，構(gòu)建了LMD-SVD降噪方法，得到的結(jié)論為：

1）LMD-SVD 降噪方法的降噪流程為：①利用LMD 方法對(duì)GNSS 測(cè)站坐標(biāo)時(shí)間序列進(jìn)行自適應(yīng)分解，并通過(guò)CMSE 局部極小值確定高低頻分量分界點(diǎn)；②利用SVD方法對(duì)包含部分有用信息的高頻分量進(jìn)行降噪，提取高頻分量中有用信息；③重構(gòu)低頻分量、余量以及降噪后高頻分量。

2）根據(jù)奇異值差分準(zhǔn)則，SVD 方法可以自適應(yīng)選擇重構(gòu)奇異值個(gè)數(shù)，實(shí)現(xiàn)高頻分量的有效重構(gòu)與降噪，由于趨勢(shì)項(xiàng)通常為低頻分量，使用奇異值差分準(zhǔn)則確定高頻分量重構(gòu)奇異值個(gè)數(shù)可有效避免趨勢(shì)項(xiàng)對(duì)奇異值個(gè)數(shù)選取的影響。

3）相較于單一LMD方法，LMD-SVD降噪方法能更加準(zhǔn)確地提取原始時(shí)間序列中的有用信息。通過(guò)對(duì)比5 個(gè)IGS 測(cè)站U 方向坐標(biāo)時(shí)間序列降噪結(jié)果發(fā)現(xiàn)，本文方法的信噪比、相關(guān)系數(shù)均明顯提高，均方根誤差明顯降低，表現(xiàn)出更好的降噪性能。