基于EEMD-SCBSS的巖石聲發(fā)射信號去噪方法

趙 奎,楊道學(xué),曾 鵬,王曉軍,鐘 文,龔 囪,閆 雷

(1.江西理工大學(xué) 資源與環(huán)境工程學(xué)院，江西 贛州 341000; 2.江西理工大學(xué) 江西省礦業(yè)工程重點(diǎn)實(shí)驗室，江西 贛州 341000)

巖石等固體材料在受到外部應(yīng)力作用時，其內(nèi)部儲存的應(yīng)變能將以彈性波形式被快速釋放，發(fā)出聲響，稱其為聲發(fā)射[1-3]。聲發(fā)射作為一種無損監(jiān)測技術(shù)，現(xiàn)在已被廣泛地應(yīng)用到巖石力學(xué)與工程領(lǐng)域中，已成為揭示巖石內(nèi)部不穩(wěn)定性變形、損傷及破壞演化規(guī)律的一種重要手段[4-10]。

AE信號中包含著巖石內(nèi)部的損傷演化信息，為了研究巖石變形破壞過程的損傷演化特征，需對包含巖石內(nèi)部損傷信息的AE信號(含噪聲)進(jìn)行濾波處理，進(jìn)而準(zhǔn)確地提取AE信號相關(guān)參數(shù)對巖石變形破壞過程中的損傷信息進(jìn)行有效地評價。目前，常用的信號濾波方法主要有：傅里葉變換去噪、小波閾值、經(jīng)驗?zāi)B(tài)分解、BSS等方法[11]。由于聲發(fā)射信號具有低信噪比、隨機(jī)性強(qiáng)、非平穩(wěn)性等特點(diǎn)，傳統(tǒng)的傅里葉變換濾波方法對周期信號具有較好的濾波效果，對非平穩(wěn)信號濾波效果不佳；小波閾值濾波是目前應(yīng)用較為廣泛的濾波方法，但其濾波效果與小波基函數(shù)類型及閾值的設(shè)置息息相關(guān)[12]，在一定程度上限制了其在AE信號濾波領(lǐng)域中的應(yīng)用；經(jīng)驗?zāi)B(tài)分解(EMD)可以較好地適用于非平穩(wěn)信號去噪，但是EMD在分解過程中容易產(chǎn)生模態(tài)混疊現(xiàn)象[13]，為了解決這個問題，Wu等提出了EEMD濾波方法，在一定程度上抑制了模態(tài)混疊問題[14]，但EEMD分解是提取信號時間序列的局部特征，而且在聲發(fā)射信號濾波方面效果不佳[15-16]； BSS技術(shù)是近十幾年發(fā)展起來的一種重要的信號濾波手段，且在非平穩(wěn)信號濾波中得到廣泛應(yīng)用[17]，但是盲信號處理過程中存在太多未知條件，而且?guī)r石內(nèi)部隨機(jī)分布的微裂紋導(dǎo)致彈性波在傳播過程中發(fā)生不同機(jī)制的衰減，因此，欠定情況也是聲發(fā)射信號中普遍存在的問題，等等這一系列問題都會導(dǎo)致直接使用BSS的濾波效果不佳[18]。

為了提高BSS技術(shù)在AE信號中的濾波效果，本文提出了EEMD及SCBSS相結(jié)合的AE濾波方法，該方法對含噪聲的AE信號進(jìn)行EEMD分解，得到高頻的噪聲信號IMF分量，與觀測信號一起構(gòu)建虛擬多通道觀測信號，再利用FastICA算法對構(gòu)建的虛擬多通道觀測信號進(jìn)行BSS，實(shí)現(xiàn)對含噪聲的AE信號進(jìn)行有效地濾波處理。

1 基于EEMD分解單通道聲發(fā)射信號的盲源分離理論

1.1 盲源分離理論

1.1.1 獨(dú)立分量分析基本理論

