李紫鵬, 紀(jì)永強(qiáng), 郭兵勇, 楊坤德

(1.西北工業(yè)大學(xué) 航海學(xué)院, 陜西 西安 710072; 2.中國(guó)船舶集團(tuán)有限公司 系統(tǒng)工程研究院, 北京 100094; 3.西北工業(yè)大學(xué) 海洋研究院, 江蘇 蘇州 215400)

水聲目標(biāo)檢測(cè)與識(shí)別是水聲探測(cè)領(lǐng)域的關(guān)鍵技術(shù),也是水聲信號(hào)處理中的重點(diǎn)和難點(diǎn)[1-2],從輻射噪聲信號(hào)中提取目標(biāo)特征是其中最常用也是最有效的方法之一[3]。

但在現(xiàn)實(shí)環(huán)境中所獲取到的輻射噪聲信號(hào)往往會(huì)受到復(fù)雜海洋環(huán)境的影響,導(dǎo)致信號(hào)具有強(qiáng)非平穩(wěn)性和極低信噪比[4]。這種情況下,從原始輻射噪聲信號(hào)中直接提取目標(biāo)特征幾乎是不可能的。在這種情況下,基于小波變換(wavelet transform, WT)[5]、經(jīng)驗(yàn)?zāi)Ｊ椒纸?empirical mode decomposition, EMD)[6]、變分模式分解(variational mode decomposition, VMD)[7]等先進(jìn)信號(hào)分析技術(shù)被廣泛應(yīng)用于輻射噪聲信號(hào)的分析與處理中。

基于WT的信號(hào)處理方法,例如離散小波、連續(xù)小波、雙樹(shù)復(fù)小波、多小波等,最先被應(yīng)用于水聲信號(hào)處理中,并取得了較好的結(jié)果。但是小波變換的核心缺點(diǎn)在于其固定的基函數(shù)只能應(yīng)對(duì)單一的水聲目標(biāo)特征[8],盡管?chē)?guó)內(nèi)外學(xué)者通過(guò)構(gòu)造包含多個(gè)小波基函數(shù)的基函數(shù)族滿足不同水聲目標(biāo)特征的提取需求,但是固定的基函數(shù)仍限制了小波方法的發(fā)展。不同于使用固定基函數(shù)的傅里葉變換與小波變換,基于EMD和VMD的方法給出了自適應(yīng)的特征提取思路,并在水聲信號(hào)處理領(lǐng)域得到了最為廣泛的應(yīng)用。然而,EMD方法的性能受制于端點(diǎn)效應(yīng)和局部極值搜索算法的影響,且作為一種遞歸方法,每一次模式分解的誤差都會(huì)向下累積,最終導(dǎo)致分解失真[9];而VMD方法的性能則受制于分解參數(shù)的選取,在實(shí)際水聲目標(biāo)識(shí)別中無(wú)法準(zhǔn)確定義相關(guān)分解參數(shù)(模式數(shù)目與帶寬),從而降低了VMD分解的準(zhǔn)確性與可靠性[10]。因此,研究者們希望尋求一種更為有效的、自組織的、非遞歸的自適應(yīng)水聲信號(hào)特征提取方法。

受經(jīng)驗(yàn)?zāi)Ｊ椒纸?、變分模式分解等方法的啟發(fā),苗永浩[11]提出了一種全新的時(shí)間序列信號(hào)自適應(yīng)分解方法,即特征模式分解(feature mode decomposition, FMD),并將其應(yīng)用于旋轉(zhuǎn)機(jī)械的故障診斷中,取得了較好的特征提取效果與運(yùn)算效率。

本文提出一種參數(shù)優(yōu)化的特征模式分解方法,通過(guò)相關(guān)峭度對(duì)特征模式分解中的參數(shù)進(jìn)行優(yōu)化以實(shí)現(xiàn)最優(yōu)分解,并根據(jù)模式之間的相似性進(jìn)行模式融合,以提高微弱特征在信號(hào)中的表征,從而實(shí)現(xiàn)復(fù)雜水聲環(huán)境下的微弱目標(biāo)特征提取。

