一種EMD改進方法及其在旋轉機械故障診斷中的應用

時培明，丁雪娟，李 庚，韓東穎

(1.燕山大學 電氣工程學院，秦皇島 066004;2.燕山大學 車輛與能源學院 ，秦皇島 066004)

旋轉機械故障診斷中，常常遇到諸如瞬變、調幅或調頻等非平穩、非線性信號，提取和分析這些信號中的特征信息是機械診斷的關鍵。EMD是Huang首先提出的一種新型的時頻分析方法［1］，能將復雜信號分解為有限個內稟模態函數(Intrinsic Mode Function，IMF)，適合處理非線性、非平穩信號。自EMD方法問世以來，它就引起了眾多學者的廣泛關注，已被應用于結構分析、設備診斷等領域［2-3］。

EMD方法的分析質量很大程度上取決于EMD分解的質量。而EMD分解由于采用三次樣條插值來獲取信號的瞬時平均，使得這種方法存在特殊的端點效應。另外，在Hilbert變換的過程中也存在端點效應，嚴重時會影響整個信號的分析結果。針對這一問題，研究人員已經提出了一些改進方法，如鏡像延拓［4］、神經網絡延拓［5］、相似極值延拓［6］、波形特征匹配延拓［7］、支持矢量回歸機［8］等。但是，這些延拓法存在一個共同的問題，就是在延拓后，信號的端點仍然不確定，這樣隨著分解過程的進行，包絡線兩端仍可能存在發散現象，并逐漸向內“污染”，導致其端點效應問題依然存在。

本文提出一種波形特征匹配延拓和余弦窗函數結合的EMD端點效應處理方法。該方法在保留了波形匹配延拓方法優勢的基礎上，通過窗函數處理，使延拓部分逐漸減小直到歸零，從而使延拓誤差減小［9-10］。通過仿真分析和不對中故障診斷實例，驗證了該方法的有效性。

1 改進的EMD時頻分析方法

1.1 EMD分解方法

EMD時頻分析方法適于處理非線性、非平穩信號。該方法包括兩個過程:經驗模態分解(EMD)和Hilbert變換。

EMD方法的分解步驟為:

(1)確定信號所有的局部極值點。

(2)利用三次樣條線將所有的局部極大值點連接起來形成上包絡線;將所有的局部極小值點連接起來形成下包絡線。上、下包絡線應該包絡所有的數據點。

(3)上、下包絡線的平均值記為m1(t)。

上述分解過程如圖1所示，圖(a)為原始信號，圖(b)中·表示極大值點，?為極小值點。圖(c)中a為上包絡線，b為下包絡線，c為上、下包絡線的均值。

圖1 EMD分解示意圖Fig.1 Diagram of the EMD decomposition

(4)信號x(t)和m1(t)的差值為第一個分量，h1(t)為:

如果h1(t)是一個IMF，那么h1(t)就是第一個 IMF分量。

(5)如果h1(t)不滿足IMF條件，把h1(t)作為原始數據，重復上述步驟，得到:

式中，m11(t)是h1(t)上下包絡的平均值。反復篩選k次后，使得h1k(t)為IMF分量，即:

令:

從原始信號中獲得的第一個IMF分量c1(t)應該包含信號最好的范圍或者最短的周期成分。

(6)從x(t)分離出c1(t)，我們得到:

將r1(t)看作原始數據重復以上步驟，得到x(t)的第2個IMF分量c2(t)。重復循環n次，得到信號x(t)的n個IMF分量。這樣就有:

當rn(t)成為一個單調函數不能再從中提取出滿足IMF條件的分量時，循環結束。這樣由式(5)和式(6)得到:

因此可以將信號分解為n個經驗模態，殘余函數rn(t)代表信號的平均趨勢。從而對每一個IMF進行Hilbert變換，得到Hilbert譜及其邊際譜。

對式(7)中的每一個內稟模態函數ci(t)作Hilbert變換得到:

