基于ICEEMDAN的蓮藕超聲探測(cè)回波HHT時(shí)頻能量分析

楊誠(chéng)，王艾倫，王計(jì)劃，蔡金

（中南大學(xué)a.輕合金研究院；b.高性能復(fù)雜制造國(guó)家重點(diǎn)實(shí)驗(yàn)室；c.機(jī)電工程學(xué)院，長(zhǎng)沙 410083）

0 引言

蓮藕主要分布于亞洲，是我國(guó)產(chǎn)量最高的水生蔬菜，占世界產(chǎn)量80%以上[1]。蓮藕是一種易種難收的農(nóng)產(chǎn)品，采挖難度大，導(dǎo)致挖藕成本約占其市場(chǎng)價(jià)格的40%。其中最難解決的問(wèn)題就是如何找到稀疏分布于泥土中的蓮藕。

蔡金[2]提出了使用超聲實(shí)現(xiàn)非接觸探測(cè)蓮藕，并得出合適的超聲探藕頻率段。王計(jì)劃[3]用等效流體模型對(duì)超聲探藕過(guò)程進(jìn)行分析，得出探測(cè)角度對(duì)探藕回波的影響。然而他們均未對(duì)信號(hào)進(jìn)行處理。蓮藕探測(cè)涉及多層介質(zhì)，導(dǎo)致蓮藕層的回波極其微弱，很容易被噪聲掩蓋。故回波電壓數(shù)值變化巨大，用該值判斷十分不準(zhǔn)確，易造成誤判。

超聲探藕的回波是非線性、非平穩(wěn)信號(hào)。傳統(tǒng)的時(shí)域分析如快速傅里葉變換、小波分析不能滿足要求。1998年，Huang等[4]提出了一種基于經(jīng)驗(yàn)?zāi)B(tài)分解（EMD）與希爾伯特頻譜分析（HT）集合的方式: 希爾伯特-黃變換（HHT）。HHT比起傅里葉變換和小波變換更合適分析局部[5]，更合適非平穩(wěn)信號(hào)分析。EMD是HHT的核心，但因其模態(tài)混疊、端點(diǎn)效應(yīng)及噪聲殘余問(wèn)題[6]，一直在研究發(fā)展至ICEEMDAN，這是目前效果最理想的方法[7]。

為減輕噪聲及雜波干擾，準(zhǔn)確又快速地從回波信號(hào)分析出蓮藕存在與否。本文將基于ICEEMDAN算法的HHT時(shí)頻能量分析用于蓮藕超聲探測(cè)中，解決從蓮藕回波信號(hào)判斷蓮藕是否存在的問(wèn)題，并用仿真實(shí)驗(yàn)回波及模擬實(shí)驗(yàn)回波來(lái)驗(yàn)證此方法的有效性。

1 仿真及實(shí)驗(yàn)信號(hào)的獲取

1.1 基本信息

超聲蓮藕探測(cè)的“水－泥土－蓮藕－空氣”多層介質(zhì)聲場(chǎng)模型如圖1所示，換能器垂直發(fā)射聲波經(jīng)過(guò)水層、泥土層、蓮藕外層到達(dá)蓮藕內(nèi)部的孔洞，遇到孔洞中的空氣反射絕大部分能量。其中探頭到泥土層0.35 m,蓮藕距離泥土表層0.15 m。

圖1 超聲蓮藕探測(cè)示意圖

1.2 聲場(chǎng)仿真獲取模擬回波數(shù)據(jù)

本文中使用COMSOL Multiphysics多物理場(chǎng)數(shù)值仿真軟件，并采用“聲—壓—電—固”多物理場(chǎng)耦合建立2D模型進(jìn)行仿真計(jì)算，得到換能器下表面的聲壓回波數(shù)據(jù)。

