李 薇，白艷萍

(中北大學 理學院， 太原 030051)

現(xiàn)今，傳感器水聲微弱信號的提取方法非常多，主要有傅里葉濾波法、小波變換法、快速獨立分量分析法(Fast ICA)、自適應濾波法和經(jīng)驗模態(tài)分解法(EMD)等[1-10]。這些方法對含噪傳感信號進行去噪會起到一定效果，但依然存在一些不足。HUANG等[1]提出一種新的時頻分析方法—經(jīng)驗模態(tài)分解法，其不足之處是分解的固有模態(tài)函數(shù)(IMF)會產(chǎn)生模態(tài)混疊的情況，導致在重構(gòu)信號時仍有大量噪聲混入。為了降低模態(tài)混疊的影響，WU等在EMD基礎(chǔ)上提出了集合經(jīng)驗模態(tài)(EEMD)方法，通過加入輔助白噪聲來降低模態(tài)混疊影響。為了消除白噪聲引起的重構(gòu)誤差，Jia Rong Yeh等提出了互補集合經(jīng)驗模態(tài)分解(CEEMD)方法，但是CEEMD方法中如果添加的白噪聲和迭代次數(shù)不合適時，會出現(xiàn)一些IMF偽分量。針對此問題研究者提出了MEEMD算法。由于單一的模態(tài)分解方法在去噪方面存在缺陷，有學者進行了模態(tài)分解方法與小波閾值去噪法的結(jié)合，降低了模態(tài)分解方法的混疊效應，比單一模態(tài)分解方法更優(yōu)越。本文在上述研究基礎(chǔ)上提出了MEEMD分解和小波軟閾值的聯(lián)合去噪方法。通過對比EEMD、CEEMD和小波軟閾值[2]的聯(lián)合去噪方法，發(fā)現(xiàn)本文算法更具優(yōu)勢。

1 幾種模態(tài)分解算法

1.1 EEMD算法

算法步驟：① 原信號中加入不一樣的高斯白噪聲；② 對目標信號進行模態(tài)分解；③ 重復上面兩步驟；④ 對分解結(jié)果取平均值，消除多次加入高斯白噪聲對固有模態(tài)的影響。

1.2 CEEMD算法

CEEMD算法在原始信號的基礎(chǔ)上，加入了n組正負成對的輔助噪聲，獲得兩套IMF集合：

其中：S表示原始信號；N表示輔助噪聲；M1、M2分別表示加入正、負成對噪聲后的信號。

算法步驟同EEMD算法：最后，得到2n個集合信號。

1.3 MEEMD算法

EEMD、CEEMD算法限制了迭代次數(shù)，會使得分解得到的IMF分量不滿足定義，所以當有異常分量出現(xiàn)后，沒必要對加入噪聲進行EMD分解。針對此，提出了對高頻和間歇信號的檢測方法——基于排列熵的信號隨機性檢測。

1.3.1 排列熵定義

排列熵(PE)[4]是一種檢測隨機數(shù)列隨機性和動力學突變的方法，抗干擾能力強，計算快，適合于非線性數(shù)據(jù)。

對長度為N的時間序列{x(i),i=1,2,…,N}進行相空間重構(gòu)，得到下面序列：

(1)

其中：m是嵌入維數(shù)；λ是時間延遲。將x(i)的m個向量X(i)={x(i),x(i+λ),…,x(i+(m-1)λ)}按照升序重新排列得：

X(i)={x(i+(j1-1)λ)≤

x(i+(j2-1)λ)≤…,x(i+(jm-1)λ)}

(2)

若存在

X(i+(ji1-1)λ)=X(i+(ji2-1)λ)

則按j值的大小進行排序，任意一個X(i)都可以得到一組序列：

S(g)={j1,j2, …,jm}

(3)

時間序列{x(i)=1,2,3,…,N}的排列熵可以按照Shannon熵的形式定義為

(4)

(5)

顯然Hp的取值范圍是0≤Hp≤1，值的大小表示隨機數(shù)列的隨機性。

1.3.2PE參數(shù)的選取

計算排列熵需要確定3個參數(shù)：時間序列長度N、切入維數(shù)m、時間延遲λ。有學者建議嵌維數(shù)取3～7，一般取m=6。當嵌入維數(shù)較小時，要求數(shù)據(jù)長度也較??；若m=6,則數(shù)據(jù)長度需大于1 024。時間延遲為λ，本文取λ=1。選取仿真信號為S(t)=0.1cos(2π·500t)，含有高斯白噪聲、隨機噪聲和脈沖噪聲。在選定上述參數(shù)的情況下，4種成分的PE值分別為0.373 9、0.969 8、0.971 4、0.614 6，可以發(fā)現(xiàn)高斯白噪聲和隨機信號的PE值較大，較為隨機，余弦信號PE值較小，不隨機。脈沖信號PE值大于0.6，相比余弦信號來說也是較為隨機，因此基于排列熵的隨機性檢測可以用于異常信號檢測。