BSS理論是近十幾年發(fā)展起來一種重要的信號處理技術(shù)，在信號處理領(lǐng)域中被廣泛地應(yīng)用[19]。其中獨(dú)立分量分析(ICA)是BSS理論中一種重要的統(tǒng)計學(xué)方法[20]，ICA的核心思想就是使用優(yōu)化算法將觀測信號分離成若干個近似的估計源信號，但需滿足源信號之間相互獨(dú)立的這一前提條件。

假設(shè)巖石變形破壞過程中產(chǎn)生的含噪聲AE信號是由n個相互獨(dú)立的源信號(AE信號和噪聲信號)經(jīng)線性混合而成，被m個接收傳感器所接收，那么每個AE接收傳感器接收含噪聲的AE信號就是這n個源信號的一個線性組合信號，第j個聲發(fā)射傳感器接收的觀測信號可表示如下

(1)

式中:si(t)為第i個非高斯分布的源信號；bji為混合參數(shù)(i=1,2,…,n，j=1,2,…,m) ；N為含噪聲AE信號長度。m個含噪聲AE信號組成觀測信號矩陣為X(t)=[x1(t),x2(t),…,xm(t)]，將n個源信號組成源信號矩陣S(t)=[s1(t),s2(t),…,sn(t)]，經(jīng)過ICA分解可以得到n個估計源信號，n個估計源信號組成源信號近似解矩陣Y(t)=[y1(t),y2(t),…,yn(t)]。由式(1)可知觀測信號矩陣X(t)與源信號矩陣S(t)存在如下關(guān)系

X(t)=BS(t)

(2)

式中，B是由混合參數(shù)bji組成的混合矩陣。BSS的核心思想就是在混合矩陣B和源信號矩陣S(t)未知情況下通過優(yōu)化算法找到一個分解矩陣W∈Rn×m，使輸出估計源信號矩陣盡可能的接近于S(t)，即：

Y(t)=WX(t)=WBS(t)=GS(t)

(3)

式中，G∈Rn×n表示全局矩陣。通過優(yōu)化迭代算法計算的全局矩陣G=I(I∈Rn×n)，由于噪聲與源信號相關(guān)性及迭代算法計算精度的影響，計算的估計源信號肯定會存在數(shù)值誤差，所以全局矩陣G的每行每列元素中有且僅有一個元素接近于1，而其它元素都接近于零，就可以認(rèn)為分離效果較好。

1.1.2 觀測信號預(yù)處理

在觀測信號進(jìn)行ICA處理之前，需要對觀測信號進(jìn)行預(yù)處理。預(yù)處理主要包括兩個步驟：中心化處理和白化處理。

中心化又稱去均值，實(shí)際上就是減掉觀測信號均值，使觀測信號矩陣均值為零。中心化是對觀測信號最基本的預(yù)處理，其處理過程如下式所示

X(t)=X(t)-E(X(t))

(4)

觀測信號白化處理的基本原理：是指對經(jīng)過中心化處理后的多維觀測信號進(jìn)行一個線性變換，使其生成各維度之間互不相關(guān)的觀測信號，觀測信號經(jīng)過白化處理可以提高分離算法的穩(wěn)定性及收斂速度。

由式(4)可知：經(jīng)中心化處理后的觀測信號X(t)自相關(guān)函數(shù)矩陣RXX=E[X(t)XT(t)]，其是一個由m個特征值組成的對角矩陣D=diag(d1,d2,…,dm)，m個特征值對應(yīng)的特征向量組成特征矩陣U。則去中心化處理后的觀測信號X(t)的白化矩陣Q如下式所示

Q=D-1/2UT

(5)

中心化處理后的觀測信號X(t)經(jīng)過白化處理可得

(6)

1.1.3 基于負(fù)熵理論構(gòu)建目標(biāo)函數(shù)

