基于CEEMDAN-MPE算法的隧道爆破地震波信號(hào)降噪方法及應(yīng)用*

黃智剛，呂虎波，林一庚，彭亞雄，吳 立，陳 勁

(1.中國地質(zhì)大學(xué)(武漢) 工程學(xué)院，武漢 430074；2.福州水務(wù)平潭引水開發(fā)有限公司，福州 350001；3.浙江省隧道工程集團(tuán)有限公司，杭州 310030；4.福州城建設(shè)計(jì)研究院有限公司，福州 350001;5.湖南科技大學(xué) 巖土工程穩(wěn)定控制與健康監(jiān)測湖南省重點(diǎn)實(shí)驗(yàn)室，湘潭 411201)

由于工程環(huán)境復(fù)雜、電磁干擾和監(jiān)測儀器誤差等因素影響，實(shí)測爆破地震波信號(hào)包含大量高頻噪聲，掩蓋了地震波真實(shí)信息，直接影響隧道爆破有害效應(yīng)分析與評價(jià)。為了準(zhǔn)確掌握隧道爆破地震波波形特征、能量特性和衰減規(guī)律，必須對實(shí)測爆破地震波信號(hào)進(jìn)行降噪處理。小波算法具備較好的時(shí)頻局域化特征，利用其進(jìn)行信號(hào)降噪處理是一種廣泛使用的方法[1]。熊正明等利用平移不變小波對爆破振動(dòng)信號(hào)進(jìn)行去噪處理[2]，消除信號(hào)的偽吉布斯現(xiàn)象，同時(shí)減小降噪后信號(hào)與原始信號(hào)的誤差。路亮等提出了基于提升小波包最優(yōu)基分解算法的爆破振動(dòng)信號(hào)的降噪和能量提取方法[3]，驗(yàn)證了方法的有效性。由于小波變換算法降噪過程中小波基函數(shù)和分解層次難以確定，使得這類方法的自適應(yīng)性不強(qiáng)，降噪效果難以保證[4,5]。經(jīng)驗(yàn)?zāi)B(tài)分解(EMD)是一種處理非平穩(wěn)信號(hào)的分解方法，對爆破振動(dòng)信號(hào)分解有較好的適應(yīng)性[6]。費(fèi)鴻祿等將改進(jìn)EMD和小波閾值算法結(jié)合進(jìn)行降噪處理[7]，較好地去除了爆破振動(dòng)信號(hào)所含噪聲。

自適應(yīng)噪聲的完全集合經(jīng)驗(yàn)?zāi)B(tài)分解算法(CEEMDAN)是一種基于EMD的改進(jìn)算法，能夠消除人為添加噪聲對原始信號(hào)完備性的影響，抑制了模態(tài)混疊問題又避免了原始信號(hào)失真[8]。本文通過對信號(hào)進(jìn)行CEEMDAN分解，利用多尺度排列熵(MPE)檢測分解得到的模態(tài)函數(shù)(IMF)的隨機(jī)性，去除噪聲IMF分量以達(dá)到信號(hào)降噪的目的，構(gòu)建了一種適合于隧道爆破地震波信號(hào)的降噪方法。將該方法應(yīng)用于福建省平潭及閩江口水資源配置工程，對實(shí)測隧道爆破振動(dòng)信號(hào)進(jìn)行降噪處理，并驗(yàn)證了方法的有效性。

1 信號(hào)降噪算法

1.1 CEEMDAN