1 參數(shù)優(yōu)化的特征模式分解算法

1.1 特征模式分解

受經(jīng)驗(yàn)?zāi)Ｊ椒纸獾膯l(fā),特征模式分解也是一種非遞歸的信號(hào)分解方法,它通過(guò)設(shè)計(jì)最優(yōu)FIR濾波器組,將原始信號(hào)自適應(yīng)地分解為不同的模式。在尋找最優(yōu)濾波器組與分解結(jié)果的過(guò)程中,特征模式分解迭代更新濾波器系數(shù),并趨向于使子信號(hào)的相關(guān)峭度最大。因此,模式分解算法可以被視為求解約束下的最優(yōu)化問(wèn)題,即:

(1)

式中:C表示相關(guān)峭度(correlative kurtosis,CK),反映模式的脈沖性和周期性[12];uk表示第k個(gè)模式;fk(l)表示第k個(gè)濾波器,其窗長(zhǎng)為L(zhǎng);Ts是估計(jì)的特征周期。

采用迭代特征值分解算法求解公式,首先將公式(1)表示為矩陣形式:

(2)

模式的相關(guān)峭度則可以表示為:

(3)

式中WM是加權(quán)相關(guān)矩陣的中間變量:

(4)

結(jié)合式(1)、(2)、(4),可以得到FMD目標(biāo)函數(shù)的最終表達(dá):

(5)

式中RXWX和RXX分別是加權(quán)相關(guān)矩陣和相關(guān)矩陣。至此,第k個(gè)濾波器的系數(shù)可以根據(jù)式(5)進(jìn)行更新,并且可以通過(guò)優(yōu)化的FIR濾波器組來(lái)獲得單分量模式。

考慮到多調(diào)制現(xiàn)象,分解得到的單分量模式可能包含相同的目標(biāo)特征,并且其中一些模式是無(wú)意義的噪聲分量。為了消除冗余模式并減少計(jì)算量,FMD設(shè)計(jì)了一種迭代算法:在每次分解之后,首先選擇具有最大相似性的2個(gè)模式,并從這2個(gè)模式中剔除CK較小的模式;其次,通過(guò)剩余的模式來(lái)重構(gòu)新的輸入并進(jìn)行FMD分解;直到殘差信號(hào)滿足終止條件,對(duì)最終保留下的模式進(jìn)行解調(diào)分析,其過(guò)程如圖1所示。

圖1 特征模式分解算法流程

1.2 參數(shù)優(yōu)化的特征模式分解

1.2.1 濾波器數(shù)量與階次優(yōu)化

基于包絡(luò)信號(hào)獲得的峭度指標(biāo)可以很好地表征目標(biāo)特性,因此可以作為濾波器數(shù)量與階次的優(yōu)化目標(biāo):首先,將濾波器組階次設(shè)置為默認(rèn)值,并使用模式的包絡(luò)譜峭度來(lái)優(yōu)化和確定濾波器數(shù)量。然后,基于優(yōu)化的濾波器數(shù)量(固定值),使用包絡(luò)譜峭度優(yōu)化濾波器組階次。最后使用優(yōu)化后的參數(shù)對(duì)原始信號(hào)進(jìn)行FMD分解。具體過(guò)程如下。

1)固定濾波器組階次為30,濾波器數(shù)量的搜索范圍固定為M∈[2,15],針對(duì)每一個(gè)濾波器數(shù)量Mi,對(duì)原始信號(hào)進(jìn)行FMD分解,獲得Mi個(gè)模式,并對(duì)這Mi個(gè)模式進(jìn)行包絡(luò)分析獲取其包絡(luò)譜;

(6)

(7)

(8)

