張慧娟 李 冬

（中國人民解放軍92124部隊(duì) 大連 116023）

1 引言

工程測量中，由于受到各種環(huán)境因素的影響，測試信號(hào)常伴隨趨勢項(xiàng)，這嚴(yán)重影響時(shí)域分析和頻域分析某些參數(shù)的估計(jì)精度。趨勢項(xiàng)會(huì)產(chǎn)生低頻峰值，淹沒高頻分量信息，降低頻譜結(jié)果的可讀性，嚴(yán)重時(shí)使譜分析喪失原有的物理意義。如何將數(shù)據(jù)中包含的趨勢項(xiàng)提取的同時(shí)保證數(shù)據(jù)的可靠性和真實(shí)性，是振動(dòng)數(shù)據(jù)處理的重要的研究方向之一。

試驗(yàn)振動(dòng)數(shù)據(jù)中常見的趨勢項(xiàng)類型可以主要概括為三種，即線性趨勢項(xiàng)、指數(shù)型趨勢項(xiàng)和多項(xiàng)式趨勢項(xiàng)。其中對于線性趨勢項(xiàng)提取處理最為簡單，但大多數(shù)振動(dòng)數(shù)據(jù)含有的趨勢項(xiàng)較為復(fù)雜，線性方法往往難以滿足需要。文獻(xiàn)［1～4］中的方法，對指數(shù)性趨勢項(xiàng)的提取較為準(zhǔn)確，但由于其對低頻分量較為敏感，處理時(shí)會(huì)將部分真實(shí)信號(hào)當(dāng)作趨勢項(xiàng)，造成真實(shí)信號(hào)的損失；對多項(xiàng)式趨勢項(xiàng)的提取，取決于多項(xiàng)式的階數(shù)，如果對多項(xiàng)式階數(shù)估計(jì)不準(zhǔn)確，會(huì)造成較大誤差。經(jīng)驗(yàn)?zāi)B(tài)分解［5～6］（Empirical Mode Decomposition，EMD）是針對非平穩(wěn)信號(hào)的一種有效的分解方法，在結(jié)構(gòu)識(shí)別、外彈道降噪等領(lǐng)域得到了廣泛應(yīng)用［7～8］。為了得到高精度的振動(dòng)數(shù)據(jù)處理結(jié)果，文中設(shè)計(jì)的方法利用經(jīng)驗(yàn)?zāi)B(tài)分解結(jié)合對信號(hào)相關(guān)性的分析判決，提出一種適用于振動(dòng)數(shù)據(jù)趨勢項(xiàng)處理的新方法，既發(fā)揮了EMD對與被提取趨勢項(xiàng)類型具有普適性的優(yōu)點(diǎn)，同時(shí)考慮被提取信號(hào)的特征，保證了信號(hào)不被污染。

2 EMD基本原理和實(shí)現(xiàn)過程

經(jīng)驗(yàn)?zāi)B(tài)分解（EMD）分解是將目標(biāo)信號(hào)分解為若干按頻率高低排列的固有模態(tài)函數(shù)（Intrinsic Mode Function，IMF）。把隨機(jī)序列分解為若干固有模態(tài)函數(shù)IMF的線性組合稱為IMFs，對于每個(gè)單分量信號(hào)IMF必須保證滿足兩個(gè)條件：1）就整個(gè)時(shí)間序列來說，極值的個(gè)數(shù)和過零點(diǎn)的個(gè)數(shù)必須相等或只相差1；2）在任何一點(diǎn)，最大值和最小值形成的上下包絡(luò)線的均值為零。EMD與小波分析［9～10］相比具有無需驗(yàn)前選定基函數(shù)等先驗(yàn)信息的優(yōu)點(diǎn)。EMD算法的基本流程為

圖1 EMD算法流程