負(fù)熵理論是基于信息論對信號隨機(jī)性評價參數(shù)熵提出的一種理論。由熵的定義可知：當(dāng)觀測信號的方差相同時，變量的隨機(jī)性越大，則該觀測信號的熵也就越大。為了使用方便，定義負(fù)熵作為度量標(biāo)準(zhǔn)，觀測信號xi(t)的負(fù)熵表示為[21]

Ji=H(xGauss)-H(xi) (i=1,2,…,m)

(7)

式中：xGauss是高斯變量。負(fù)熵理論就是使負(fù)熵最大化，即觀測信號xi(t)熵遠(yuǎn)離xGauss的微分熵，進(jìn)而實(shí)現(xiàn)從觀測信號中分離出聲發(fā)射源信號。

則基于負(fù)熵理論建立的目標(biāo)函數(shù)如下式所示

(8)

1.1.4 快速不動點(diǎn)優(yōu)化算法

快速不動點(diǎn)優(yōu)化算法又被稱為FastICA算法，是一種基于牛頓迭代法的數(shù)值優(yōu)化算法。其核心思想就是通過數(shù)值計算解決基于負(fù)熵理論構(gòu)建的目標(biāo)函數(shù)非高斯性最大值尋優(yōu)問題。因為負(fù)熵值計算過程需使用理論源信號的分量，因此在實(shí)際應(yīng)用時需要對其進(jìn)行近似處理，即：

JG(W)=(E(G(WTυ))-E(G(υ)))2

(9)

式中，υ為一高斯變量，函數(shù)G表示非二次光滑函數(shù)，常用的G函數(shù)如下表示

(10)

(11)

(12)

FastICA算法的迭代過程。

(13)

(14)

由于篇幅原因，本文簡單地敘述FastICA算法的迭代過程，詳細(xì)地介紹參見文獻(xiàn)[21]。

1.2 EEMD信號分解理論

Wu等提出了一種信號處理新算法-希爾伯特黃變換(HHT)理論，該理論中包含了經(jīng)驗?zāi)B(tài)分解(EMD)。EMD是一種新的信號分解方法，其原理就是依據(jù)信號的局部時間特性，將信號分解成一系列單量信號分量，這種分量被稱為本征模函數(shù)[22]。由于EMD具有良好的時頻分析能力，在非平穩(wěn)聲發(fā)射信號處理等許多領(lǐng)域被廣泛的重視。但隨著國內(nèi)外學(xué)者對EMD深入的研究，發(fā)現(xiàn)EMD算法在分解信號過程中經(jīng)常會出現(xiàn)模態(tài)混疊問題。為了解決這一問題，Zheng等[23]對EMD算法進(jìn)行了一系列的改進(jìn)，提出了EEMD。與EMD不同之處在于在分解觀測信號之前，需向觀測信號中引入一任意白噪聲，就可以有效地抑制EMD分解的本征模函數(shù)中的模態(tài)混疊現(xiàn)象。

EEMD算法具體分解過程如下：

將觀測信號x(t)中加入隨機(jī)正態(tài)分布的白噪聲υ(t)

(15)

(16)

由于篇幅原因，本文只是簡單地敘述了EEMD算法分解過程，EEMD算法詳細(xì)推導(dǎo)過程參見文獻(xiàn)[17]。

1.3 EEMD-SCBSS理論