構造解析信號:

于是得到幅值函數:

和相位函數:

進一步可以求出瞬時頻率:

這樣，可以得到:

這里省略了殘余分量rn，RP表示取實部。展開式(13)稱為Hilbert譜，記作:

再定義Hilbert邊際譜:

式中，T為信號的總長度。

1.2 EMD分解中的端點效應

EMD分解中由于無法保證端點處的極值點，導致求包絡平均過程中，會在樣條插值時產生數據的擬合誤差。并且隨著分解的進行，誤差不停積累，由端點處向內逐漸傳播，嚴重時會使分解的數據失去意義。

圖2為仿真信號x(t)的EMD分解結果，信號x(t)由兩個正弦信號和一個調幅信號組成，其表達式為

分解時沒有對原始數據進行任何處理。由圖2所示，信號x(t)分解出來的3階IMF都有嚴重的端點效應。

圖2 未經過處理的EMD分解結果Fig.2 The decomposition result of simulative signal without being processed

1.3 基于波形特征匹配延拓與余弦窗函數的EMD方法改進

基于波形特征匹配的數據延拓方法是通過采用信號內部和邊緣處變化趨勢最為相似的子波來對端點處數據進行延拓，它是一種自適應的方法。在具體實現中，通過計算波形匹配來量化兩端波形的變化趨勢。

如圖3所示，以左邊界第一個極值點為極大值為例，Mi、Ni(n=1，2，3，…)分別為曲線的極大值、極小值點，分別對應時間tmi、tni，S1為左端點，波形特征匹配延拓法以S1-M1-N1為邊界特征波形，在全部數據中找到與S1-M1-N1構成的三角形最接近的波形為匹配波形，如Si-Mi-Ni，從Si的前一點(右邊界為后一點)數據開始，向前(右邊界向后)延拓波形數據，使延拓數據符合信號的自然走向。

圖3 波形特征匹配原理圖Fig.3 Diagram of wave characteristic matching

具體步驟為:

(1)根據S1-M1-N1邊界三角形的時間坐標，尋找與S1對應的點Si，其對應的時間由式(17)求得:

圖4 信號x(t)的波形匹配延拓結果Fig.4 Results of wave characteristic matching

由于得到的tsi不一定在采樣點上，采用線性插值求其精確值Si。

(2)計算匹配誤差Ei

式中，最后一項為匹配波形的趨勢項，反映特征波形在曲線中相對極值點的位置。

(3)取Ei最小的波形為匹配波形，若存在多個相等的Ei，為獲得足夠的延拓數據供選擇，取與起始點距離最遠的波形為匹配波。

(4)自匹配波形與Si距離最近的數據點的前一點開始，將實際波形數據延拓到S1前，延拓點數據根據需要選擇，若信號中Si前數據點個數小于需要延拓的點數，可反復延拓此段波形。

(5)采用同樣的原理可以對右邊界數據進行延拓。

(6)采用延拓后的數據完成EMD分解。信號包絡插值計算采用延拓后的全部數據，而在EMD分解結束條件的計算中，仍采用有效數據。

(7)對EMD分解得到的IMF作Hilbert變換，此時數據長度為延拓后數據總長，按原始信號長度及位置取分析結果，得到有效數據的EMD分析結果。

在上述延拓過程中兼顧了數據的極值點及非極值點的波形數據，使得延拓數據與原信號交界處光滑過渡，避免了邊界處瞬時頻率的跳躍。

圖4是信號x(t)基于波形匹配延拓的結果。由圖可以看出延拓部分雖然與實際信號比較接近，但還是存在誤差。對延拓數據加余弦窗函數處理，控制端點效應向內“污染”，得到更準確的IMF。

定義余弦窗函數:

式中，L為信號延拓后的長度，A為信號兩端延拓中較長的延拓長度。

圖5 余弦窗函數Fig.5 The shape of cosine window

