楊志飛，崔春生，杜桂云，劉雙峰，趙海霞

(1 中北大學(xué)省部共建動態(tài)測試技術(shù)國家重點(diǎn)實(shí)驗(yàn)室，山西 太原 030051；2 中北大學(xué)電氣與控制工程學(xué)院，山西 太原 030051；3 中北大學(xué)環(huán)境與安全工程學(xué)院，山西 太原 030051)

0 引言

爆炸沖擊波信號在軍事、航天等領(lǐng)域具有重要意義[1-3]。在采集爆炸沖擊波超壓信號的過程中,由于高溫環(huán)境、磁場干擾、測試設(shè)備誤差等多種因素的影響,實(shí)測沖擊波信號中往往混有大量環(huán)境噪聲和系統(tǒng)自身的高頻噪聲。為了提高信號處理的準(zhǔn)確性,人們對爆炸沖擊波信號的降噪方法進(jìn)行了廣泛研究[4-18]。

張衍芳等[6]針對不同類型的典型實(shí)測沖擊波信號進(jìn)行了濾波處理,發(fā)現(xiàn)傅里葉變換等傳統(tǒng)濾波方法無法完全消除信號中的毛刺。劉嘉慧[7]對實(shí)測數(shù)據(jù)曲線進(jìn)行了貝塞爾濾波與巴特沃斯濾波處理,結(jié)果表明,截止頻率為40 kHz的6階貝塞爾濾波器在沖擊波信號中去噪效果較好。由于傳統(tǒng)濾波方法在處理非平穩(wěn)信號時容易造成信號失真現(xiàn)象,如傅里葉變換只能在頻域內(nèi)進(jìn)行分析,若信號在時域上某處發(fā)生突變,則無法分辨信號的尖峰是突變還是噪聲導(dǎo)致,限制了其在爆炸沖擊波信號降噪中的應(yīng)用[8]。近年來,基于經(jīng)驗(yàn)?zāi)B(tài)分解(EMD)的信號處理方法在信號降噪領(lǐng)域得到了廣泛關(guān)注。自適應(yīng)白噪聲完整集合經(jīng)驗(yàn)?zāi)B(tài)分解(CEEMDAN)是經(jīng)驗(yàn)?zāi)B(tài)分解(EMD)的一種改進(jìn)方法,通過引入白噪聲信號,增強(qiáng)了模態(tài)分解的魯棒性,能夠更準(zhǔn)確地分解非線性、非平穩(wěn)信號,但EMD在去噪過程中分解出的固有模態(tài)函數(shù)分量之間出現(xiàn)模態(tài)混疊現(xiàn)象[15],通過經(jīng)驗(yàn)?zāi)B(tài)分解去除高頻模態(tài)分量的方式實(shí)現(xiàn)去噪,然而這種“一刀切”式的降噪方式會在降噪的同時消除信號中的有用成分,造成信號失真,對去噪效果會產(chǎn)生影響。梁晶等[16]采用CEEMDAN-小波閾值聯(lián)合降噪算法對含噪模型和實(shí)測數(shù)據(jù)進(jìn)行處理,發(fā)現(xiàn)該方法能夠有效去除沖擊波信號中的高頻噪聲,提高信號質(zhì)量。張冉等[9]采用CEEMDAN-SG聯(lián)合降噪算法解決了其他傳統(tǒng)算法在信號重構(gòu)后局部信號與原始信號差異過大的問題。張健等[10]采用CEEMDAN-PCA-CEEMDAN算法,利用PCA對CEEMDAN算法分解后的IMF分量正交化,以抑制模態(tài)混疊現(xiàn)象。

盡管上述方法在一定程度上提高了爆炸沖擊波信號的降噪效果,但在信號失真、模態(tài)混疊等方面仍存在一定問題。因此,文中提出了一種基于自適應(yīng)白噪聲完整集合經(jīng)驗(yàn)?zāi)B(tài)分解(CEEMDAN)-高斯混合模型(GMM)和小波閾值方法的爆炸沖擊波降噪算法。首先,通過CEEMDAN對原始爆炸沖擊波信號進(jìn)行分解,得到多個本征模態(tài)分量(IMF)。然后,采用GMM聚類算法對IMF進(jìn)行分類,并利用赤池信息準(zhǔn)則(AIC)和輪廓系數(shù)(SC)確定最佳分類類別。剔除高頻噪聲類別的IMF分量后,對受噪聲干擾的IMF分量采用小波閾值方法進(jìn)行濾波去除噪聲。最后,將濾波后的IMF分量與未受噪聲干擾的IMF分量一起重構(gòu),從而獲得降噪后的爆炸沖擊波信號。