巖石力學(xué)實(shí)驗室試驗過程中需要采集應(yīng)力、軸向應(yīng)變、側(cè)向應(yīng)變及聲發(fā)射信息，由于巖石試件表面積有限，不能提供足夠的空間同時安裝多個聲發(fā)射探頭，而且?guī)r石并不是各向同性介質(zhì)，聲發(fā)射接收傳感器所在位置的不同，也會導(dǎo)致接收的聲發(fā)射信號也存在著一定的差異性，巖石在塑性階段的聲發(fā)射事件率較高，分離出不同接收傳感器接受的同一個聲發(fā)射事件難度較大。因此，研究SCBSS技術(shù)對實(shí)驗室內(nèi)巖石力學(xué)聲發(fā)射試驗過程中的AE信號濾波具有實(shí)際意義，從方程組的角度出發(fā)SCBSS是一個病態(tài)方程，直接求解這類方程組難度較大，故本文將通過EEMD分解觀測信號，經(jīng)分解得到一系列IMF分量，AE信號降噪過程中，通常把高頻的IMF分量作為噪聲信號，利用噪聲信號與觀測信號構(gòu)建虛擬多通道觀測信號，然后再使用BSS算法對虛擬多通道觀測信號進(jìn)行濾波處理。

因此，基于EEMD-SCBSS算法，具體流程如下：

步驟1 利用EEMD將單通道的AE信號x(t)分解成多個本征模函數(shù)IMFi(t)(i=1,2,…,h)，即：

IMF=[IMF1(t) IMF2(t) … IMFh(t)]

(17)

步驟2 使用選取的本征模函數(shù)高頻信號矩陣IMF1(t)與觀測信號x(t)組成虛擬多通道觀測信號，即：

IMF1=[x(t) IMF1(t)]

(18)

步驟3 使用FastICA算法對虛擬多通道觀測信號IMF1進(jìn)行BSS，進(jìn)而得到AE信號與噪聲信號。

1.4 分離算法性能評價標(biāo)準(zhǔn)

本文還需要對EEMD-SCBSS濾波效果進(jìn)行評價，目前常用的信號濾波效果評價參數(shù)有：信噪比、峭度、相關(guān)系數(shù)、信號標(biāo)準(zhǔn)差等，本文將使用信噪比這一參數(shù)對該算法的濾波效果進(jìn)行評價。

信噪比定義如下

(19)

式中：x(t)為模擬原始AE信號；d(t)為降噪后的AE信號；T為AE信號采樣長度；RSN為信噪比，信噪比是描述含噪聲AE信號濾波效果的量化評價參數(shù)，信噪比越高表明信號中的真實(shí)AE信號的信息量越大，濾波效果越好。

2 仿真分析

為了驗證本文提出的基于EEMD-SCBSS的AE信號濾波方法的可行性及有效性，利用含噪聲模擬AE信號在MATLAB計算平臺上進(jìn)行濾波實(shí)驗，將與小波閾值濾波后的模擬AE信號時域波形及頻譜圖進(jìn)行比較。

(20)

式中:Ai表示第i個疊加信號的振幅；Qi表示第i個疊加信號的衰減因子；fi表示第i疊加信號的主頻；ti表示第i個疊加信號的延遲時間；n表示模擬聲發(fā)射信號中疊加信號的個數(shù)。本文模擬AE信號選用的參數(shù)如下所示：n=2，A1=A1=2，Q1=6.14×108，Q2=1.48×108，t1=4×10-4s，t1=5×10-4s，f1=80 kHz，f2=50 kHz。模擬AE信號采樣率設(shè)置為1MSPS，在巖石力學(xué)AE試驗過程中會受到各種背景噪聲對采集的AE信號產(chǎn)生干擾，為了更真實(shí)地模擬巖石力學(xué)AE試驗過程中采集的AE信號，在模擬AE信號中加入一組白噪聲。為了分析白噪聲對模擬AE信號頻域信息的影響，利用快速傅里葉變換方法對模擬AE信號及含噪聲模擬AE信號的頻譜特征進(jìn)行分析。模擬AE信號的時域波形及頻譜圖與含噪聲模擬AE信號的時域波形及頻譜圖如1所示。

(a) 模擬原始聲發(fā)射信號(b) 模擬原始聲發(fā)射信號頻譜圖