換能器定義為壓電材料PZT-5H。對(duì)于聲學(xué)仿真，為了計(jì)算足夠的精度，至少要滿足一個(gè)波長(zhǎng)內(nèi) 有4 ～6 個(gè) 以上的網(wǎng)格[8]。仿真具體的詳細(xì)參數(shù)見(jiàn)表1。

表1 超聲蓮藕探測(cè)仿真參數(shù)表

根據(jù)聲波在的傳播速度，以脈沖激勵(lì)開(kāi)始為時(shí)間零點(diǎn)，換能器將在0.47 ms接收到聲波在土層的回波，在0.67 ms接收到聲波在蓮藕及其孔洞處的回波。超聲蓮藕探測(cè)的仿真時(shí)間選取0～1.5 ms。仿真計(jì)算利用換能器表面的聲壓來(lái)表征實(shí)驗(yàn)中的回波電壓。

同時(shí)計(jì)算有藕和無(wú)藕時(shí)的探測(cè)回波，得到的換能器表面聲壓值隨時(shí)間變化的圖像為了排除始波的干擾，只截取0.4～1.5 ms的信號(hào)進(jìn)行研究，如圖2所示。

圖2 仿真有無(wú)蓮藕聲壓對(duì)比

有藕無(wú)藕的仿真信號(hào)差異小，難以發(fā)現(xiàn)蓮藕存在的回波特征。本文將采用ICCEMDAN找出信號(hào)的特征模態(tài)進(jìn)行分析。根據(jù)對(duì)超聲探藕的信號(hào)傳播過(guò)程分析，蓮藕存在時(shí)聲波經(jīng)過(guò)孔洞，與空氣界面接觸會(huì)發(fā)生反射和散射，產(chǎn)生的回波會(huì)使該時(shí)間接收到的回波有所增強(qiáng)，如圖3所示。而蓮藕不存在時(shí)，則不會(huì)有該過(guò)程。所以本文采用對(duì)比信號(hào)特征模態(tài)的瞬時(shí)能量譜進(jìn)行蓮藕存在與否的識(shí)別。

圖3 聲波經(jīng)過(guò)蓮藕孔洞發(fā)生反射及散射

1.3 超聲蓮藕探測(cè)實(shí)驗(yàn)平臺(tái)獲取實(shí)驗(yàn)回波數(shù)據(jù)

本文實(shí)驗(yàn)選用DYSTM32F4型號(hào)的超聲探測(cè)開(kāi)發(fā)板，中心頻率為400 kHz。實(shí)驗(yàn)采用的水聲換能器型號(hào)為DW-400-ELJZ。實(shí)驗(yàn)平臺(tái)采用了一個(gè)亞克力板制成的90 cm×90 cm×80 cm的水缸。在水缸中央預(yù)埋蓮藕至泥土下15 cm深，加水至50 cm處，用換能器放置在水下15 cm，正對(duì)蓮身探測(cè)，將探測(cè)信號(hào)導(dǎo)出進(jìn)行信號(hào)處理，如圖4所示。

圖4 模擬平臺(tái)實(shí)驗(yàn)

有藕無(wú)藕實(shí)測(cè)所得信號(hào)對(duì)比圖如圖5所示，由于本次實(shí)驗(yàn)的水聲換能器只有單層的壓電片，只能接收正壓電，所以結(jié)果與仿真略有不同，有大量有用成分集中在低頻分量。

圖5 實(shí)驗(yàn)有無(wú)藕回波電壓對(duì)比(0.4～1.5 ms)

2 ICEEMDAN分解HHT分析

根據(jù)上述分析，識(shí)別蓮藕存在與否的關(guān)鍵是去除噪聲及低頻雜波干擾，分析蓮藕的時(shí)域能量特征，就能判別蓮藕存在與否。下面對(duì)兩個(gè)信號(hào)進(jìn)行ICEEMAN分解，去除噪聲和雜波的干擾，然后進(jìn)行時(shí)域能量分析。

2.1 仿真計(jì)算信號(hào)的分析