4)固定濾波器數(shù)量Mi,濾波器數(shù)量的搜索范圍固定為N∈[20,50],針對(duì)每一個(gè)濾波器數(shù)量Ni,對(duì)原始信號(hào)進(jìn)行FMD分解,獲得Ni個(gè)模式,并對(duì)這Ni個(gè)模式進(jìn)行包絡(luò)分析獲取其包絡(luò)譜;

1.2.2 模式融合與優(yōu)選

考慮到原始FMD算法實(shí)際上剔除較多了包含相似目標(biāo)特征的模式,因此在獲取最有分解參數(shù)后,僅對(duì)原始信號(hào)進(jìn)行一次FMD分解以獲得單組分模式,不再進(jìn)行相似模式篩除、信號(hào)重構(gòu)、再分解等步驟。

獲得若干單組分模式后,分別計(jì)算這些模式的包絡(luò)譜之間的皮爾遜相關(guān)系數(shù)[14-15]：

(9)

1)子信號(hào)必須包含相關(guān)峭度最大的模式;

2)如果有多個(gè)子信號(hào)滿足條件(1),則從中選取由最多模式融合得到的子信號(hào)。

1.2.3 模式融合與優(yōu)選

完整參數(shù)優(yōu)化的特征模式分解算法流程,如圖2所示。

圖2 參數(shù)優(yōu)化的特征模式分解算法流程

2 目標(biāo)特征提取方法仿真分析

在仿真信號(hào)中,假設(shè)輻射噪聲信號(hào)是由若干調(diào)幅調(diào)頻信號(hào)與高斯白噪聲混合而來(lái)[16]：

x(t)=∑iAie-2πfiβ(t-ti)cos(2πfi(t-ti))+n(t)

(10)

式中:Ai為振幅;β為振動(dòng)特征的衰減特性;n(t)為噪聲。特征頻率設(shè)置為10 Hz。測(cè)試信號(hào)的采樣率為5 000 Hz,采樣時(shí)長(zhǎng)為1 s。仿真信號(hào)如圖3(a)所示。加噪后的仿真信號(hào)信噪比為-15 dB,噪聲信號(hào)如圖3(b)所示。噪聲信號(hào)的頻譜和包絡(luò)頻譜如圖所示。分別為圖3(c)和圖3(d)。

從圖3(d)中可以看出,由于強(qiáng)噪聲的干擾,信號(hào)的包絡(luò)譜中僅出現(xiàn)了一個(gè)頻率為10 Hz譜峰,且強(qiáng)度較低,無(wú)法判定10 Hz就是目標(biāo)的特征頻率。

因此利用本文提出的方法對(duì)加噪后的仿真信號(hào)進(jìn)行處理,先設(shè)置濾波器組階次為30,當(dāng)濾波器數(shù)量M為13時(shí),模式包絡(luò)譜的峭度取得最大值,為3.23;然后固定濾波器數(shù)量N為13,當(dāng)濾波器數(shù)量階次N為41時(shí),模式包絡(luò)譜的峭度取得最大值,為3.46。根據(jù)上述獲得的最優(yōu)參數(shù)對(duì)仿真信號(hào)進(jìn)行分解,獲得13個(gè)模式。

然后計(jì)算這13個(gè)模式兩兩之間的相關(guān)系數(shù),選取大于0.65的模式進(jìn)行合并,獲得的子信號(hào)如表1所示,其中模式3為CK最大的模式,因此只考慮包含模式3的子信號(hào)。

表1 模式合并后的子信號(hào)

根據(jù)子信號(hào)的篩選規(guī)則,選擇子信號(hào)7進(jìn)行希爾伯特解調(diào)分析,結(jié)果如圖4所示?？梢钥闯?圖中目標(biāo)特征10 Hz及其高階諧波均十分明顯,且頻譜中的地毯值被控制在了極低的水平,說(shuō)明參數(shù)優(yōu)化的FMD方法能夠有效提取并增強(qiáng)噪聲信號(hào)中的微弱特征。