(c) 含噪聲模擬聲發(fā)射信號(d) 含噪聲模擬聲發(fā)射信號頻譜圖圖1 模擬聲發(fā)射信號的時域波形及頻譜圖與含噪聲模擬聲發(fā)射信號的時域波形及頻譜圖Fig.1 Time domain waveform and spectrum diagram of simulated acoustic emission signal and time domainwaveform and spectrum diagram of simulated acoustic emission signal with noise

由圖1可知：含噪聲模擬AE信號的頻譜也出現(xiàn)了雙峰現(xiàn)象，主頻所對的振幅譜值低于模擬原始AE信號的振幅譜值，高頻部分出現(xiàn)了明顯的振幅譜值“振蕩”現(xiàn)象；噪聲影響了模擬AE信號到時、振鈴計數(shù)、上升時間、持續(xù)時間及能量等參數(shù)的提取，對AE參數(shù)在定性分析巖石內(nèi)部微裂紋演化特征中的應(yīng)用產(chǎn)生嚴(yán)重地干擾。

對圖1中的含噪聲模擬AE信號進(jìn)行EEMD自適應(yīng)分解，得到IMF分量以頻率的高低依次排列，由于篇幅原因，在這里只畫出了前4個IMF分量時域波形圖，如圖2所示。

依據(jù)EEMD-SCBSS濾波理論，利用IMF1的高頻噪聲信號與含噪聲模擬AE信號構(gòu)建虛擬多通道觀測信號，利用BSS算法對構(gòu)建的虛擬多通道觀測信號進(jìn)行去噪，進(jìn)而得到噪聲信號及模擬AE信號，其時域波形及頻譜圖如圖3所示。由于巖石內(nèi)部微裂紋非均勻分布及各向異性，在受載荷過程中產(chǎn)生數(shù)以萬計的AE信號也將具有較大的差異性，因此在小波閾值對AE信號濾波時，使用默認(rèn)設(shè)置，以便實(shí)現(xiàn)對AE信號自動化處理。含噪聲模擬聲發(fā)射信號經(jīng)小波閾值濾波后的時域波形及頻譜圖，如圖4所示。應(yīng)用式(18)計算含噪聲模擬AE信號經(jīng)EEMD-SCBSS和小波閾值濾波前后的信噪比，計算結(jié)果如表1所示。

(c) IMF3(d) IMF4圖2 含噪聲模擬聲發(fā)射信號EEMD分解結(jié)果Fig.2 EEMD decomposition results of simulated acoustic emission signals with noise

(a) 盲源分離后的噪聲信號(b) 盲源分離后的噪聲信號頻譜圖

(c) 盲源分離去噪后的模擬聲發(fā)射信號(d) 盲源分離去噪后的模擬聲發(fā)射信號頻譜圖圖3 EEMD-SCBSS濾波后的信號時域波形及其頻譜圖Fig.3 Time domain waveform and spectrum diagram of signal after EEMD-SCBSS denoising

(a) 小波閾值去噪后的模擬聲發(fā)射信號

表1 不同濾波方法的濾波效果對比

由表1可知：模擬原始AE信號中加入白噪聲后的信噪比為1.172 dB，此時模擬AE信號嚴(yán)重失真，小波閾值與EEMD-SCBSS濾波后的模擬AE信號信噪比相差較小，基本維持在5 dB左右。但由圖3與圖4可知：EEMD-SCBSS濾波后的模擬AE信號主頻為50 kHz，次主頻為80 kHz，次主頻特征顯著；小波閾值濾波后的模擬AE信號主頻為50 kHz，但是其次主頻特征不明顯，說明小波閾值濾波干擾了模擬AE信號的次主頻特征；模擬原始AE信號的最大振幅值在400 μs，EEMD-SCBSS濾波后的模擬AE信號最大振幅值在400 μs，而小波閾值濾波后的模擬AE信號最大振幅值在500 μs左右，濾波后的時域波形也發(fā)生了變形失真，影響了利用AE信號上升時間，能量等參數(shù)對巖石內(nèi)部損傷演化規(guī)律研究的準(zhǔn)確性。