圖6 信號延拓加窗處理結果Fig.6 The results processed with the proposed method

余弦窗的形狀如圖5所示。

余弦函數窗將延拓誤差“控制”在信號兩端，使其無法(或以較慢速度)向數據內部發展，保證信號中部數據的正確分解。首先，將延拓信號x(t)與余弦窗函數進行內積運算，得到信號y(t)=〈x(t)，ω(t)〉。

然后，對處理后的信號y(t)進行EMD分解，再將分解得出的IMF的兩端去掉相應的延拓部分A。最后對減去延拓部分后的IMF進行邊界譜分析。該方法既考慮到了延拓誤差又考慮到了信號的完整性。

圖6為信號x(t)延拓后再加余弦窗處理的結果。從圖中可以看出延拓的部分逐漸減小直到歸零，從而使延拓部分誤差減小，為得到更準確的IMF提供了可靠條件。

2 仿真分析

對經過處理的信號x(t)進行EMD分解得到3階IMF，并與實際信號作比較如圖7。通過與圖2對比發現經過延拓加窗處理得到的IMF明顯更符合實際值，說明該方法對抑制端點效應有良好的效果。

圖7 延拓加窗后的分解結果Fig.7 The decomposition result of simulative signal with the proposed method

圖8、圖9分別是信號x(t)未經過處理和用改進方法得到的Hilbert譜與邊界譜。由圖8可以看出，未處理的Hilbert譜在信號兩端有比較嚴重的發散現象，而通過本文方法處理后的Hilbert譜效果有明顯的改善。圖9中(a)的幅值顯示明顯有很多微弱振蕩，尤其是在高頻部分尤為明顯。通過端點效應處理后得到的邊界譜(b)中，幅值振蕩基本消除，得到了更好的處理效果。

3 EMD改進方法在不對中故障診斷中的應用研究

圖8 Hilbert時頻譜對比圖Fig.8 The composition of Hilbert spectrum

圖9 邊際譜對比圖Fig.9 The composition of marginal spectrum

將本文提出的EMD改進方法應用于轉子不對中故障信號的特征提取及診斷。圖10是采集信號的時域波形，轉速為900 r/min，采樣頻率為768 Hz。圖11是故障信號的EMD分解結果，圖12、圖13分別是改進前后信號的Hilbert譜及邊際譜。

圖10 不對中故障的時域波形Fig.10 Time domain waveform of misalignment fault

圖11 故障信號的EMD分解結果Fig.11 The decomposition result of fault signal

圖11中的內稟模態函數IMF1，IMF2，IMF3被依次分解出來，分別對應著多倍頻、2倍頻和基頻振動模態，但是因為噪聲的存在，一定程度地影響了EMD的分解精度［11-12］。由Hilbert時頻譜圖可觀察到，基頻分量和2倍頻左右的分量在分析的時間內一直穩定存在，除此之外還有多倍頻分量，但并不穩定;比較圖12，改進后的Hilbert譜的端點發散現象得到了明顯的改善，尤其是信號的左端更為明顯。在Hilbert邊際譜的頻幅譜中觀察到，除了基頻還有其他倍頻信息存在，基頻分量和2倍頻分量占主導地位，且2倍頻分量強度并未超過基頻分量;比較圖13，改進后的邊際譜幅值微弱振蕩基本消除，且2倍頻分量得到了突出。根據以上分析診斷該轉子的故障為不對中故障。

圖12 故障信號的Hilbert譜Fig.12 The Hilbert spectrum of fault signal

圖13 故障信號的Hilbert邊際譜Fig.13 The marginal spectrum of fault signal

4 結論

(1)提出了一種抑制端點效應的新方法，首先利用波形特征匹配延拓對信號兩端進行延拓，然后根據延拓情況對信號加余弦窗進行處理，把延拓誤差控制在兩邊。最后在EMD分解后，只取信號的有效長度，提高了分解的精度。