經(jīng)驗(yàn)?zāi)B(tài)分解算法(EMD)根據(jù)信號(hào)的時(shí)標(biāo)特性，將多分量信號(hào)分解為一系列固有模態(tài)函數(shù)分量和剩余分量，并按瞬時(shí)頻率由高到低的順序排列，具有良好的適應(yīng)性、完備性和正交性[9]。然而該方法在處理含有不連續(xù)、脈沖和噪聲的信號(hào)時(shí)存在模態(tài)混合問題。Torres等人對EMD算法進(jìn)行改進(jìn)[10]，提出了自適應(yīng)噪聲的完全集合經(jīng)驗(yàn)?zāi)Ｊ椒纸馑惴?CEEMDAN)。該算法在EMD分解各階段自適應(yīng)添加白噪聲，計(jì)算唯一的殘差信號(hào)以獲取固有模態(tài)函數(shù)(IMF)，能夠在集成次數(shù)較少的情況下，使得重構(gòu)誤差幾乎為零，重構(gòu)信號(hào)與原信號(hào)幾乎完全相同，并在一定程度上解決了EMD算法的模態(tài)混疊現(xiàn)象[11]。CEEMDAN的主要步驟如下。

在原始信號(hào)x(t)中添加不同幅值的白噪聲nj(t)，可表示為x(t)+ε0nj(t)，其中ε0為噪聲系數(shù)。利用EMD對加噪信號(hào)進(jìn)行I次分解，通過集成平均得到第一個(gè)IMF分量。IMF分量和殘差分量如下所示

(1)

r1(t)=x(t)-IMF1(t)

(2)

定義EMDj(·)是EMD分解的第j個(gè)模態(tài)函數(shù)。對加噪信號(hào)r1(t)+ε1·EMD1[nj(t)]進(jìn)行分解I次分解，得到第二個(gè)IMF分量

(3)

計(jì)算k階殘差分量

rk(t)=rk-1(t)-IMFk(t)

(4)

從r1(t)+ε1·EMD1[nj(t)]中提取第一個(gè)IMF，得到IMFk+1。

(5)

重復(fù)上述計(jì)算直到殘差分量不能繼續(xù)分解，得到所有的IMF分量。

(6)

則原始信號(hào)x(t)可以表示為

(7)

CEEMDAN算法利用了噪聲輔助分析技術(shù)，能夠完整地重構(gòu)原始信號(hào)。針對不同類型信號(hào)，利用噪聲系數(shù)ε加入不同信噪比的白噪音，能夠有效地提高分解效果。

1.2 多尺度排列熵

多尺度排列熵(MPE)是一種檢測信號(hào)隨機(jī)性和動(dòng)力突變的方法，將時(shí)間序列進(jìn)行多尺度粗?；?，進(jìn)而計(jì)算其排列熵[12]。具體步驟如下：

①對時(shí)間序列X={x1,x2,...,xL}進(jìn)行多尺度粗?；幚?/p>

(8)

(9)

式中：τ為時(shí)間延遲；m為嵌入維數(shù)。

(10)

(11)

(12)

⑤對上述計(jì)算的排列熵進(jìn)行歸一化處理

(13)

CEEMDAN-MPE算法是對原始信號(hào)進(jìn)行CEEMDAN分解得到IMF分量，對各IMF分量進(jìn)行多尺度排列熵的隨機(jī)性檢測，計(jì)算得到各分量信號(hào)的MPE平均值。當(dāng)MPE平均值大于設(shè)定的熵值時(shí)，則被認(rèn)為是異?；蛟肼暢煞郑瑢⑦@些成分從原始信號(hào)中剔除，達(dá)到降噪的目的。采用該算法對信號(hào)進(jìn)行降噪處理，避免了不必要的集成平均，減小了計(jì)算量和添加白噪聲的重構(gòu)誤差，保證了分解算法的完備性；有效地去除了原始信號(hào)中的噪聲成分，能夠獲得較好的降噪效果。

1.3 降噪效果評價(jià)指標(biāo)

為研討爆破地震波信號(hào)的降噪效果，采用信噪比ξ、降噪后信號(hào)和原始信號(hào)的均方根誤差ε作為評價(jià)指標(biāo)[13]，如下所示：

(1)信噪比ξ

(14)

(2)均方根誤差ε

(15)