設(shè)X(t)是信號(hào)的時(shí)間序列，EMD方法對信號(hào)的分解模型可概括為：通過計(jì)算，得到測試信號(hào)X(t)的所有局部極值；用三次樣條函數(shù)擬合方法［11］求解測試信號(hào)中所有的極大值點(diǎn)構(gòu)成的上包絡(luò)線與所有極小值點(diǎn)構(gòu)成的下包絡(luò)線，分別記作emax(t)和emin(t)；上、下包絡(luò)的均值為原信號(hào)的平均包絡(luò)線：m(t)=[emin(t)+emax(t)]，原信號(hào)減去均值X(t)-m(t)=h(1)(t)得到新的信號(hào)h(1)(t)；判斷新信號(hào)是否滿足IMF的兩條性質(zhì)。若滿足，則h(1)(t)為IMF。否則，將h(1)(t)作為X(t)，重復(fù)以上過程，假定經(jīng)過k次分解后得到的h(k)(t)滿足IMF性質(zhì)，則X(t)的第一階IMF分量就是c1(t)=h(k)(t)；接下來記X1(t)=X(t)-c1(t)為原始信號(hào)，依次重復(fù)前述過程，得到第二個(gè)IMF分量，記作c2(t)，余項(xiàng)為X2(t)=X1(t)-c2(t)。重復(fù)上述步驟，直到的余項(xiàng)Xm(t)是一個(gè)單調(diào)信號(hào)或Xm(t)的值無法分解為止；分解完成后，可得到m個(gè)IMFs，c1(t)，c2(t)，…，cm(t)，最終的殘余項(xiàng)記為Xm(t)。由于實(shí)際數(shù)據(jù)中常含有不可消除的系統(tǒng)誤差，采用無偏估計(jì)時(shí)應(yīng)注意到，對于IMFs上、下包絡(luò)線均值為零的要求，顯然太過嚴(yán)格，結(jié)合振動(dòng)數(shù)據(jù)處理特性，在實(shí)際處理時(shí)，可適當(dāng)放寬，只要滿足包絡(luò)線均值小于某個(gè)閾值ε，即可認(rèn)為h(k)(t)滿足IMF性質(zhì)。其中ε的值，根據(jù)實(shí)際工程背景，取為ε∈[0.2,0.3]。

3 相關(guān)系數(shù)判決的EMD趨勢項(xiàng)提取方法

相關(guān)系數(shù)又稱皮氏積矩相關(guān)系數(shù)，是說明兩個(gè)變量之間相關(guān)關(guān)系密切程度的統(tǒng)計(jì)分析指標(biāo)。振動(dòng)信號(hào)中包含非常豐富信息，其中真實(shí)信號(hào)與原始信號(hào)具有強(qiáng)相關(guān)性。大量的工程實(shí)踐統(tǒng)計(jì)結(jié)果證明，通常在作EMD分解時(shí)，IMF分量與原信號(hào)之間的相關(guān)系數(shù) ||ζ＜0.1，認(rèn)為IMF分量與真實(shí)信號(hào)基本不再具有相關(guān)性，分解得到的余項(xiàng)或者最低1-2階即是信號(hào)中的趨勢項(xiàng)，下面以仿真信號(hào)s(t)=sin(2p×20t)+sin(2p×40t)為處理對象，說明利用相關(guān)系數(shù)判決的EMD趨勢項(xiàng)提取過程。

3.1 線性趨勢項(xiàng)

將線性信號(hào)作為趨勢項(xiàng)，疊加到處理對象中，構(gòu)造新的待分解信號(hào)［12］，得到的新信號(hào)為s(t)=sin(2p×20t)+sin(2p×40t)+x。

圖2 含有線性趨勢項(xiàng)的待分解信號(hào)

將混疊有線性趨勢項(xiàng)的信號(hào)作為新的處理對象，使用EMD方法進(jìn)行分解，原信號(hào)被分解為3個(gè)固有模態(tài)函數(shù)分量及余項(xiàng)，如圖3所示。

圖3 EMD法分解帶有線性趨勢項(xiàng)的信號(hào)

對EMD分解后得到的每個(gè)IMF分量分別求解其相關(guān)系數(shù)矩陣，結(jié)果為

相關(guān)系數(shù)矩陣對角線上的數(shù)，即為各個(gè)IMF分量與原信號(hào)之間的相關(guān)系數(shù)，得到的結(jié)果分別是：ζi1=0.691，ζi2=0.693，ζi3=-0.026；現(xiàn)認(rèn)為當(dāng)相關(guān)系數(shù) ||ζ≥0.1時(shí)，該IMF分量與原信號(hào)具有相關(guān)性。據(jù)此可得固有模態(tài)函數(shù)IMF1和IMF2與原信號(hào)具有較強(qiáng)的相關(guān)性，而IMF3與原信號(hào)之間幾乎不具有相關(guān)性。因此，可認(rèn)為IMF3與余項(xiàng)疊加后得到的信號(hào)分量就是趨勢項(xiàng)。