為進(jìn)一步評(píng)估所提出方法的有效性,將本文所提出方法與EEMD、VMD和固定參數(shù)的FMD方法進(jìn)行比較,將原加噪沖擊信號(hào)分別進(jìn)行EEMD、VMD和固定參數(shù)的FMD分解,其中EEMD選用的白噪聲標(biāo)準(zhǔn)差為0.2,混淆噪聲分解的迭代次數(shù)為50次,即取仿真信號(hào)進(jìn)行50加噪與EMD分解,將獲得模態(tài)的平均值作為最終分解結(jié)果;VMD分解中模式數(shù)目K和帶寬α分別5和2 000,即將原始仿真信號(hào)分解為5個(gè)帶寬為2 000的模式,且這些模式能夠最大程度上重構(gòu)原始信號(hào);FMD中的濾波器數(shù)目和濾波器組階次分別取5和30。然后以最大峭度為依據(jù)選取敏感模式進(jìn)行包絡(luò)譜分析,結(jié)果如圖5所示。

圖5 最優(yōu)子信號(hào)的包絡(luò)譜

EEMD、VMD和固定參數(shù)的FMD均能夠在不同程度上識(shí)別出特征分量,但是3種方法特征頻率與本文方法相比不明顯,且均無(wú)法識(shí)別更高次的諧波分量。其中EEMD方法的效果最差,其在于EEMD方法的前幾個(gè)寬帶模式包含了大量噪聲組分;VMD和固定參數(shù)的FMD方法中存在大量無(wú)意義的干擾,這些干擾的幅值甚至大于目標(biāo)特征頻率。因此,本文方法在微弱特征提取及表達(dá)方面具有更好的效果。為了進(jìn)一步驗(yàn)證該方法的有效性和工程應(yīng)用價(jià)值,下一節(jié)將采用該方法對(duì)真實(shí)海試數(shù)據(jù)進(jìn)行分析。

3 實(shí)驗(yàn)驗(yàn)證與分析

課題組于南海某處設(shè)計(jì)并開(kāi)展了一次水聲目標(biāo)輻射噪聲信號(hào)探測(cè)和采集實(shí)驗(yàn)。帶有3個(gè)水聽(tīng)器的測(cè)量船和目標(biāo)船如圖6所示。測(cè)量船的發(fā)動(dòng)機(jī)熄火,目標(biāo)船以勻速航行。在本次實(shí)驗(yàn)中共獲得了61組輻射噪聲信號(hào),這些信號(hào)的采頻為20 kHz,大多數(shù)輻射噪聲信號(hào)可以很容易地進(jìn)行分析,所包含的船舶軸系特征頻率也很明顯。

圖6 水聲目標(biāo)輻射噪聲信號(hào)探測(cè)和采集實(shí)驗(yàn)

為了驗(yàn)證本文所提出方法的有效性,仍然采用參數(shù)優(yōu)化的FMD方法、EEMD、VMD和固定參數(shù)的FMD4種方法分別處理所有的輻射噪聲信號(hào)。每種方法的性能如表2所示。

表2 不同信號(hào)分解方法的目標(biāo)特征提取能力

隨后,選取一條目標(biāo)特征微弱的輻射噪聲信號(hào)進(jìn)行分析。首先對(duì)其進(jìn)行降采樣,將采樣頻率降至5 000 Hz。降采樣后的信號(hào)時(shí)域波形如圖7(a)所示,其頻譜和包絡(luò)譜分別如圖7(b)和7(c)所示。圖7(a)中無(wú)明顯的周期性震蕩,圖7(b)和7(c)中也無(wú)顯著的特征譜線。

圖7 實(shí)測(cè)輻射噪聲信號(hào)時(shí)域波形、頻譜及包絡(luò)譜