對(duì)仿真計(jì)算得到信號(hào)回波進(jìn)行ICCEMDAN分解，得到各階的IMF信號(hào)。利用相關(guān)性和方差貢獻(xiàn)率兩個(gè)指標(biāo)來(lái)選取合適的IMFs。

相關(guān)性體現(xiàn)了模態(tài)分量與原始信號(hào)的相似程度[9]。計(jì)算公式為

為了更好地利用相關(guān)性和方差貢獻(xiàn)率濾掉高頻噪聲和無(wú)關(guān)的低頻分量，構(gòu)建一個(gè)評(píng)價(jià)IMFs性能的綜合評(píng)價(jià)指標(biāo)。選擇相關(guān)性權(quán)重為α=0.4，方差貢獻(xiàn)率權(quán)重為1-α，得到各IMFs的綜合性能指標(biāo)。

圖6為有藕時(shí)仿真計(jì)算經(jīng)過(guò)ICEEMDAN分解后特征模態(tài)分量的相關(guān)性、方差貢獻(xiàn)率及其綜合性能指標(biāo)圖。可以看出IMF3和IMF4的相關(guān)性和方差貢獻(xiàn)率都占據(jù)絕對(duì)優(yōu)勢(shì)（相關(guān)性系數(shù)＞0.4,方差貢獻(xiàn)率＞10%，綜合性能＞0.2），所以在本次分析中選用IMF3 和IMF4進(jìn)行分析。

圖6 仿真有藕IMFs評(píng)價(jià)指標(biāo)

信號(hào)經(jīng)過(guò)ICEEMDAN分解成為IMFs,經(jīng)過(guò)希爾伯特變換，可以得到每個(gè)IMFs分量的振幅、瞬時(shí)頻率和時(shí)間的關(guān)系H(ω,t)。對(duì)Hilbert譜中振幅的平方對(duì)頻率進(jìn)行積分就可得出瞬時(shí)能量譜IE(t)：

式中：ai(t)為信號(hào)實(shí)部；ωi(t)為信號(hào)虛部。

將IMF3和IMF4重構(gòu)得到的信號(hào)進(jìn)行時(shí)域能量分析，求得時(shí)頻能量譜。同理將無(wú)藕仿真的聲壓信號(hào)圖也經(jīng)過(guò)上述處理，選用IMF3和IMF4進(jìn)行重構(gòu)，得到瞬時(shí)能量譜，將有無(wú)蓮藕的瞬時(shí)能量譜進(jìn)行對(duì)比，如圖7所示?？梢钥闯鲈?.46～0.80 ms、1.1 ～1.2 ms 處，2 個(gè)信號(hào)有明顯的能量差異，而在1.3 ms之后，能量差異不再明顯。

圖7 仿真有藕無(wú)藕聲壓瞬時(shí)能量譜對(duì)比

為了能更簡(jiǎn)單直接地從信號(hào)中判定蓮藕是否存在，滿足智能化機(jī)械的需求，我們引入一個(gè)值：能量差值百分比Dif，如式（5）所示，用于判斷蓮藕是否存在。若該值大于設(shè)定的閾值，則蓮藕存在，反之則蓮藕不存在。

式中：e1為探測(cè)信號(hào)的平均瞬時(shí)能量；e2為無(wú)藕信號(hào)的平均瞬時(shí)能量。

考慮到不同水域超聲探測(cè)回波的差異、蓮藕大小及其氣孔位置的差別、水域中的雜物如石子雜草等的影響，經(jīng)過(guò)多次實(shí)驗(yàn)驗(yàn)證，區(qū)分蓮藕存在與否的閾值設(shè)定為12%比較合適。當(dāng)瞬時(shí)能量差值百分比大于12%時(shí)判斷為蓮藕存在，反之則不存在。