將使用EMD法提取出的趨勢項(xiàng)與仿真信號(hào)中加入的趨勢項(xiàng)進(jìn)行比較，其中用實(shí)線表示前者，而用虛線表示后者。比較的結(jié)果如圖4所示。

圖4 比較提取結(jié)果與加入的線性擾動(dòng)

使用EMD法提取出的趨勢項(xiàng)與仿真信號(hào)中加入的線性信號(hào)作差，得到的殘差如圖5所示。

圖5 提取趨勢項(xiàng)與真實(shí)趨勢項(xiàng)殘差

3.2 指數(shù)型趨勢項(xiàng)

仿真得到指數(shù)型趨勢項(xiàng)，將該趨勢項(xiàng)疊加到仿真信號(hào)中，新的信號(hào)作為待分解信號(hào)，如圖6所示。

圖6 含有指數(shù)型趨勢項(xiàng)的待分解信號(hào)

使用EMD法對新的信號(hào)進(jìn)行分解，分解后得到四個(gè)由高頻至低頻排列的IMF分量及余項(xiàng)，如圖7所示。

分別對原信號(hào)分解得到的各個(gè)IMF分量求解其相關(guān)系數(shù)矩陣，如下：

圖7 分解指數(shù)型趨勢項(xiàng)

由相關(guān)系數(shù)矩陣可知，各IMF分量與原信號(hào)之間 的 相 關(guān) 系 數(shù) 為 ：ζi1=0.690，ζi2=0.691，ζi3=-0.034，ζi4=0.083；根據(jù)相關(guān)系數(shù)判定趨勢項(xiàng)的方法可知，固有模態(tài)函數(shù)分量IMF1與IMF2與原信號(hào)的相關(guān)性較強(qiáng)，認(rèn)為是信號(hào)中的有用信息，而后兩項(xiàng)與原信號(hào)之間相關(guān)性較弱，所以后兩項(xiàng)與余項(xiàng)疊加后即認(rèn)為是仿真信號(hào)中帶有的趨勢項(xiàng)，分離出的趨勢項(xiàng)如圖8所示。

圖8 比較提取結(jié)果與加入的指數(shù)型擾動(dòng)

使用EMD法提取出的趨勢項(xiàng)與仿真信號(hào)中加入的指數(shù)型趨勢項(xiàng)作差，得到的殘差如圖9所示。

圖9 提取結(jié)果與真實(shí)值殘差

4 新方法提取實(shí)測數(shù)據(jù)趨勢項(xiàng)算例分析

利用某試驗(yàn)中獲取的一段加速度數(shù)據(jù)，對文中提出的方法加以驗(yàn)證，原信號(hào)的時(shí)間序列如圖10所示。

圖10 某試驗(yàn)加速度數(shù)據(jù)片段

該信號(hào)帶有未知的趨勢項(xiàng)。首先，將對該加速度信號(hào)進(jìn)行EMD分解，分解的結(jié)果如圖11所示，該信號(hào)被分解為5個(gè)IMF分量和余項(xiàng)?，F(xiàn)設(shè)定IMF分量與信號(hào)的相關(guān)系數(shù)閾值為0.1，認(rèn)為當(dāng)閾值小于0.1時(shí)，該分量與信號(hào)不具有相關(guān)性，應(yīng)將其作為趨勢項(xiàng)處理。所得到的5個(gè)IMF分量與信號(hào)的相關(guān)系數(shù)見表1。

表1 IMF分量與信號(hào)的相關(guān)系數(shù)表

由表1可知固有模態(tài)函數(shù)IMF5與待分解信號(hào)幾乎不相關(guān)，因此，將IMF5以及余項(xiàng)疊加后，即得到趨勢項(xiàng)。通過EMD分解，對加速度原始數(shù)據(jù)各IMF分量進(jìn)行提取，本例經(jīng)過5次分解成功提取出滿足限定條件的頻帶及趨勢項(xiàng)，如圖12所示。實(shí)例表明，EMD分解是目前提取數(shù)據(jù)序列趨勢項(xiàng)的十分有效的方法。

圖12 趨勢項(xiàng)提取結(jié)果

圖13 去除趨勢項(xiàng)后的真實(shí)信號(hào)