(2)通過仿真信號的分析，證明了該方法能有效地抑制端點效應，為得到準確地邊際譜和Hilbert譜提供了保障。

(3)將改進方法應用到旋轉機械故障診斷中，通過對含有不對中故障的信號進行分析，證明了該方法能從非線性故障信號中得到真實有用的故障信息，實現旋轉機械故障的有效診斷。

［1］Huang N E.The empirical mode decomposition and the hilbert spectrum for nonlinear and non-stationary time series analysis［J］. Proc. R. Soc. Lond. A，1998，454:903-995.

［2］ Cheng J S，Yu D J，Tang J S，et al.Local rub-impact fault dianosis of the rotor systems based on EMD［J］.Mechanism and Machine Theory，2009，44:784-791.

［3］Yang B，Suh C S.Interpretation of crack-induced rotor nonlinear response using instantaneous frequency［J］.Mechanical System and Signal Processing，2004，18(3):491-513.

［4］韓建平，錢 炯，董小軍.采用鏡像延拓和RBF神經網絡處理EMD中端點效應［J］.振動、測試與診斷，2010，30(4):414-417.

HAN Jian-ping，QIAN Jiong，DONG Xiao-jun.Suppression of end-effect in empirical mode decomposition by mirror extension and radial basis function neural network prediction［J］.Journal of Vibration，Measurement＆Diagnosis，2010，30(4):414-417.

［5］胡勁松，楊世錫.EMD方法基于徑向基神經網絡預測的數據延拓與應用［J］.機械強度，2007，19(6):894-899.

HU Jin-song，YANG Shi-xi.Application of EMD method with data extension technique based on RBF neural network to time-frequency analysis［J］.Journal of Mechanical Strength，2007，19(6):894-899.

［6］沈 路，周曉軍，張志剛，等.Hilbert-Huang變換中的一種端點延拓方法［J］.振動與沖擊，2009，28(8):168-174.

SHEN Lu， ZHOU Xiao-jun， ZHANG Zhi-gang， et al.Boundary-extension method in Hilbert-Huang transform［J］.Journal of Vibration and Shock，2009，28(8):168-174.

［7］胡愛軍，安連鎖，唐貴基.HILBERT-HUANG變換端點效應處理新方法［J］.機械工程學報，2008，44(4):154-158.

HU Ai-jun， AN Lian-suo， TANG Gui-ji. New process method for end effects of Hilbert-Huang transform［J］.Journal of Mechanical Engineering，2008，44(4):154-158.

［8］程軍圣，于德介，楊 宇.基于支持矢量回歸機的Hilbert-Huang變換端點效應問題的處理方法［J］.機械工程學報，2006，42(4):23-31.

CHENG Jun-sheng，YU De-jie，YANG Yu.Process method for end effects of Hilbert-Huang transform based on support vector regression machines［J］.JournalofMechanical Engineering，2006，42(4):23-31.

［9］Qi K Y，He Z J，Zi Y Y.Cosine window-based boundary processing method for EMD and its application in rubbing faultdiagnosis［J］. MechanicalSystems and Signal Processing，2007，21:2750-2760.

［10］于慧君，陳章位，王慶豐.一種加窗重疊信號平滑連接方法及其在振動信號預處理中的應用［J］.振動與沖擊，2007，26(8):39-40.

YU Hui-jun，CHEN Zhang-wei，WANG Qing-feng.OLA signal smooth linking method and its application in vibration signal pre-processing［J］.Journal of Vibration and Shock，2007，26(8):39-40.

［11］ Jalan A K，Mohanty A R.Model based fault diagnosis of a rotor-bearing system for misalignment and unbalance under steady-state condition［J］.Journal of Sound and Vibration，2009，327(3-5):604-622.

［12］ Rybczyński J.The possibility of evaluating turbo-set bearing misalignment defects on the basis of bearing trajectory features［J］.Mechanical Systems and Signal Processing，2011，25(2):521-536.