2 小波軟閾值去噪

小波閾值去噪的基本原理是設(shè)置一個臨界閾值λ。若小波系數(shù)小于λ，則系數(shù)主要由噪聲產(chǎn)生，將這部分系數(shù)去掉；若小波系數(shù)大于λ，則系數(shù)主要由信號產(chǎn)生，則留下這部分系數(shù)。然后，利用小波反變換對小波系數(shù)進行處理得到去噪后的信號[1]。小波去噪基本原理見圖1。

圖1 小波去噪流程

在進行小波閾值去噪時首先選擇小波基，然后確定合適的分解層數(shù)以及閾值，最后選擇合適的閾值函數(shù)。

小波閾值去噪方法的步驟：① 將原始信號x(t)進行小波變換到小波域，于是得到一組分解系數(shù)；② 在小波域進行閾值處理后得到包含大量隨機噪聲的較小的小波系數(shù)；③ 利用處理后得到的小波系數(shù)進行信號重構(gòu)，得到去噪后的信號。軟閾值函數(shù)為

其中：sgn(·)為符號函數(shù)；λ為閾值。

3 基于EEMD、CEEMD、MEEMD的小波軟閾值

將原始信號分別進行EEMD、CEEMD、MEEMD分解，得到一組固有模態(tài)分量。由于前幾層IMF中仍含有少量細節(jié)信號，所以采用小波閾值函數(shù)對前幾層IMF提取細節(jié)信息，得到新的分量后與剩余分量進行重構(gòu)。

假設(shè)噪聲信號為

y(t)=x(t)+η(t)

其中：x(t) 為無噪聲信號；η(t) 為有限振幅的獨立噪聲。噪聲信號y(t)首先被分解為固有模態(tài)IMF分量，選用小波基為db7，分解層數(shù)為6。選擇含噪高的分量采用軟閾值進行系數(shù)提取，該方法有一個可調(diào)用函數(shù)，選用全局閾值。在EEMD、CEEMD、MEEMD分解完成后，選取合適的模態(tài)應用小波軟閾值函數(shù)進行信號提取，提取后進行重構(gòu)。小波軟閾值函數(shù)如下：

去噪信號的重構(gòu)如下：

4 仿真和實測

4.1 仿真

對所提出的3種算法進行仿真比較。根據(jù)中北大學國防重點實驗室在汾河進行的汾機實測數(shù)據(jù)分析，發(fā)現(xiàn)該信號為一單頻正弦信號序列。信號在發(fā)射和傳輸過程中，不僅受機器本身所產(chǎn)生干擾和漂移影響，還受傳播過程中環(huán)境噪聲的影響。因此，本文選用的仿真實驗信號為s(t)=0.1cos(2π·500t)，振幅值為0.1，頻率為500 Hz。實驗所用軟件為Matlab(2014a)。向該信號中加入隨機噪聲、脈沖噪聲，使得仿真實驗更逼真。選用信噪比和均方差作為性能指標，其計算公式如下：

選取信噪比為7.75 dB的含躁信號進行去噪，圖2～4分別是含噪信號經(jīng)過EEMD、CEEMD、MEEMD方法分解后得到的各層固有模態(tài)。從圖2、3可以看出：原始信號主要集中在imf2、imf3上，圖4的信號則主要集中在imf1、imf2上。

對含信號較多的層分別采用小波軟閾值處理進行細節(jié)信號的提取，得到新的分量后與剩余未處理過的分量進行重構(gòu)。得到3種方法的去噪效果后，取其中200～400個點放大觀察，得到圖5～7。

從圖5～7可以看出：本文提出的方法效果最好。仿真中選取5組信噪比不同的含噪信號進行去噪，分別得到基于EEMD、CEEMD、MEEMD方法的軟閾值仿真去噪的性能指標。表1是去噪性能指標對比。

圖2 含躁信號經(jīng)過EEMD方法分解后得到的各層固有模態(tài)

圖3 含噪信號經(jīng)過CEEMD方法分解后得到的各層固有模態(tài)

圖4 含噪信號經(jīng)過MEEMD方法分解后得到的各層固有模態(tài)

圖5 EEMD方法去噪效果采樣

圖6 CEEMD方法去噪效果采樣

圖7 MEEMD方法去噪效果采樣

指標去噪前EEMD小波軟閾值去噪CEEMD小波軟閾值去噪MEEMD小波軟閾值去噪SNR1MSE1 -4.481.42×10-26.281.18×10-37.109.75×10-411.963.18×10-4SNR2 MSE21.323.69×10-311.473.56×10-412.722.67×10-418.147.68×10-5SNR3 MSE34.261.88×10-314.101.95×10-415.191.51×10-420.514.44×10-5SNR4MSE47.758.39×10-416.291.17×10-417.419.07×10-522.342.92×10-5SNR5 MSE515.131.54×10-423.062.47×10-524.871.63×10-530.334.64×10-6