EMD方法在無需預(yù)設(shè)趨勢項(xiàng)類型的前提下，通過自然分解的方式，從實(shí)測信號(hào)中提取出趨勢項(xiàng)，再將剩余的部分疊加即可得到原始信號(hào)，該加速度的原始信號(hào)如圖13所示，EMD方法得到的趨勢項(xiàng)真實(shí)反映了原始測量信號(hào)的實(shí)際變化趨勢。

5 結(jié)語

文中對振動(dòng)數(shù)據(jù)中存在趨勢項(xiàng)的可能性及原因進(jìn)行了分析，并對利用經(jīng)驗(yàn)?zāi)B(tài)分解提取振動(dòng)數(shù)據(jù)中趨勢項(xiàng)的方法進(jìn)行了研究，仿真和實(shí)測數(shù)據(jù)計(jì)算結(jié)果證明，屬于真實(shí)信號(hào)部分的IMF分量與真實(shí)信號(hào)之間具有較大相關(guān)性，而屬于趨勢項(xiàng)的IMF分量與真實(shí)信號(hào)相關(guān)性較小。因此，可設(shè)定閾值，當(dāng)相關(guān)系數(shù)大于閾值時(shí)判定為真實(shí)信號(hào)分量，反之則為趨勢項(xiàng)部分。最后將判定為趨勢項(xiàng)的IMF分量相加得到最終趨勢項(xiàng)。該方法能夠有效地對信號(hào)中的趨勢項(xiàng)進(jìn)行提取，并盡可能保留原信號(hào)中的有用信息不受損失，保證了信號(hào)的完整性。根據(jù)算例中的信號(hào)特性，該方法可應(yīng)用于外彈道測量數(shù)據(jù)的趨勢項(xiàng)提取和遙測測量中的振動(dòng)信號(hào)處理工作。

參考文獻(xiàn)

［1］陳雋，李杰.振動(dòng)信號(hào)趨勢項(xiàng)提取的幾種方法及其比較［J］.福州大學(xué)學(xué)報(bào)，2005，33（s1）：42-45.

［2］李慧浩，許寶杰，左云波，等.基于小波變換和EMD方法提取趨勢項(xiàng)對比研究［J］.儀器儀表與分析監(jiān)測，2013（3）：28-30.

［3］郭淑卿.EMD的頻帶濾波篩分方法［J］.中國民航大學(xué)學(xué)報(bào)，2008，26（4）：30-33.

［4］高品賢.趨勢項(xiàng)對時(shí)域參數(shù)識(shí)別的影響及消除［J］.振動(dòng)、測試與診斷，1994（2）：20-26.

［5］Huang N E，Shen Z，Long S R.A NEW VIEW OF NON?LINEAR WATER WAVES：The Hilbert Spectrum1［J］.Annual Review of Fluid Mechanics，1999，31（1）：417-457.

［6］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］.Proceedings Math?ematical Physical&Engineering Sciences，1998，454（1971）：903-995.

［7］郭小紅，徐小輝，趙樹強(qiáng).基于經(jīng)驗(yàn)?zāi)B(tài)分解的外彈道降噪方法及應(yīng)用［J］. 宇航學(xué)報(bào)，2008，29（4）：1272-1275.

［8］李振興.結(jié)合經(jīng)驗(yàn)?zāi)B(tài)分解的振動(dòng)信號(hào)趨勢項(xiàng)提取方法［J］.飛行器測控學(xué)報(bào)，2011，30（1）：56-60.

［9］秦品樂，林焰，陳明.基于平移不變小波閾值算法的經(jīng)驗(yàn)?zāi)B(tài)分解方法［J］. 儀器儀表學(xué)報(bào)，2008，29（12）：2637-2641.

［10］馬社祥，劉貴忠，曾召華.基于小波分析的非平穩(wěn)時(shí)間序列分析與預(yù)測［J］.系統(tǒng)工程學(xué)報(bào)，2000，51（4）：305-311.

［11］Deering R，Kaiser J F.The use of a masking signal to im?prove empirical mode decomposition［C］//IEEE Interna?tional Conference on Acoustics，Speech，and Signal Pro?cessing，2005.Proceedings.IEEE，2005（4）：485-488.

［12］費(fèi)業(yè)泰.誤差理論與數(shù)據(jù)處理［M］.機(jī)械工業(yè)出版社，2015.