趙海峰 張 亞 李世中 郭 燕
1.中北大學(xué),太原,030051 2.南京信息職業(yè)技術(shù)學(xué)院,南京,2100233.渥太華大學(xué),渥太華,K1N 6N5
侵徹彈體頻率特性分析及過(guò)載信號(hào)處理
趙海峰1,2,3張亞1李世中1郭燕2
1.中北大學(xué),太原,0300512.南京信息職業(yè)技術(shù)學(xué)院,南京,2100233.渥太華大學(xué),渥太華,K1N 6N5
為解決硬目標(biāo)侵徹過(guò)載信號(hào)降噪問(wèn)題,提出融合總體經(jīng)驗(yàn)?zāi)B(tài)分解(EEMD)和小波變換(WT)的聯(lián)合濾波方法。首先對(duì)實(shí)測(cè)信號(hào)進(jìn)行總體經(jīng)驗(yàn)?zāi)B(tài)分解,獲得信號(hào)的本征模態(tài)函數(shù)(IMF)分量,然后計(jì)算各分量功率譜并與原信號(hào)比較,得出信號(hào)的有效分解尺度和彈體的過(guò)載響應(yīng)頻率,接著對(duì)高頻IMF分量采用小波閾值降噪,最后將降噪后的高頻分量與分解后的低頻分量組合重構(gòu)獲得侵徹特征信號(hào)。實(shí)驗(yàn)證明,這一方法可以有效提取彈體響應(yīng)頻率,消除侵徹過(guò)程中彈體的高頻振動(dòng)信號(hào)和外部噪聲,且處理后的加速度曲線具有更高的信噪比,積分所得速度和位移時(shí)程曲線也與實(shí)驗(yàn)結(jié)果相近。
侵徹過(guò)載;總體經(jīng)驗(yàn)?zāi)B(tài)分解;小波變換;本征模態(tài)函數(shù);功率譜
彈丸侵徹硬目標(biāo)的過(guò)程是一個(gè)非常復(fù)雜的動(dòng)力學(xué)問(wèn)題,其信號(hào)包含軸向加速度、彈體振動(dòng)信號(hào)以及外部噪聲信號(hào),屬于典型的非平穩(wěn)隨機(jī)振動(dòng)過(guò)程。通過(guò)對(duì)實(shí)測(cè)過(guò)載信號(hào)進(jìn)行分析處理,可以得出彈體侵徹過(guò)程的諸多重要參數(shù)(如位移、速度、時(shí)間等),從而為彈體的結(jié)構(gòu)強(qiáng)度設(shè)計(jì)、侵徹參數(shù)預(yù)測(cè)提供重要依據(jù)。彈體在侵徹硬目標(biāo)時(shí),所測(cè)過(guò)載信號(hào)主要包含三種信號(hào)成分:一是彈體在侵徹目標(biāo)介質(zhì)過(guò)程中遇到的侵徹阻力所形成的減加速度信號(hào);二是彈體在侵徹過(guò)程中彈體及測(cè)試結(jié)構(gòu)產(chǎn)生的振動(dòng)信號(hào)[1];三是測(cè)試環(huán)境的外部噪聲和干擾信號(hào)。對(duì)于侵徹過(guò)載信號(hào)的處理,關(guān)鍵在于找到侵徹彈體的剛體過(guò)載響應(yīng)頻率,將侵徹過(guò)程中所產(chǎn)生的彈體振動(dòng)信號(hào),以及外部噪聲和干擾信號(hào)剔除,保留侵徹阻力形成的彈體剛體加速度信號(hào)。
國(guó)內(nèi)外學(xué)者對(duì)硬目標(biāo)侵徹信號(hào)的處理技術(shù)開(kāi)展了諸多研究[2-6],但這些信號(hào)處理方法都屬于固定閾值濾波,其特點(diǎn)都是基于特定條件找出與彈體相關(guān)的某一固定頻率閾值,當(dāng)頻域變換大于該閾值時(shí)保留原值,否則置零。這種濾波方法在一定條件下可以濾除彈體運(yùn)動(dòng)時(shí)所含高頻分量,對(duì)所得加速度曲線積分可以求得彈體運(yùn)動(dòng)的速度和深度等參數(shù),是目前侵徹過(guò)載信號(hào)處理的常用方法。1995年,Donoho[7]提出了小波閾值消噪方法。侵徹信號(hào)的小波分析一般通過(guò)選擇合適的小波對(duì)測(cè)試信號(hào)進(jìn)行分解,然后對(duì)分解后的高頻信號(hào)采用設(shè)定的閾值函數(shù)進(jìn)行濾波,最后對(duì)處理后的信號(hào)進(jìn)行重構(gòu)獲得消噪信號(hào)。其處理過(guò)程能夠去除信號(hào)攜帶的部分高頻噪聲,獲得較高的信噪比,但在實(shí)際濾波過(guò)程中忽視了侵徹信號(hào)本身的頻率分析,所以容易去除侵徹阻力本身形成的侵徹高頻信號(hào)。
本文在上述固定閥值濾波和小波閾值濾波兩種方法的基礎(chǔ)上,提出了融合總體經(jīng)驗(yàn)?zāi)B(tài)分解(ensemble empirical mode decomposition,EEMD)和小波變換(wavelet transform,WT)的侵徹信號(hào)處理方法。該方法首先對(duì)侵徹信號(hào)進(jìn)行總體經(jīng)驗(yàn)?zāi)B(tài)分解,然后計(jì)算分解后的各本征模態(tài)函數(shù)(intrinsic mode function,IMF)功率譜,并與原信號(hào)功率譜進(jìn)行比較,找出侵徹彈體的過(guò)載響應(yīng)頻率,然后對(duì)高頻IMF分量進(jìn)行小波閾值降噪,最后將降噪后的高頻IMF和反映信號(hào)特征的低頻IMF進(jìn)行組合以重構(gòu)信號(hào)。
1.1基于EEMD的本征模態(tài)分解
經(jīng)驗(yàn)?zāi)B(tài)分解法(empirical mode decomposition,EMD)是一種自適應(yīng)信號(hào)時(shí)頻處理方法,能夠很好地處理非平穩(wěn)、非線性信號(hào)[8]。EMD本質(zhì)是利用信號(hào)的特征時(shí)間尺度來(lái)獲得信號(hào)的固有波動(dòng)模式,據(jù)此對(duì)數(shù)據(jù)進(jìn)行逐級(jí)分解[9]。其處理過(guò)程如下:
(1)確定測(cè)試信號(hào)S(t)序列的所有極值點(diǎn),分別記為Smax序列和Smin序列。
(2)依據(jù)Smax和Smin序列用三次樣條插值函數(shù)擬合形成信號(hào)S(t)的上下包絡(luò)線,同時(shí)求得上下包絡(luò)線的均值M(t)。
(3)記H1(t)=S(t)-M(t),若原信號(hào)序列S(t)減去包絡(luò)均值M(t)后的新數(shù)據(jù)Hl不存在負(fù)的局部極大值和正的局部極小值,則將其作為S(t)的第一個(gè)IMF分量,若不滿足條件,說(shuō)明它還不是一個(gè)本征模態(tài)函數(shù),需要繼續(xù)進(jìn)行“篩選”。此時(shí)重新賦值,使H1(t)=S(t),重復(fù)上述步驟,直至得到第一個(gè)基本模態(tài)分量C1。
(4)將C1從原數(shù)據(jù)序列內(nèi)分離,得到殘余項(xiàng)R1=S(t)-C1,然后將殘余項(xiàng)R1作為新的數(shù)據(jù),按上述步驟分解,可以得到若干不可分解的固有模態(tài)函數(shù)Ci,直到Cn為止。
這樣原始信號(hào)序列就可以表征為
(1)
其中,Ci就是信號(hào)的本征模態(tài)分量,代表了第i個(gè)固有模態(tài)函數(shù),Rn為信號(hào)分解后的殘余分量,代表了信號(hào)中低頻干擾噪聲,在信號(hào)重構(gòu)時(shí)可以將其去除。
在上述分解中,得到合理的Ci取決于信號(hào)極值點(diǎn)的分布情況,如果信號(hào)極值點(diǎn)分布不均勻,分解的信號(hào)就會(huì)出現(xiàn)頻率混疊的情況。EEMD方法為在EMD的基礎(chǔ)上加入高斯白噪聲的信號(hào)處理方法。EEMD方法主體分解過(guò)程同上,不同之處在于在EMD分解前將N組不同的高斯白噪聲加入被分析信號(hào)中,在形成新的總體信號(hào)序列S1(t)、S2(t)、…、Sn(t)之后對(duì)新的信號(hào)序列進(jìn)行EMD分解,得到分解后的不同固有模態(tài)函數(shù)Ci1、Ci2、…、Cin和殘余項(xiàng)Ri1、Ri2、…、Rin,最后分別將經(jīng)過(guò)多次分解的結(jié)果平均,求得信號(hào)總體的C1、C2、…、Cn。由于EEMD加入了零均值噪聲的特性,經(jīng)過(guò)多次分解的結(jié)果平均后,信號(hào)噪聲將相互抵消[10],集成均值的結(jié)果就可作為最終結(jié)果。EEMD解決了EMD會(huì)出現(xiàn)頻域混疊的問(wèn)題,可以清楚地分解出信號(hào)的各階模態(tài),其分解流程如圖1所示。
圖1 信號(hào)的EEMD分解流程
1.2侵徹信號(hào)的EEMD分解
將設(shè)計(jì)的彈載信號(hào)測(cè)試裝置(圖2)加裝于直徑100 mm的彈體底部,用于無(wú)限混凝土板靶的侵徹實(shí)驗(yàn),實(shí)驗(yàn)測(cè)得的加速度信號(hào)如圖3所示。從圖3可以看出,加速度信號(hào)中含有較高的頻率分量,要想獲取反映彈體侵徹阻力的信號(hào),必須對(duì)測(cè)試信號(hào)進(jìn)行處理,即濾除測(cè)試中彈體的高頻振動(dòng)和噪聲信號(hào)。
利用圖1所示的侵徹信號(hào)EEMD分解流程對(duì)實(shí)驗(yàn)測(cè)得的信號(hào)進(jìn)行處理[11],分解后的各IMF分量如圖4所示。
圖2 侵測(cè)信號(hào)測(cè)試裝置
圖3 侵徹加速度曲線
圖4 實(shí)測(cè)信號(hào)EEMD分解后的IMF分量
從圖4可以看出,分解過(guò)程突出了信號(hào)的高頻振動(dòng)分量,有利于提取侵徹加速度信號(hào)中疊加的高頻振動(dòng)噪聲。信號(hào)的高頻部分主要集中在C1、C2、C3、C4分量中,從C5開(kāi)始,分解后的信號(hào)頻率大幅降低,且C6、C7總體保持一致,說(shuō)明6階分解后,所測(cè)信號(hào)已經(jīng)達(dá)到趨勢(shì)項(xiàng),無(wú)需繼續(xù)分解。但C6是否為信號(hào)的EEMD最優(yōu)分解仍需進(jìn)一步分析。
圖5 加速度信號(hào)功率譜
圖6 IMF分量的功率譜
本文提出結(jié)合信號(hào)功率譜分析的方法來(lái)進(jìn)一步確定信號(hào)的EEMD最優(yōu)分解層數(shù),因?yàn)楣β首V反映了信號(hào)對(duì)應(yīng)頻率的功率大小。對(duì)于侵徹信號(hào)而言,侵徹阻力所形成的減加速度信號(hào)屬于測(cè)試信號(hào)的主體,其頻率對(duì)應(yīng)的功率值最大,在功率譜圖中就會(huì)形成較大峰值。具體分析過(guò)程可以分為兩步:①計(jì)算信號(hào)的總體功率譜和EEMD分解后的IMF分量功率譜;②兩者進(jìn)行比較,若分解后的分量功率譜與總體功率譜一致,則表明EEMD分解已經(jīng)提取了侵徹時(shí)阻力形成的主體信號(hào)。圖5所示為侵徹信號(hào)的整體功率譜,可以看出信號(hào)在1083Hz和3813Hz處存在明顯峰值。將其與分解后的IMF分量功率譜比較(圖6),分析發(fā)現(xiàn)分解后的二階IMF功率譜最大能量在3813 Hz附近,四階IMF功率譜最大能量在1071 Hz處,此后五階IMF功率譜能量開(kāi)始急劇下降,沒(méi)有明顯能量集中。表明侵徹測(cè)試信號(hào)經(jīng)EEMD四階分解,已經(jīng)提取了侵徹過(guò)程中彈體阻力形成的主要信號(hào),所以C4為EEMD的最優(yōu)分解。這一結(jié)果與傳統(tǒng)傅里葉頻域分析(圖7)所得的主體頻率1071 Hz和3836 Hz非常接近,并且經(jīng)EEMD分解所得的主體頻率較傅里葉譜分析更加清晰、直觀。此時(shí)可以認(rèn)為1083 Hz為彈體侵徹過(guò)載時(shí)的固有頻率,3813 Hz可能是彈體的二次或多次諧振頻率。且由于EEMD分解時(shí)加入了噪聲零均值的特性,侵徹過(guò)程中環(huán)境和外部引入的噪聲經(jīng)EEMD分解平均后,噪聲將相互抵消。
圖7 侵徹信號(hào)頻域分析
為了比較本文方法和傳統(tǒng)小波分析在侵徹信號(hào)處理上的不同,對(duì)測(cè)得侵徹加速度信號(hào)進(jìn)行小波分解,其步驟如下: ①根據(jù)侵徹過(guò)載信號(hào)的特點(diǎn)和前人經(jīng)驗(yàn)選擇Daubechies II小波作為侵徹信號(hào)分析的母小波進(jìn)行分解[12],由分解后的小波系數(shù)(圖8)確定3層小波分解為侵徹信號(hào)的最佳分解;②對(duì)小波分解后信號(hào)的高頻噪聲進(jìn)行強(qiáng)制閾值濾波和默認(rèn)閾值濾波;③利用濾波后的高頻系數(shù)和第3層低頻系數(shù)進(jìn)行信號(hào)重構(gòu)。
分解后不同閾值濾波重構(gòu)所得的信號(hào)分別如圖9所示。其降噪后的效果由信噪比(signal to noise ratio,SNR)來(lái)衡量,SNR越大,降噪效果越好。SNR定義如下:
(2)
其中,RSNR為信噪比;P1為測(cè)試信號(hào)的功率,P2為噪聲信號(hào)的功率。由濾波后信噪比結(jié)果(表1)可知,小波分解默認(rèn)閾值濾波優(yōu)于強(qiáng)制閾值濾波。
(a)低頻系數(shù)
(b)高頻系數(shù)圖8 dB2小波分解后的系數(shù)
雖然上述信號(hào)的小波分解方法在時(shí)域和頻域具有局部化特征的特點(diǎn),能夠分離信號(hào)的高低頻分量,但對(duì)于侵徹信號(hào)而言,由于沒(méi)有太多的實(shí)驗(yàn)比對(duì)和大量的經(jīng)驗(yàn)數(shù)據(jù)支撐,分解時(shí)小波基的選擇只能按照少量實(shí)例作為參照,分解層數(shù)的確定也沒(méi)有與彈體固有頻率或功率結(jié)合來(lái)確定,所以對(duì)侵徹信號(hào)直接進(jìn)行單一小波分解,可能會(huì)濾去彈體侵徹過(guò)程中的高頻本征振動(dòng)信號(hào),獲得較差的濾波效果。
(a)強(qiáng)制閾值降噪
(b)默認(rèn)閾值降噪圖9 不同閾值的小波去噪
濾波方法強(qiáng)制消噪默認(rèn)閾值信噪比RSNR29.686645.2798
結(jié)合以上分析,針對(duì)侵徹非平穩(wěn)過(guò)載信號(hào),本文提出對(duì)EEMD分解后的高頻信號(hào)采用小波濾波進(jìn)行處理。首先對(duì)實(shí)測(cè)加速度信號(hào)采用EEMD分解,得到信號(hào)的IMF分量,然后求得各分量的功率譜,再通過(guò)與原始信號(hào)的功率譜比較,得出信號(hào)的EEMD分解尺度,確定彈體過(guò)載的響應(yīng)頻率。其次對(duì)EEMD分解后的高頻IMF分量進(jìn)行小波閾值降噪,濾除信號(hào)中包含的高頻振動(dòng)噪聲。最后將反映信號(hào)特征的低頻IMF分量與小波降噪后的高頻IMF分量進(jìn)行信號(hào)重構(gòu)。具體濾波流程如圖10所示。
根據(jù)EEMD分解流程,結(jié)合1.2節(jié)中侵徹信號(hào)的EEMD分解過(guò)程,按最優(yōu)分解層數(shù)(4層)分解,然后對(duì)分解所得的C1、C2、C3高頻信號(hào)進(jìn)行小波強(qiáng)制閾值濾波和默認(rèn)閾值濾波,最后將濾波后的高頻IMF信號(hào)和反映彈體侵徹的C4低頻分量進(jìn)行信號(hào)重構(gòu),得到EEMD和WT聯(lián)合濾波后的信號(hào)(圖11),消噪后信號(hào)的信噪比如表2所示??梢钥闯?,本文提出的EEMD聯(lián)合WT的濾波方法,能夠很好地確定侵徹信號(hào)的分解尺度,在保持侵徹過(guò)程高頻細(xì)節(jié)(如圖9a、圖11a的①、②處)的同時(shí),舍棄了信號(hào)低頻噪聲信號(hào)Rn,獲得了比小波濾波方法更高的信噪比。
圖10 EEMD和WT聯(lián)合濾波流程
(a)強(qiáng)制閾值降噪(b)默認(rèn)閾值降噪圖11 EEMD和WT聯(lián)合濾波信號(hào)
聯(lián)合濾波方法強(qiáng)制閾值默認(rèn)閾值信噪比RSNR31.132854.7859
圖12 彈體侵徹的速度時(shí)程曲線
圖13 彈體侵徹的位移時(shí)程曲線
對(duì)濾波后的加速度信號(hào)進(jìn)行積分,求得彈體侵徹過(guò)程的速度和位移時(shí)程曲線,見(jiàn)圖12、圖13。由圖12、圖13可知,彈體濾波后侵徹入靶速度為424 m/s,終止速度為0,彈體侵徹深度為1.114 m,與實(shí)驗(yàn)測(cè)得參數(shù)(表3)相近,其誤差在11%范圍內(nèi)。
表3 實(shí)驗(yàn)和本文方法所得侵徹參數(shù)及誤差
上述分析證明融合總體經(jīng)驗(yàn)?zāi)B(tài)分解(EEMD)和小波變換(WT)的侵徹引信信號(hào)處理方法能夠有效確定侵徹信號(hào)的分解層數(shù),得出信號(hào)侵徹過(guò)程中彈體的固有振動(dòng)頻率,在保留侵徹本身高頻分量的情況下,消除侵徹過(guò)程中頻域內(nèi)疊加的外部噪聲,處理后的信號(hào)能夠獲得較高的信噪比,積分后所得參數(shù)(速度、深度)都具有較高的準(zhǔn)確度,為侵徹實(shí)驗(yàn)測(cè)試過(guò)程中非平穩(wěn)隨機(jī)振動(dòng)信號(hào)的處理提供了一種較好的濾波方法。
[1]徐鵬.高g值沖擊測(cè)試及彈載存儲(chǔ)測(cè)試裝置本征特性研究[D].太原:中北大學(xué),2006.
[2]Franco R J,Platzbecker M R.Miniature Penetrator(MINPEN) Acceleration Recorder Development Test[R].New Mexico:Sandia National Laboratories,1998.[3]Forrestal M J,Frew D J,Hickerson J P,et al.Penetration of Concrete Targets with Deceleration-time Measurements[J].International Journal of Impact Engineering,2003,28(5):479-497.
[4]Rothacher T.Giger B:High G Ballistic Flight Data Recorder[C]//18th International Symposium on Ballistic.San Antonia,1999:379-386.
[5]張志安.硬目標(biāo)侵徹引信半實(shí)物仿真技術(shù)研究[D].南京:南京理工大學(xué),2007.
[6]王成華,陳佩銀,徐孝誠(chéng).侵徹過(guò)載實(shí)測(cè)數(shù)據(jù)的濾波及彈體侵徹剛體過(guò)載的確定[J].爆炸與沖擊,2007,27(5):416-419.Wang Chenghua,Chen Peiyin,Xu Xiaocheng.Filtering of Penetration Deceleration Data and Determining of Penetration Deceleration on the Rigid-body[J].Explosion and Shock Waves,2007,27(5):416-419.[7]Donoho D L.De-noising by Soft-thresholding[J].J. IEEE Transactions on Information Theory,1995,41(3):6123627.
[8]Wu Zhaohua,Huang N E.Ensemble Empirical Mode Decomposition:a Noise Assisted Data Analysis Method[J].Advances in Adaptive Data Analysis,2009,1(1):1-41.
[9]吳虎勝,呂建新,吳廬山,等.基于EMD 和SVM 的柴油機(jī)氣閥機(jī)構(gòu)故障診斷[J].中國(guó)機(jī)械工程,2010,21(22):2710-2713.
Wu Husheng,Lü Jianxin,Wu Lushan,et al.Fault Diagnosis for Diesel Valve Train Based on SVM and EMD[J].China Mechanical Engineering,2010,21(22):2710-2713.
[10]梁曲生,張西寧,沈玉綈.機(jī)械故障診斷理論與方法[M].西安:西安交通大學(xué)出版社,2009.
[11]郝慧艷,李曉峰,孫運(yùn)強(qiáng),等.侵徹過(guò)程彈體結(jié)構(gòu)響應(yīng)頻率特性的分析方法[J]. 振動(dòng),測(cè)試與診斷.2013,33(2):307-311.
Hao Huiyan, Li Xiaofeng, Sun Yunqiang,et al.Projectile Structural Response Frequency Characteristics Analysis Method in Penetration Process[J].Journal of Vibration, Measurement & Diagnosis,2013,33(2):307-311.
[12]朱麗, 婁國(guó)偉.自適應(yīng)閾值的小波去噪研究[J].制導(dǎo)與引信,2003,24(1):13-16.
Zhu Li,Lou Guowei.Wavelet De-noising Study with Adaptive Threshold Value[J].J. Guidance & Fuze,2003,24(1):13-16.
(編輯王艷麗)
Frequency Characteristics Analyses of Penetrating Missile and Penetration Overload Signal Processing
Zhao Haifeng1,2,3Zhang Ya1Li Shizhong1Guo Yan2
1.North University of China,Taiyuan,030051 2.Nanjing College of Information Technology,Nanjing,210023 3.University of Ottawa,Ottawa,K1N 6N5
In order to solve the hard target penetration overload signals de-noising problem,this paper proposed a joint filtering method based on EEMD and WT. Firstly,the measured signals were decomposed by EEMD method, and the IMF components could be got.Secondly,the original signals EEMD decomposition scale could be drawn through comparing the power spectrum of each components with the original signals’.Then,the high-frequency components of IMF were filtered based on the WT threshold.Finally,the signals were reconstructed by using the low-frequecy IMF components and the filtered high-frequency IMF components. Experiments show that proposed method can effectively extract the response frequency of missile body, eliminate high-frequency vibration and noise in the penetration. The results of the proposed method can get better signal to noise ratio(SNR) than that of WT.And the integrated velocity and displacement time-history curves are close to the experiments.
penetration overload;ensemble empirical mode decomposition(EEMD);wavelet transform(WT);intrinsic mode function(IMF);power spectrum
2015-03-24
國(guó)家自然科學(xué)基金資助項(xiàng)目(51275547);江蘇省第二批中青年骨干教師和校長(zhǎng)境外研修計(jì)劃資助項(xiàng)目(2012-13)
TN911.6DOI:10.3969/j.issn.1004-132X.2015.22.009
趙海峰,男,1981年生。中北大學(xué)機(jī)電工程學(xué)院博士研究生,南京信息職業(yè)技術(shù)學(xué)院機(jī)電學(xué)院講師,加拿大渥太華大學(xué)訪問(wèn)學(xué)者。主要研究方向?yàn)闄C(jī)電一體化、目標(biāo)識(shí)別與探測(cè)。發(fā)表論文10余篇。張亞,男,1964年生。中北大學(xué)機(jī)電工程學(xué)院教授、博士研究生導(dǎo)師。李世中,男,1969年生。中北大學(xué)機(jī)電工程學(xué)院教授。郭燕,女,1982年生。南京信息職業(yè)技術(shù)學(xué)院機(jī)電學(xué)院講師。