與文獻(xiàn)[13,15]相比,最大的不同便是文中在獲取了模態(tài)分解之后采用了GMM聚類算法對其進(jìn)行分類并對此進(jìn)行小波閾值濾波、重構(gòu)等。

1 信號降噪算法實(shí)現(xiàn)

1.1 CEEMDAN-GMM

CEEMDAN-GMM算法是GMM算法[11]針對CEEMDAN降噪時容易造成信號失真所改進(jìn)的算法,雖然CEEMDAN方法在信號降噪過程中減輕了EMD方法的模態(tài)混疊問題和避免了殘余噪聲影響,但該方法容易造成部分有效信號信息缺失,而CEEMDAN-GMM降噪算法充分利用了CEEMDAN在處理非平穩(wěn)信號方面的優(yōu)勢,結(jié)合GMM聚類方法在數(shù)據(jù)分析領(lǐng)域的強(qiáng)大能力,而且能充分發(fā)揮小波閾值去噪法在不同尺度上分別處理信號中的噪聲,從而更好地保留信號中有用信息的能力,實(shí)現(xiàn)對爆炸沖擊波信號的高效降噪。

CEEMDAN-GMM算法流程如圖1所示。其主要步驟為:

圖1 CEEMDAN-GMM算法流程圖Fig.1 CEEMDAN-GMM algorithm flowchart

1)在獲取的原始爆炸沖擊波信號x(t)中加入噪聲序列,得到新信號x′(t),采用CEEMDAN算法分解新信號x′(t)獲得m個不同的分量,每個IMF分量包含n個樣本點(diǎn)IMFi(t);

2)利用高斯混合模型(GMM)[12]對IMF進(jìn)行分類。首先,初始化高斯混合模型的參數(shù),包括每個高斯分布的均值μk、協(xié)方差σk和權(quán)重πk。

3)使用EM算法[13]迭代更新GMM的參數(shù),包括E-step和M-step兩步。其中,E-step用于計算每個樣本屬于各高斯分布的后驗(yàn)概率γ(Zik);M-step根據(jù)后驗(yàn)概率更新高斯分布的均值μk、協(xié)方差σk和權(quán)重πk。其數(shù)學(xué)表達(dá)式分別為:

(1)

式中:i=1,2,…,n;k=1,2,…,m;N(xi|μk,σk)為高斯分布的概率密度函數(shù)。

(2)

(3)

(4)

4)反復(fù)執(zhí)行步驟3),當(dāng)滿足收斂條件時停止迭代,完成聚類過程。根據(jù)后驗(yàn)概率將不同IMF分量劃分至不同類別。

對于每個IMF分量IMFi(t),計算其屬于每個高斯分布k的后驗(yàn)概率γ(Zik);

找出具有最大后驗(yàn)概率的高斯分布k*,將IMF分量IMFi(t)劃分至對應(yīng)的類別k*。

k*=arg maxkγ(zik)

(5)

式中:arg max表示使后驗(yàn)概率γ(Zik)取最大值的高斯分布。通過將每個IMF分量劃分至具有最大后驗(yàn)概率的高斯分布,可以將不同IMF分量劃分至不同類別。

5)計算不同IMF分量的方差貢獻(xiàn)率校核聚類結(jié)果,剔除噪聲類別的IMF分量,將剩余IMF分量重構(gòu)獲得降噪后的爆炸沖擊波信號。

1.2 降噪效果評價指標(biāo)

文中以信噪比RSNR和均方根誤差MRMSE[15]作為評價指標(biāo)來衡量信號降噪效果。信噪比用于衡量信號與噪聲的比值,均方根誤差用于衡量過濾后的信號與原始信號的緊密程度。其計算公式分別為:

(6)

(7)

式中:RSNR為信噪比;MRMSE為均方根誤差;X為原始信號;Y為降噪后的信號;n為信號的序列長度。

一般情況下,降噪后信號RSNR越大,MRMSE越小,降噪效果越好。但降噪效果不僅要考慮評價指標(biāo),還要考慮降噪算法對信號超壓峰值、正壓時間等信號特征的影響較小。

1.3 相關(guān)算法信息

1.3.1 赤池信息AIC

為了在降噪過程中更好地區(qū)分信號與噪聲,從而提高降噪算法的性能,文中引入赤池信息準(zhǔn)則AIC和輪廓系數(shù)SSC作為評價指標(biāo)來確定最佳聚類類別,從而有效地剔除噪聲類別的IMF分量,保留信號中的有用信息,從而達(dá)到較好的降噪效果。

AIC[18]是評估統(tǒng)計模型的復(fù)雜度和衡量統(tǒng)計模型“擬合”資料之優(yōu)良性的一種標(biāo)準(zhǔn)。

在一般情況下,AIC量化值NAIC可以表示為:

NAIC=2k-2lnL

(8)