ξ反映了原始信號(hào)和噪聲的能量關(guān)系，ξ越大表明降噪后信號(hào)更好地保留了原始信號(hào)含有的信息與特征。ε反映了噪聲的平均能量值，體現(xiàn)了降噪后信號(hào)與原始信號(hào)的相似程度，通常ε越小降噪效果越好。此外，除了采用定量參數(shù)客觀評價(jià)降噪效果，還應(yīng)該分析降噪前后信號(hào)的波形特征，確保特征波形的一致性和明顯噪點(diǎn)已經(jīng)去除干凈。

2 工程應(yīng)用

2.1 工程簡介

福建省平潭及閩江口水資源配置工程是一項(xiàng)跨區(qū)域的重大水利工程，屬于國務(wù)院推進(jìn)建設(shè)的172項(xiàng)節(jié)水供水工程之一。工程第4標(biāo)段(大樟溪～石溪輸水線路)由主洞和多條支洞組成，隧洞累計(jì)長度高達(dá)42078 m。隧洞區(qū)沿線分布的地層巖性主要有流紋巖、凝灰?guī)r、凝灰質(zhì)砂礫巖、凝灰質(zhì)砂巖等，埋深一般在70～180 m，最大埋深520 m。

輸水隧洞采用光面爆破開挖，現(xiàn)場爆破監(jiān)測采用TC-4850型測振儀。選取主洞爆破開挖的一條實(shí)測典型地震波信號(hào)為研究對象(如圖1)，信號(hào)采樣頻率為4000 sps，根據(jù)Nyquist采樣定理，實(shí)測信號(hào)的Nyquist頻率為2000 Hz，采用時(shí)間為1s，共采集4000個(gè)采樣點(diǎn)。

2.2 降噪處理與分析

對現(xiàn)場監(jiān)測的地震波信號(hào)進(jìn)行CEEMDAN分解，分解過程中加入了200組信號(hào)標(biāo)準(zhǔn)差為0.2的高斯白噪聲,分解得到的各IMF分量如圖2所示。

由圖2可知，原始信號(hào)經(jīng)過CEEMDAN分解后共得到12個(gè)IMF分量，IMF1～I(xiàn)MF12的中心頻率逐漸降低，高頻噪聲對IMF分量的影響逐漸減弱，IMF分量所含真實(shí)信號(hào)成分不斷增加。根據(jù)IMF分量的波形和中心頻率變化，可以推斷IMF1～I(xiàn)MF5可能為高頻噪聲分量，IMF6～I(xiàn)MF12則為地震波真實(shí)信息。

為了準(zhǔn)確確定真實(shí)信號(hào)成分和噪聲，利用多尺度排列熵方法計(jì)算各IMF分量的MPE值。計(jì)算過程中，需要選取合適的嵌入維數(shù)m、時(shí)間延遲τ和尺度因子s，經(jīng)過多次試算取m=6，τ=1，s=5。計(jì)算得到各IMF分量的MPE平均值如表1所示。

表1 IMF分量的MPE平均值Table 1 Mean MPE of IMF

由表1可知IMF1～I(xiàn)MF12，MPE平均值是逐漸減小的，說明噪聲成分逐漸減少，說明噪聲對不同IMF分量影響不同，與上述波形分析結(jié)果一致。對于爆破地震波信號(hào)[14]，通常有效信號(hào)成分的MPE閾值為0.6，IMF1～I(xiàn)MF5的MPE平均值大于閾值為噪聲信號(hào)成分，需要將其從原始信號(hào)中除去。因此，得到降噪后爆破地震波信號(hào)如圖3所示。采用AOK時(shí)頻技術(shù)[15]分別對原始信號(hào)和降噪后信號(hào)進(jìn)行處理，得到二者的時(shí)頻譜如圖3～圖4所示，圖中X為峰值能量，Y為主頻。