將仿真計(jì)算得到的有藕重構(gòu)信號(hào)平均瞬時(shí)能量代入e1，無(wú)藕的重構(gòu)信號(hào)平均瞬時(shí)能量代入e2，計(jì)算得Dif=20.82%。數(shù)值較大，有足夠的區(qū)分度來(lái)排除仿真實(shí)驗(yàn)中誤差的影響，足以作為判斷蓮藕存在與否的依據(jù)。當(dāng)不采用基于信號(hào)相關(guān)性和方差貢獻(xiàn)率的篩選方法，即使用全部的IMFs進(jìn)行時(shí)域分析，并求能量差值百分比，得Dif=15.28%。區(qū)分度減小，說(shuō)明進(jìn)行IMFs篩選是有顯著成效的。

2.2 模擬平臺(tái)實(shí)驗(yàn)信號(hào)分析

和上文分析仿真數(shù)據(jù)步驟一致。有藕實(shí)驗(yàn)回波數(shù)據(jù)經(jīng)過(guò)ICEEMDAN分解得到的12個(gè)IMF。根據(jù)2個(gè)信號(hào)基于相關(guān)性和方差貢獻(xiàn)率的綜合性能評(píng)價(jià)（如圖8）。有藕信號(hào)的IMF3、IMF8、IMF11、IMF12綜合性能評(píng)價(jià)大于0.2，無(wú)藕信號(hào)的IMF3、IMF4、IMF9、IMF13的信號(hào)綜合性能評(píng)價(jià)較高。含有原信號(hào)的信息量大，可以用于重構(gòu)信號(hào)進(jìn)行分析。

圖8 實(shí)驗(yàn)有藕無(wú)藕IMFs綜合評(píng)價(jià)

將上述分析中得到的有藕重構(gòu)信號(hào)和無(wú)藕重構(gòu)信號(hào)分別進(jìn)行時(shí)域能量分析，得到瞬時(shí)能量譜，對(duì)比如圖9所示?？梢?jiàn)在0.47～0.80 ms、1.0 ～1.2 ms處有藕信號(hào)的瞬時(shí)能量明顯較無(wú)藕信號(hào)大，驗(yàn)證了上文中的仿真計(jì)算是能夠符合實(shí)際的。

圖9 實(shí)驗(yàn)有藕無(wú)藕瞬時(shí)能量譜對(duì)比

對(duì)2 個(gè)信號(hào)的瞬時(shí)能量譜積分，進(jìn)一步計(jì)算實(shí)驗(yàn)信號(hào)的能量差值百分比，得到差值百分比Dif=41.8%，區(qū)分蓮藕是否存在效果十分明顯。而不經(jīng)過(guò)相關(guān)性及方差貢獻(xiàn)率篩選，采用全部IMF進(jìn)行計(jì)算得能量差值百分比僅為Dif=26.1%，表明文中的IMF篩選方法效果明顯。

3 結(jié)論

1）本文通過(guò)ICEEMDAN-HHT的時(shí)頻能量分析方法，可以簡(jiǎn)單地通過(guò)回波信號(hào)判斷是否有蓮藕存在。通過(guò)仿真和實(shí)驗(yàn)對(duì)比顯示，該方法能夠有效判別蓮藕存在與否。該方法具有可行性、易操作性，為后續(xù)蓮藕采挖的智能機(jī)械化生產(chǎn)有一定的借鑒作用。

2）ICEEMDAN-HHT作為新型的非平穩(wěn)信號(hào)分析方法，對(duì)非平穩(wěn)信號(hào)分析效果明顯，對(duì)信號(hào)的沖擊部分非常敏感，能夠提取出超聲蓮藕探測(cè)信號(hào)的回波特征，有效減少噪聲及雜波的干擾。使用基于相關(guān)性和方差貢獻(xiàn)率的綜合性能指標(biāo)可以有效地挑選信號(hào)的占據(jù)主要信息的IMFs進(jìn)行重構(gòu)，減少噪聲和無(wú)關(guān)的低頻分量的影響，提高分析的精度。