式中:k是參數(shù)的數(shù)量,L是似然函數(shù)。NAIC鼓勵數(shù)據(jù)擬合的優(yōu)良性,但是盡量避免出現(xiàn)過度擬合的情況。因此優(yōu)先考慮的模型應(yīng)是NAIC值最小的那一個。

1.3.2 輪廓系數(shù)SSC

輪廓系數(shù)SSC[17]使用數(shù)據(jù)集中對象之間的相似性度量來評估聚類的質(zhì)量,是簇的密集與分散程度的評價指標(biāo)。其計算公式為:

(9)

式中:a為該數(shù)據(jù)幀與簇內(nèi)其他數(shù)據(jù)幀的平均距離;b為該數(shù)據(jù)幀與距離它最近的另一個簇中樣本的平均距離。

輪廓系數(shù)SSC的取值范圍為[-1,1],其值越大,聚類結(jié)果越合理。如果輪廓系數(shù)SSC=-1,則表明該數(shù)據(jù)幀應(yīng)該被分到其他的類;如果SSC接近0,則表明該數(shù)據(jù)幀在兩個類的交叉點(diǎn)上。對全部樣本的輪廓系數(shù)求平均值,可以得到聚類結(jié)果整體輪廓系數(shù)STCS,即:

(10)

STSC的取值范圍為[-1,1],STSC值越大,聚類效果越好。

1.3.3 方差貢獻(xiàn)率η

為保證GMM算法聚類分析結(jié)果的準(zhǔn)確性,采用方差貢獻(xiàn)率η[14]進(jìn)行校核,方差貢獻(xiàn)率計算公式為:

(11)

由式(11),根據(jù)方差貢獻(xiàn)率η來判斷類別數(shù)據(jù)受高頻噪聲干擾的高低,方差貢獻(xiàn)率η越低說明此類別受到的噪聲越大,反之為有效數(shù)據(jù)。

2 實(shí)測數(shù)據(jù)處理及分析

2.1 實(shí)驗(yàn)數(shù)據(jù)與預(yù)處理

實(shí)測爆炸沖擊波信號來自60 kg級TNT爆炸沖擊波試驗(yàn),實(shí)測爆炸沖擊波信號曲線如圖2所示。采用CEEMDAN算法分解原始爆炸沖擊波信號,分解后的曲線如圖3所示。其中,高斯白噪聲標(biāo)準(zhǔn)差與原始信號標(biāo)準(zhǔn)差之比為0.2,加入噪聲的次數(shù)為100,最大迭代次數(shù)為800,獲得17個不同頻率的IMF分量和1個殘余分量。根據(jù)圖3可知,噪聲對不同分量的影響逐漸降低,其中,IMF1～I(xiàn)MF8分量具備噪聲分量的特征,IMF9～I(xiàn)MF10分量含有部分信號特征,但是也受到了噪聲干擾,IMF11之后的分量具有和原始信號相似的波形。

圖2 實(shí)測爆炸沖擊波信號Fig.2 Measured explosion shock wave signal

圖3 實(shí)測爆炸沖擊波信號IMF分量Fig.3 IMF component of measured explosion shock wave signal

不同聚類類別的NAIC變化如圖4所示,隨著聚類類別個數(shù)不斷增加,NAIC逐漸減小。由圖4可知,NAIC曲線的轉(zhuǎn)折聚類類別為5,在聚類類別為5之后,NAIC曲線下降曲線逐漸平緩,NAIC可以對模型復(fù)雜性和擬合優(yōu)度進(jìn)行權(quán)衡,以防止模型過于復(fù)雜而導(dǎo)致的過擬合。不同聚類類別的SSC變化如圖5所示,SSC衡量了類別間的分離度和類別內(nèi)部的緊密度。SSC值范圍為-1～1,值越接近1,表示類別間的分離度越高,類別內(nèi)部的緊密度越好。通過比較不同類別數(shù)量下的SSC值,可以找到一個最佳的類別數(shù)量,使得類別間的分離度和類別內(nèi)部的緊密度都達(dá)到最優(yōu)。同樣地,當(dāng)聚類類別為5時,SSC值最接近1。因此為了獲得更好的聚類效果,將聚類類別設(shè)置為5。在CEEMDAN-GMM算法中,使用赤池信息準(zhǔn)則(AIC)和輪廓系數(shù)(SSC)確定最佳聚類類別的意義在于:為高斯混合模型(GMM)提供一個合適的類別數(shù)量,以便在降噪過程中更好地區(qū)分信號與噪聲。這兩個評價指標(biāo)有助于評估聚類結(jié)果的質(zhì)量,并找到最佳的聚類劃分。

圖4 不同聚類類別AIC變化Fig.4 Changes in AIC of different cluster categories

圖5 不同聚類類別SC變化Fig.5 Changes in SC of different cluster categories