對比圖3和圖4，與實(shí)測爆破地震波原始信號(hào)相比，降噪后信號(hào)的噪聲成分明顯減少，更好地反映了地震波波形特征。由頻譜圖可知，通過降噪處理去除了信號(hào)的高頻成分，對信號(hào)的主頻沒有影響，峰值能量也僅降低了0.5。說明CEEMDAN-MPE算法不僅能成功地去除高頻噪聲能量，而且對地震波信號(hào)所含主要信息的影響極小。

2.3 降噪效果對比

為驗(yàn)證CEEMDAN-MPE算法的有效性，采用EEMD-MPE算法、CEEMDAN閾值算法對上述實(shí)測地震波信號(hào)降噪處理。原始信號(hào)與降噪后信號(hào)如圖5所示。計(jì)算信噪比ξ、降噪后信號(hào)和原始信號(hào)的均方根誤差ε，如表2所示。

表2 爆破振動(dòng)信號(hào)降噪效果指標(biāo)Table 2 Denoised effect index of blasting vibration signals

由表2可知，CEEMDAN-MPE算法的信噪比ξ為23.49 dB，均大于EEMD-MPE和CEEMDAN算法，表明該算法得到的降噪后信號(hào)更好地保留了原始信號(hào)含有的信息與特征；CEEMDAN-MPE算法的均方根誤差ε最小，說明降噪后信號(hào)與原始信號(hào)有更高的相似度。表明CEEMDAN-MPE算法在處理爆破地震波信號(hào)中具有更好降噪效果。由圖5可以看出，CEEMDAN-MPE算法將爆破地震波所含噪聲成分基本去除干凈，能夠很好的展現(xiàn)其波形特征；而EEMD-MPE降噪后的信號(hào)仍有明顯的噪聲，CEEMDAN降噪后的信號(hào)則在峰值點(diǎn)處有較為明顯噪聲。其主要原因是，EEMD分解過程中加入了白噪聲，以減少了模態(tài)混疊現(xiàn)象，但由于加入的白噪聲無法消除，導(dǎo)致降噪效果不理想；CEEMDAN方法成對加入白噪聲，消除了白噪聲的影響。通過計(jì)算IMF分量的MPE平均值，能夠更好的判斷各分量所含噪聲成分，去除高頻噪聲成分，提高了降噪效果。因此，由波形分析和降噪效果指標(biāo)可知，CEEMDAN-MPE算法的降噪效果優(yōu)于EEMD-MPE和CEEMDAN算法。

3 結(jié)論

由于工程環(huán)境和監(jiān)測設(shè)備的影響，隧道爆破實(shí)測地震波信號(hào)中不可避免的存在大量噪聲，掩蓋了真實(shí)信號(hào)所包含信息，不利于爆破振動(dòng)效應(yīng)分析與控制。針對這一問題提出了CEEMDAN和MPE相結(jié)合的算法用于地震波信號(hào)降噪處理。主要研究結(jié)論如下：

(1)利用CEEMDAN算法對隧道爆破地震波信號(hào)進(jìn)行分解，得到不同頻帶的IMF分量，對各IMF分量進(jìn)行多尺度排列熵的隨機(jī)性檢測，利用MPE熵值去除噪聲IMF分量，達(dá)到降噪的目的。工程應(yīng)用表明該方法達(dá)到較好地去除高頻噪聲的目的。

(2)利用AOK時(shí)頻分析技術(shù)，對比分析降噪前后地震波信號(hào)的時(shí)頻特征，說明通過CEEMDAN-MPE算法的降噪處理去除了信號(hào)的高頻成分，對地震波信號(hào)所含主要信息的影響極小。

(3)將EEMD-MPE、CEEMDAN和CEEMDAN-MPE三種算法的降噪效果進(jìn)行對比分析，三種方法均具有一定的降噪效果。波形分析和降噪效果指標(biāo)均表明CEEMDAN-MPE算法的降噪效果最優(yōu)，驗(yàn)證了基于該方法的有效性，對隧道爆破地震波信號(hào)降噪及分析具有指導(dǎo)意義。