使用本文提出的參數(shù)優(yōu)化FMD對(duì)原始信號(hào)進(jìn)行處理,當(dāng)濾波器數(shù)量M為10時(shí),模式包絡(luò)譜的峭度取得最大值,為3.78;然后固定濾波器數(shù)量N為10,當(dāng)濾波器數(shù)量階次N為37時(shí),模式包絡(luò)譜的峭度取得最大值,為4.59。根據(jù)上述獲得的最優(yōu)參數(shù)對(duì)仿真信號(hào)進(jìn)行分解,獲得10個(gè)模式。

然后計(jì)算這10個(gè)模式兩兩之間的相關(guān)系數(shù),選取大于0.65的模式進(jìn)行合并,其中模式6為相關(guān)峭度(CK)最大的模式,因此只考慮包含模式6的子信號(hào),即由模式3、6、7構(gòu)造的子信號(hào)是最優(yōu)子信號(hào)。對(duì)其進(jìn)行包絡(luò)解調(diào),其結(jié)果如圖8所示。

圖8 最優(yōu)子信號(hào)的包絡(luò)譜

圖8中有較為明顯的目標(biāo)軸系特征頻率及其2倍頻,作為譜線特征進(jìn)行目標(biāo)的檢測(cè)與識(shí)別。同樣地,將該方法與EEMD、VMD和固定參數(shù)的FMD方法進(jìn)行比較,其中EEMD的噪聲標(biāo)準(zhǔn)差為0.2,分解次數(shù)為50;VMD分解中模式數(shù)目K和帶寬α分別8和2 000;FMD中的濾波器數(shù)目和濾波器組階次分別取7和25。然后以最大峭度為依據(jù)選取敏感模式進(jìn)行包絡(luò)譜分析,結(jié)果如圖9所示。

圖9 3種方法的最優(yōu)子信號(hào)包絡(luò)譜

從圖9中可以看出,EEMD僅能提取出目標(biāo)軸系特征頻率,且幅值較低,譜線不明顯;VMD方法可以提取出目標(biāo)軸系特征頻率及其二倍頻,但是譜線不突出,且地毯值較高;固定參數(shù)的FMD也僅能提取出目標(biāo)軸系特征頻率,且極為模糊。因此,進(jìn)一步驗(yàn)證了本文方法在水聲目標(biāo)特征提取方面比EEMD、VMD和固定參數(shù)的FMD方法具有更好的效果。

4 結(jié)論

1) 本文將FMD這一自適應(yīng)信號(hào)分解方法應(yīng)用到水聲目標(biāo)檢測(cè)與識(shí)別領(lǐng)域,并克服了FMD的參數(shù)選擇問(wèn)題,能夠自適應(yīng)地獲取與待分解信號(hào)相匹配的濾波器數(shù)目M和濾波器組階次N。

2) 引入基于相似理論的模式融合方法,增強(qiáng)了微弱目標(biāo)特征的表達(dá),同時(shí)給出了融合后的子信號(hào)優(yōu)選準(zhǔn)則,降低了算法的計(jì)算復(fù)雜度。

3) 采用該方法分別對(duì)仿真和實(shí)驗(yàn)信號(hào)進(jìn)行分析,均成功識(shí)別出信號(hào)中隱藏的微弱目標(biāo)特征。通過(guò)與EEMD、VMD和固定參數(shù)的FMD方法的比較,進(jìn)一步驗(yàn)證了該方法在水聲目標(biāo)特征提取方面的優(yōu)勢(shì)。

該方法對(duì)水聲目標(biāo)的檢測(cè)與識(shí)別具有一定的工程應(yīng)用價(jià)值,在滿足輻射噪聲采集、存儲(chǔ)及算法運(yùn)行時(shí)間要求的前提下,該方法可用于水聲目標(biāo)的在線檢測(cè)與識(shí)別。下一步將針對(duì)該方法的適用范圍進(jìn)行研究,以確定該方法在不同強(qiáng)度、不同信噪比的輻射噪聲下的適用性。