根據(jù)GMM算法聚類分析結(jié)果,IMF1～I(xiàn)MF6分量劃分為類別一;IMF7～I(xiàn)MF10分量劃分為類別二;IMF11～I(xiàn)MF12分量劃分為類別三;IMF13、IMF16、IMF17劃分為類別四;IMF14、IMF15、IMF18屬于類別五。類別一重構(gòu)信號如圖6所示。與原始爆炸沖擊波信號圖2相比,類別一IMF分量重構(gòu)波形相似度較差且具有噪聲信息特征,其他類別IMF分量重構(gòu)波形相似度較好且包含較多原始信號信息。因此,可將類別一IMF分量重構(gòu)波形初步定義為高頻噪聲序列,其他類別IMF分量重構(gòu)波形定義為有效信息序列。

圖6 類別一重構(gòu)信號Fig.6 Reconstructed signal of category 1

根據(jù)式(11)計算不同類別的方差貢獻(xiàn)率如表1。

表1 不同類別的方差貢獻(xiàn)率Table 1 Variance contribution rates of different categories

由表1可知,類別五的方差貢獻(xiàn)率最高,其值為75.89%,類別一的方差貢獻(xiàn)率最低,僅有0.46%,說明方差貢獻(xiàn)率的分析結(jié)果和GMM算法聚類的分析結(jié)論相同。因此,可以認(rèn)定類別一位高頻噪聲序列,類別二～類別四是受噪聲干擾的序列,需要降噪處理,類別五為有效信息序列。類別二～類別四經(jīng)小波閾值去噪重構(gòu)后的波形如圖7所示。

圖7 類別二～類別四經(jīng)小波閾值去噪重構(gòu)后的波形Fig.7 Waveforms reconstructed by wavelet threshold denoising for categories No.2 to No.4

將類別一高頻噪聲去除,重構(gòu)類別二～類別四小波閾值降噪信號以及類別五IMF分量得到降噪后的純凈爆炸沖擊波信號,同時與CEEMDAN降噪算法、CEEMDAN-小波閾值降噪算法相比的時域和頻域圖形如圖8所示。

圖8 爆炸沖擊波實(shí)測信號降噪方法對比圖Fig.8 Comparison of noise reduction methods for measured signals of explosion shock waves

2.2 降噪效果評價結(jié)果

根據(jù)式(6)～式(7)可以計算方法處理后各模型的RSNR指標(biāo)和MRMSE指標(biāo)如表2所示。

表2 降噪效果評價指標(biāo)Table 2 Noise reduction effect evaluation indicators

由表2可知:三種方法均具有良好的去噪效果,但是,CEEMDAN-GMM和小波閾值算法的信號降噪效果最好,分別比CEEMDAN算法和CEEMDAN-小波閾值算法分別提高了10.19 dB和2.99 dB。由表3可知,該方法對毀傷評估指標(biāo)的影響最小,相比于其他兩種方法分別降低了0.42×10-4和0.14×10-4,降噪后更好地保留了原始信號信息,與原始信號具有較高的相似度。

表3 不同降噪方法對毀傷評估參數(shù)的影響Table 3 Impact of different noise reduction methods on damage assessment parameters

3 結(jié)論

文中所提出的融合CEEMDAN-GMM和小波閾值方法的爆炸沖擊波降噪算法充分利用了CEEMDAN在非平穩(wěn)信號處理中的優(yōu)勢,并結(jié)合了GMM聚類方法在數(shù)據(jù)挖掘領(lǐng)域的數(shù)據(jù)分析能力。與其他降噪方法相比,文中提出的CEEMDAN-GMM與小波閾值方法相結(jié)合的爆炸沖擊波信號降噪算法在去除高頻噪聲方面具有較好的效果。通過與實(shí)測數(shù)據(jù)對比,CEEMDAN-GMM與小波閾值算法SNR值和RMSE值分別比CEEMDAN算法和CEEMDAN-小波閾值算法提高了10.19 dB和2.99 dB,降低了0.42×10-4和0.14×10-4。同時該算法能夠在降低信號失真程度的同時保留信號中的關(guān)鍵特征。尚需完善的地方有:

1)文中采用的實(shí)測爆炸沖擊波數(shù)據(jù)集數(shù)量有限,可能無法充分反映所有可能的噪聲情況和信號特性。因此,需要擴(kuò)大數(shù)據(jù)集,包括更多類型和來源的爆炸沖擊波信號,以提高算法的通用性和穩(wěn)健性。

2)文中固定了GMM聚類的參數(shù)和小波閾值方法的閾值參數(shù),然而,這些參數(shù)可能會對降噪效果產(chǎn)生顯著影響。需通過參數(shù)優(yōu)選和策略調(diào)整來提高降噪性能。