通過對圖5～7和表1的性能指標對比發(fā)現(xiàn)：本文提出的基于MEEMD的小波軟閾值去噪方法的效果要優(yōu)于基于CEEMD的小波軟閾值[1]去噪方法和基于EEMD的小波軟閾值去噪方法。表明本文中提出的去噪方法可在提取水聽器水聲微弱信號去噪中使用。

4.2 實測

本次實測為中北大學國防重點實驗室研究人員在汾河二庫進行的MEMS矢量水聽器的湖試實驗。實驗中采用二元MEMS矢量水聽器線陣，將其固定于船舷一側(cè)，陣元之間相距0.5 m，置于水下10 m，基陣上有羅經(jīng)實時檢測基陣姿態(tài)，保持基陣水平，每個陣元輸出聲壓和兩路振速信號[2]。發(fā)射換能器被放置于基陣的90°方位，分別發(fā)射186、270、331、500、800、1 000、1 500 Hz等連續(xù)單頻信號，采樣頻率為10 kHz，采集信號時1號水聽器為3、4路信號，2號水聽器為1、2路信號。選取186、331、500 Hz這3種實測數(shù)據(jù)包來觀察實測情況，不同頻率實測數(shù)據(jù)如圖8所示。

實測數(shù)據(jù)選取fenji 331 Hz數(shù)據(jù)包的4路信號(見圖8)。對331 Hz的4路陣元信號截取50 001～51 000的1 000個點，得出去噪前的4路實測信號，如圖9所示。圖9的頻譜顯示，信號中摻雜高頻噪聲和低頻的一些干擾，所以在仿真時加入相似的噪聲?；贛EEMD改進小波軟閾值的去噪結(jié)果如圖10所示，從圖中可以看出：該算法在去噪方面有較好的效果，但其幅值和實測信號有出入。這是由于原信號在低頻處含有噪聲，造成實測信號幅值偏移，而本文算法在去噪的過程中較好地還原了原信號。

圖8 不同頻率實測數(shù)據(jù)

圖9 去噪前4路實測信號數(shù)據(jù)

圖10 基于MEEMD改進小波軟閾值的去噪數(shù)據(jù)

5 結(jié)束語

針對MEMS水聽器在接收信號時混入噪聲的情況，提出了一種基于MEEMD的小波閾值去噪方法。在仿真實驗中通過對比EEMD、CEEMD的小波軟閾值去噪方法的性能指標，發(fā)現(xiàn)基于MEEMD的小波軟閾值去噪方法要優(yōu)于本文的其余方法，因此可用于對實測信號數(shù)據(jù)進行去噪。實測結(jié)果表明：新的去噪方法效果明顯，且原信號未出現(xiàn)失真，說明本文提出的方法具有一定的參考價值。

[1] HUANG N E,SHEN Z,LONG S R,et al.The Empirical Mode Decomposition and the Hilbert Spectrum for Nonlinear and Non-stationary Time Series Analysis[J].Proc.R.Soc.Lond.A[J].1998,454:903-995.

[2] 赫彬,白艷萍,張雅婷.基于不同經(jīng)驗模態(tài)分解的小波軟閾值ECG去噪研究[J].數(shù)學的實踐與認識,2016,46(6):136-144.

[3] 王鵬.基于MEMS矢量水聽器陣列的聲目標定向定位技術(shù)研究[D].太原：中北大學，2013.

[4] 楊向林，嚴紅，許志.基于Hilbert-huang變換的ECG消噪[J].電子學報,2011,39(4):819-824.

[5] 鄭近德，程軍圣，楊宇.改進的EEMD算法及其應用研究[J].振動與沖擊，2013,32(21):21-26.

[6] 赫彬，白艷萍.基于CEEMDAN小波包的MEMS水聽器信號去噪[J].數(shù)學的實踐與認識.2016,46(13):139-147.

[7] LU Jingyi,LIN Hong,YE Dong,et al.A New Wavelet Threshold Function and Denoising Application[J].Mathematication Problems in Engineering,2016.

[8] KOPSINIS Y,MCLAUGHLIN S.Emprical Mode Decomposition Based Soft-Thresholding[C]//In Proc.16thEur.Signal Process.Conf.(EUSIPCO).Switzerland:[s.n.],2008.

[9] DONO D L.Denoising by soft threshold[J].TEEE Transactions on Information Theory,1995,41(3):613-627.

[10] CAO Junhong,WEI Zhuobin.Independent component analysis in frequency domain and its application in stractural vibrition signal seperation[J].Procedia Enginering,2011,16(3):511-517.