劉海波,玄志武
(91550部隊94分隊,遼寧 大連 116023)
Hilbert-Huang變換[1]為對非平穩、非線性信號進行分析的有效方法,且具有自適應特點,能揭示非平穩信號的局部特性,已廣泛應用于機械故障檢測[2-3]、結構識別[4]、外彈道去噪[5]等領域。但由于 EMD 分解由實踐中總結而來,理論基礎尚不完善。實際應用中難以避免出現模態混疊現象,即一個IMF分量中包含不同時間尺度組分情況,削弱了EMD算法的頻帶分解能力,導致分解出的本征模態函數失去原有物理意義。
Li等[6]用小波對信號進行預處理,濾除高頻間斷事件,一定程度上改善了分解結果,但小波基的選擇需人工干預,削弱了EMD最突出的自適應性特點。Ryan等[7]在EMD分解前加入掩膜信號,對區分頻率較接近的頻帶有一定效果,但選取掩膜信號需信號的先驗信息,影響EMD算法的自適應性。
本文提出改進的EMD算法,先將提高原始信號頻帶分解能力的微分信號進行EMD分解,再將所得IMFs進行積分處理獲得原始信號有效IMFs。仿真計算與實際工程應用表明,該方法優于傳統的EMD算法,可有效分離頻率較接近的本征模態函數,增強EMD方法的頻帶分解能力。
Rilling等[8]對EMD算法的信號頻帶分解能力給出初步解答。Feldman[9]從理論層面進行深入研究得出,假設兩不同頻帶信號為s1,s2,其振幅與頻率分別為A1,A2與f1,f2,以頻率比f2/f1為橫坐標,振幅比A2/A1為縱坐標,得EMD算法頻帶分解能力見圖1。圖1中,1區表示A2/A1≤(f2/f1)-2,或f2/f1≤1.5無法用EMD算法進行分解;2 區表示(f2/f1)-2≤A2/A1≤2.4(f2/f1)-1.75,需幾次平移迭代才能將信號分解;3區表示A2/A1≥2.4(f2/f1)-1.75,一次迭代即能將信號分解。
由上述得,EMD算法的頻帶分解能力與信號振幅比A2/A1及頻率比f2/f1有關。若能提高信號的振幅比,使其從1區進入2區或3區,則可用EMD進行頻帶分解。

圖1 EMD算法頻帶分解能力Fig.1 Frequency decomposition ability of EMD
設輸入信號x(t)由兩種模態構成,即:

其中:f2>f1,A2<A1。信號振幅比y1可表示為:

對式(1)兩邊進行一階微分得:

其振幅比為:

由于f2>f1故y2>y1,即通過微分增大信號振幅比,可將1區部分信號“提升”到2區或3區(如圖1中箭頭所示),增強應用EMD算法進行頻帶分解能力,圖1演化為圖2。

圖2 經微分算子處理的EMD頻帶分解能力Fig.2 Frequency decomposition ability of EMD after differentiation
圖2中,① 區表示A2/A1≤(f2/f1)-3,或f2/f1≤1.5無法用微分算子處理后EMD算法將其分解;② 區表示(f2/f1)-3≤A2/A1≤2.4(f2/f1)-2.75,需幾次平移迭代,才能用微分算子處理后的EMD將信號分解;③區表示A2/A1≥2.4(f2/f1)-2.75,一次迭代即能用微分算子處理后EMD算法將信號分解。
圖2與圖1相比,①區部分顯著變小,說明EMD分解的盲區減小,即EMD的頻帶分解能力得以增強。
為定量分析EMD算法的頻帶分解能力,引入指標Error,表示原始信號成分與分解結果間在時域上的誤差:

其中:IMFi(t)表示第i階固有模態函數,n為EMD分解階數。
對同一待分解信號,采用不同分解模式,當Error趨于0時,說明本征函數IMFi(t)與信號si(t)較接近,頻帶分解正確。
實際工程應用中,可用差分方法代替微分過程,用累積求和方法代替積分過程。再次用EMD方法提取出累積求和后信號的高頻部分作為新的IMF,避免由替代方法所致誤差,而使累積求和后信號包絡均值不為零的情況發生。

圖3 算法流程圖Fig.3 Algorithm flowchart
正弦合成信號為:

采樣個數 210,圖 4(a)為s(t)高頻部分 0.2 sin(4πt),圖4(b)為s(t)低頻部分 sin(2πt)。s(t)振幅比A2/A1=0.2,頻率比f2/f1=2,因此A2/A1≤(f2/f1)-2,即圖1中信號s(t)處于1區,無法用EMD算法分離頻帶,已由圖 5結果得以證實,此時 Error為24.88。
對s(t)微分得:

工程應用中,信號往往含噪聲信息,直接進行微分處理會擴大噪聲,影響信號分解效果。當信噪比WS-NR為8.12時,信號s(t)經微分再經EMD分解的結果見圖8。對圖8各頻帶積分結果見圖9。圖8第三、四頻帶與圖6相似;圖9第三、四頻帶與圖7相似但存差異,說明信號s(t)中加入噪聲會對改進的EMD分解算法產生影響。
為進一步考察疊加噪聲的強度對算法影響,在合成信號s(t)中疊加不同強度白噪聲,利用本文改進的EMD算法進行處理。結果見表1,其中 WSNR表示s(t)疊加白噪聲后的信噪比。信噪比由13.64降到2.89低4/5,而 Error只增加2.4倍,說明該算法對噪聲適應性較好。

表1 疊加不同強度噪聲對算法結果影響Tab.1 Different noise strength influence on result

圖4 信號高低頻部分Fig.4 High and low frequency of the signal

圖5 直接EMD分解Fig.5 EMD decomposition directly

圖6 EMD微分后Fig.6 EMD decomposition after differentiation

圖7 積分Fig.7 Integration

圖8 加入噪聲的信號經微分再經EMD分解Fig.8 EMD after noising signal differentiation

圖9 圖8頻帶積分Fig.9 Frequencies Integration on figure
由于受日、月等天體影響,地球自轉不規律導致地球定向參數存在多種短、長周期變化。以分析地球定向參數測量數據[10]頻帶周期性規律為例,檢驗本文算法的有效性。

圖10 日長變化數據前4 096點Fig.10 First 4 096 points of LOD data
地球定向參數space2000_midnight.eop中的 LOD(length of day:日長變化)數據,記錄了從1976年9月28日至2001年1月6日共8 867個數據點。取前4 096點(圖10)進行EMD分解,結果見圖11。IFM1~IFM7每階固有模態函數均有不同振蕩周期,代表歲差、章動、極移等因素對LOD影響。

圖11 LOD數據EMD分解結果Fig.11 EMD Result of LOD data

圖12 用本文算法重新分解LOD數據Fig.12 Decomposition again using new method
由圖11虛線方框看出,IFM5中至少包含兩種尺度分量,并出現模態混疊現象,使IFM4亦受到影響。用本文改進的EMD算法重新對LOD數據進行分解,結果見圖12。與圖11虛線框內對比看出,該算法將原混疊模態成功分離,使各階IMF周期規律更明顯:IFM1為半月周期,IFM2為整月周期,IFM3為季度周期,IFM4為半年周期,IFM5為整年周期等。
由仿真分析與實際工程應用得出,本文所提算法可增強EMD頻帶分解能力。當滿足(f2/f1)-3≤A2/A1,且f2/f1>1.5條件時,可成功對混疊模態進行頻帶分解。
[1] Huang N E,Shen Z,Long S R,et al.The empirical mode decomposition and the hilbert spectrum for nonlinear and nonstationary time series analysis[J].Proceedings of the Royal Society of London,1998,454(1971):903-995.
[2]沈 路,楊富春,周曉軍,等.基于改進 EMD與形態濾波的齒輪故障特征提取[J].振動與沖擊,2010,29(3):154-157.
SHEN Lu,YANG Fu-chun,ZHOU Xiao-jun,et al.Gear faultfeature extraction based on improved EMD and morphological filter[J].Journal of Vibration and Shock,2010,29(3):154-157.
[3]李 慧,劉小峰,夏雨峰.基于諧波小波包變換的齒輪箱包絡解調分析[J].振動與沖擊,2012,31(12):129-134.
LI Hui,LIU Xiao-feng,XIA Yu-feng.A method of modal parameters identification using ambient vibration measurements for a large-scale bridge[J]. Journal of Vibration and Shock,2012,31(12):129-134.
[4]付 春,姜紹飛,杜 權.基于改進EMD的結構模態參數識別方法[J].武漢理工大學學報,2010,32(9):280-285.
FU Chun,JIANG Shao-fei,DU Quan.Structural modal parameters identification method based on improved EMD[J].Journal of Wuhan University of Technology,2010,32(9):280-285.
[5]郭小紅,徐小輝,趙樹強.基于經驗模態分解的外彈道降噪方法及應用[J].宇航學報,2008,29(4):1272-1275.
GUO Xiao-hong,XU Xiao-hui,ZHAO Shu-qiang.Noising reduction algorithm based on empirical mode decomposition(EMD)and applicationsin calculating trajectoriesof spacecraft[J].Journal of Astronautics,2008,29(4):1272-1275.
[6]Li H,Yang L,Huang D.The study of the intermittency test filtering character of Hilbert-Huang transform [J].Mathematics and Computers in Simulation,2005,70:22-32.
[7] Deering R,Kaiser J F.The use of a masking signal to improve empirical mode decomposition[C].Philadelphia,PA,USA:Proc. of IEEE International Conference on Acoustics,Speech,and Signal Processing,2005.
[8] Rilling G,Flandrin P.One or two frequency the empirical mode decomposition answers[J].IEEE Transaction on Signal Processing,2008,56(1):85-95.
[9]Feldman M.Analytical basis of the EMD:two harmonics decomposition[J]. Mechanical Systems and Signals Processing,2009,23(7):2059-2071.
[10] Gross R S.Combinations of earth orientation measurements:SPACE2000,COMB2000,and POLE2000[M].California:JPL Publication,2001.