3 實(shí)測數(shù)據(jù)分析

3.1 巖石單軸壓縮聲發(fā)射試驗

本文進(jìn)行了單軸壓縮聲發(fā)射試驗，試驗過程中的紅砂巖試件都是取自中國贛州地區(qū)，試驗中的紅砂巖試件都取自于同一塊質(zhì)地均勻的紅砂巖巖體。依照國際巖石力學(xué)學(xué)會(ISRM)試驗規(guī)程在實(shí)驗室內(nèi)制作了Φ50 mm×100 mm標(biāo)準(zhǔn)圓柱狀紅砂巖試件，如圖5(a)所示。AE接收傳感器采用UT1000，AE前置放大器閾值設(shè)置為35 dB，采樣率設(shè)置為1MSPS。

圖5(b)和(c)展示本次試驗使用的試驗設(shè)備，其中包括力學(xué)加載系統(tǒng)、AE采集系統(tǒng)。加載試驗在RMT-150C巖石力學(xué)試驗系統(tǒng)中完成，本次試驗加載方式為位移控制方式，位移加載速率為0.002 mm/s，直至紅砂巖試件發(fā)生破壞。AE采集系統(tǒng)是美國聲學(xué)公司生產(chǎn)開發(fā)的Micro-II Digtial AE System硬件和AEwin軟件。

圖5 試驗試件及設(shè)備Fig.5 Test specimens and equipment

3.2 AE信號去噪

AE信號中包含著巖石內(nèi)部的損傷演化信息，為了研究巖石在變形破壞過程中的損傷演化特征，需對包含巖石內(nèi)部損傷信息的AE 信號進(jìn)行濾波處理，小波閾值與EEMD-SCBSS濾波后的AE信號及頻譜圖如圖6所示。

(a) 含噪聲聲發(fā)射信號(b) 含噪聲聲發(fā)射信號頻譜圖

(c) 盲源分離去噪后的聲發(fā)射信號(d) 盲源分離去噪后的聲發(fā)射信號頻譜圖

(e) 小波閾值去噪后的聲發(fā)射信號(f) 小波閾值去噪后的聲發(fā)射信號頻譜圖圖6 聲發(fā)射信號經(jīng)不同濾波方法處理后的時域波形及頻譜圖Fig.6 Time domain waveform and spectrum diagram of acoustic emission signals processed by different denoising methods

由圖6可知，基于EEMD-SCBSS濾波后的AE信號高頻噪聲基本完全濾除，其頻譜特征與含噪聲AE信號的頻譜特征較為相似，經(jīng)EEMD-SCBSS濾波后的AE信號頻域特征基本保持不變；但是經(jīng)小波閾值濾波后的AE信號頻譜特征完全改變，濾波后的AE信號頻域波形發(fā)生了失真變形。上述試驗結(jié)果表明：本文方法可以高效地對實(shí)際AE信號進(jìn)行濾波處理，能夠較好地抑制AE信號中的非平穩(wěn)隨機(jī)噪聲，并能夠有效地保護(hù)AE信號中的頻譜信息。

4 結(jié) 論

針對AE信號的低信噪比、隨機(jī)性強(qiáng)、非平穩(wěn)性等特點(diǎn)，本文提出了一種基于總體經(jīng)驗?zāi)B(tài)分解及SCBSS濾波算法。通過對含噪聲的模擬聲發(fā)射信號及AE信號頻譜分析可知，小波閾值濾波方法會導(dǎo)致濾波后的AE信號頻域信息失真變形，影響濾波后的AE信號上升時間，能量等參數(shù)識別；EEMD-SCBSS可以對含噪聲AE信號進(jìn)行有效地濾波處理，能夠較好地抑制AE信號中的非平穩(wěn)隨機(jī)噪聲，并且能夠保護(hù)濾波后的AE信號頻域信息。