基于二次EEMD的Wigner-Ville分布旋轉(zhuǎn)機(jī)械故障信號分析及試驗研究

郭 奇，劉卜瑜，史立波，李 波

(1.燕山大學(xué) 環(huán)境與化學(xué)工程學(xué)院，秦皇島;2.中國石油吉林石化公司精細(xì)化學(xué)品廠，吉林)

經(jīng)過100多年的發(fā)展，F(xiàn)ourier變換不但已成為一個重要的數(shù)學(xué)分支，而且還是信號分析與處理的重要工具。然而，F(xiàn)ourier變換也存在不足，即缺乏時間和頻率的定位功能，在分辨率上的局限性以及對非平穩(wěn)信號分析的局限性［1］。工程上旋轉(zhuǎn)設(shè)備故障信號大多是非平穩(wěn)信號，故障信號分析的關(guān)鍵便成為如何對采集到的非平穩(wěn)信號進(jìn)行分析。因此大量的時頻分析方法應(yīng)運而生。其中，由Ville在1947年提出的Wigner-Ville分布憑借其優(yōu)良的時變特性和較高的時頻分辨率一直被應(yīng)用于信號的分析和處理領(lǐng)域。然而，由于其不滿足可加性而引起的交叉項卻成了不可忽視的缺陷，并使真實頻率成分之間產(chǎn)生虛假頻率成分，嚴(yán)重干擾對信號的準(zhǔn)確分析。在一定程度上限制了Wigner-Ville分布在實際中的發(fā)展應(yīng)用。

從交叉項產(chǎn)生的原因出發(fā)，可得到如下啟發(fā)，首先將疊加的信號分解成單個分量，分別對各分量進(jìn)行WVD計算，然后對計算結(jié)果求和，即可避免WVD在計算過程中不滿足可加性的缺陷。對信號進(jìn)行單分量分解有多種方法，其中近年來發(fā)展起來的經(jīng)驗?zāi)B(tài)分解(Empirical Mode Decomposition，EMD)和聚合經(jīng)驗?zāi)B(tài)分解(Ensemble EmpiricalMode Decomposition，EEMD)方法很受青睞。然而，這兩種方法也并非完美，EMD方法在分解過程中容易出現(xiàn)頻率混疊，影響分解的效果，不能得到真正的單一頻率成分，而EEMD分解雖然較EMD方法在抗混方面有所改善，但在消除Wigner-Ville分布交叉項方面還是不太理想。對此，本文提出利用二次EEMD消除Wigner-Ville分布交叉項方法，先用EEMD將信號分解成若干近似單一頻率的分量，對這些分量再次用EEMD方法進(jìn)行二次分解，得到效果較好的單一頻率成分，最后對新得到的單一分量分別進(jìn)行Wigner－Vill時頻分析，并對分析的結(jié)果求和，這樣便很好地消除了交叉項的影響。

1 Wigner-Ville分布

Wigner-Ville分布是一種二次型分布，由信號x(t)本身或其頻譜X(ω)定義為如下兩種等價的方式:

其中:* 表示復(fù)數(shù)共軛。

式(1)中的x*(t－τ/2)x(t+τ/2)稱為信號的瞬時相關(guān)函數(shù)，因此Wigner-Ville分布實質(zhì)上是對信號的瞬時相關(guān)函數(shù)的Fourier變換，其結(jié)果能夠反映信號的時頻特征。

WVD有很多優(yōu)點，但卻不滿足可加性，考慮如下信號:

將其代入式(1)，信號x(t)的WVD可寫為:

其中:

WVDx1x2(t，ω)與 WVDx2x1(t，ω)稱為互 Wigner-Ville分布，是復(fù)值的，并且可看出:

因此，WVDx1x2(t，ω)+WVDx2x1(t，ω)是實值的。式(4)可簡寫為:

由此可看出，兩個信號和的Wigner-Ville分布并不是簡單的兩個信號各自的Wigner-Ville分布之和，附加項2Re［WVDx1x2(t，ω)］通常稱為交叉項［1］。

2 EEMD基本原理

經(jīng)驗?zāi)B(tài)分解(Empirical Mode Decomposition，EMD)方法是由Huang等［3］提出的一種直觀的、直接的、后驗的和自適應(yīng)的分解方法。該方法通過一定的篩選條件將信號分解成為若干個相對平穩(wěn)的、互不相關(guān)的本征模函數(shù)(Intrinsic Mode Function，IMF)。但當(dāng)信號中存在異常事件(如間歇性成分、脈沖干擾等)時，分解得到的IMF就會存在模態(tài)混疊現(xiàn)象。由文獻(xiàn)［3］的研究成果知，模態(tài)混疊的產(chǎn)生與極值點選擇有關(guān)，利用極值點信息通過三次樣條插值生成的包絡(luò)為異常事件的局部包絡(luò)與真實信號包絡(luò)的組合，經(jīng)該包絡(luò)計算出均值并篩選出的IMF包含了信號固有模式和與異常事件模式，即模態(tài)混疊現(xiàn)象［2］。為改善該現(xiàn)象，Wu等［4］提出了 EEMD方法，通過加入高斯白噪聲解決EMD分解過程中出現(xiàn)的模態(tài)混疊問題。EEMD的計算方法為:

(1)初始化所要分解的IMF個數(shù)N和總的分解次數(shù)M(第一輪分解取m=1);

(2)對所輸入的信號加入隨機(jī)正態(tài)分布的白噪聲序列，并對被加入白噪聲的信號做歸一化處理;

(3)應(yīng)用EMD對歸一化的信號進(jìn)行分解，得到N個IMF;

(4)當(dāng)m＜M時，重復(fù)(2)、(3)步驟，并使m=m+1;

(5)將分解得到的M組IMFs對應(yīng)集成平均，當(dāng)M足夠大時，對應(yīng)剔除添加的白噪聲序列，而得到輸入信號的本征模函數(shù)組合，見式(9)、式(10):

其中:r(t)為余項［3－5］。

3 利用二次EEMD消除交叉項

WVD的交叉項是因被分析信號中含有多個頻率成分在計算過程中不滿足可加性而造成的。如果將被分析信號分解成為一系列單一頻率成分的信號，對其信號分別進(jìn)行WVD計算，再將計算結(jié)果求和，便可消除交叉項的影響［6］。

3.1 二次EEMD思想基礎(chǔ)

EMD分解IMF的過程是基于Huang的理想假設(shè)，但在實際信號分解過程中還存在比較棘手的問題。試驗證明，實際多分量信號經(jīng)EMD分解時容易產(chǎn)生虛假的IMF或者分解不出IMF(即模態(tài)混疊現(xiàn)象)［7］，因而分解出的IMF只能被稱為近似單頻率成分，對此文獻(xiàn)［7－8］均給出了詳細(xì)的論證。因模態(tài)混疊的存在，Wu等［5］提出了EEMD的思想以改善其分解效果。由于受求取信號包絡(luò)和篩選條件的限制，EEMD方法只是起到了改善分解效果的作用。即，經(jīng)過EEMD分解出的“單頻”成分還可能存在模態(tài)混疊(“單頻”成分中仍然可能存在兩個或兩個以上的頻率成分，在3.2節(jié)仿真結(jié)果中驗證)。因此本文認(rèn)為有必要對首次分解出的IMFs進(jìn)行第二輪的EEMD分解。由于EEMD分解過程依靠隨機(jī)加入均勻分布的白噪聲來改善其分解效果，并且信噪比可調(diào)，使二次分解后得到效果較好的IMFs成為可能。另外，由EEMD計算程序可知EEMD分解過程中強制規(guī)定了所要分解IMFs的數(shù)量，即TNM=fix(log2(N))－1(TNM為IMF數(shù)量，N為采樣點數(shù))，這對有些信號成分可能會造成過分解或者分解不完全?；诖?，在第二輪分解結(jié)束以后，本文對分解得到的IMFs以原信號為基準(zhǔn)進(jìn)行相關(guān)系數(shù)的篩選［1，9］，該方法可以有效地剔除由采樣率不足、樣條插值以及過分解等引起的偽分量，從而得到質(zhì)量較高的單頻率成分。至此可得出二次EEMD消除交叉項算法如下:

(1)用EEMD方法將被分析信號分解成為一系列滿足式(10)的單一頻率成分;

(2)利用EEMD方法對步驟(1)中所得到的IMF進(jìn)行二次分解，從而獲得新的IMFs;

(3)對步驟(2)中所得到的IMFs以原信號為基準(zhǔn)進(jìn)行相關(guān)系數(shù)的篩選以獲得質(zhì)量較高的單一頻率成分;

(4)對步驟(3)中所得到的一系列IMF分別進(jìn)行WVD計算:

(5)將各分量的WVD計算結(jié)果求和即得計算結(jié)果:

3.2 仿真信號模擬驗證

合成仿真信號:

可以看出信號中包含頻率分別為10 Hz的正弦信號，20 Hz的余弦信號，40 Hz的正弦信號。信號時域波形如圖1所示。圖2是經(jīng)過EEMD分解得到的IMFs，其數(shù)量由EEMD的計算程序強制規(guī)定(TNM=fix(log2(N))－1，N為采樣點數(shù))，第一條曲線和最后一條曲線分別為原始信號和余量。圖3為再次利用EEMD方法對圖2中第一個IMF分量分解后經(jīng)過相關(guān)系數(shù)法篩選出的IMFs，第一條和最后一條曲線分別為原始信號和余量。圖3的結(jié)果證實了經(jīng)過EEMD方法分解后的部分IMF是可以再分的，同時也證明了二次EEMD方法的可行性。圖4是對信號直接做WVD分析后的結(jié)果，從圖中看出，除在10 Hz、20 Hz及40 Hz處有真實的頻率成分外，在15 Hz、25 Hz及30 Hz附近也有頻率成分。而這些頻率成分在信號中并不存在，此即為WVD計算過程中出現(xiàn)的交叉項。為消除交叉項，采用3.1節(jié)中方法對信號進(jìn)行處理，圖5為一次EEMD處理后的結(jié)果，可以看出25 Hz處的虛假成分已被消除，但在15 Hz與30 Hz處虛假成分依然存在;圖6圖7為二次EEMD處理后的灰度圖和三維圖，在圖中可清楚地看到大部分虛假成分已被消除，用二次EEMD方法消除WVD交叉項效果很好。

圖4 未經(jīng)EEMD分解的Wigner-Ville分布Fig.4 The WVD of the simulated signal not based on EEMD

圖5 經(jīng)過EEMD分解的Wigner-Ville分布Fig.5 The WVD of the simulated signal based on EEMD

圖6 經(jīng)過二次EEMD分解的Wigner-Ville分布Fig.6 The 2 － D WVD of the simulated signal based on dual EEMD

圖7 經(jīng)過二次EEMD分解的Wigner-Ville分布三維圖Fig.7 The 3 － D WVD of the simulated signal based on dual EEMD

4 碰摩故障試驗研究

本次試驗采用DGZH－A型多功能轉(zhuǎn)子實驗臺模擬汽輪機(jī)轉(zhuǎn)子碰摩故障。將試驗臺左側(cè)第一個輪盤施加適當(dāng)?shù)妮d荷，使轉(zhuǎn)子旋轉(zhuǎn)過程中產(chǎn)生不平衡力，從而在軸承處產(chǎn)生摩擦，如圖8，轉(zhuǎn)子轉(zhuǎn)速設(shè)定為3000 r/min，采樣頻率為1280 Hz，采樣點數(shù)為1024，使用加速度型傳感器采集數(shù)據(jù)。

圖8 DGZH－A型多功能轉(zhuǎn)子實驗臺Fig.8 DGZH-A Multi-function rotor experiment bench

由動靜件局部碰摩的故障特征可知，局部摩擦引起的振動頻率中包含有不平衡引起的轉(zhuǎn)速頻率，即ω=3000/60=50 Hz，同時由于摩擦振動是非線性振動，所以還會包含有2ω、3ω、…一些高次諧波［10］。

為驗證二次EEMD方法分析的準(zhǔn)確度，采用京航公司HG8904C多通道數(shù)據(jù)采集故障診斷系統(tǒng)進(jìn)行故障數(shù)據(jù)采集與分析，圖9為HG8904C分析結(jié)果的時域圖和頻域圖的截圖?？梢钥闯鰧＜蚁到y(tǒng)提供的頻域圖中出現(xiàn)了 50 Hz、100 Hz、150 Hz、200 Hz 的頻率，此為動靜件碰摩故障特征頻率，可確定試驗臺此時存在動靜件碰摩故障。

圖9 HG8904C分析結(jié)果的時域圖和頻域圖Fig.9 The time-domain chart and frequency-domain chart of the fault signal givened by HG8904C

圖10為用Matlab還原的采集信號的時域圖，由圖11～圖14可以看出采用二次EEMD方法處理過的實際信號的Wigner-Ville分布中的交叉項已基本被消除，體現(xiàn)故障特征的轉(zhuǎn)速頻率及其倍頻(50 Hz、100 Hz、150 Hz、200 Hz)清晰可見，此結(jié)果與專家系統(tǒng)分析結(jié)果基本一致。由此證明了利用二次EEMD消除Wigner-Ville分布交叉項的方法是可行的。

圖10 信號的時域圖Fig.10 The time-domain chart of the fault signal

圖11 原始故障信號的Wgner-Ville分布灰度圖Fig.11 The WVD of the original fault signal

圖12 經(jīng)過EEMD分解的Wgner-Ville分布灰度圖Fig.12 The WVD of the fault signal based on EEMD

圖13 經(jīng)過二次EEMD分解的Wgner-Ville分布灰度圖Fig.13 The 2 － D WVD of the fault signal based on dual EEMD

圖14 經(jīng)過二次EEMD分解的Wgner-Ville分布的三維圖Fig.14 The 3 － D WVD of the fault signal based on dual EEMD

5 與EMD方法消除WVD交叉項比較

為克服WVD交叉項的干擾，對第4節(jié)中采集的數(shù)據(jù)利用EMD消除交叉項的方法進(jìn)行分析，見圖15，分解后經(jīng)過相關(guān)系數(shù)法篩選得到的IMFs(該方法可以得到質(zhì)量較好的IMFs，見文獻(xiàn)［9］)，第一條曲線和最后一條曲線分別為原始信號和余量。從圖中可以看出EMD方法成功分解出四個單一頻率成分，與故障特征頻率的數(shù)量相吻合。圖16、圖17，分別為經(jīng)EMD分解后WVD的灰度圖和三維圖，由圖看出雖然50 Hz、100 Hz、150 Hz、200 Hz的頻率成分清晰可見，但在每單一頻率成分之間出現(xiàn)許多不連續(xù)的干擾，由文獻(xiàn)［7－8］可知，此為EMD方法分解多分量信號時產(chǎn)生的虛假頻率成分。通過圖15、圖16與圖13、圖14對比，可以發(fā)現(xiàn)利用二次EEMD方法分解比利用EMD方法分解效果更好。

圖15 EMD分解后經(jīng)過篩選的IMFsFig.15 The IMFs given by EMD and correlation coefficient

6 結(jié)論

本文針對WVD分析時產(chǎn)生的交叉項問題，提出利用二次EEMD消除Wigner-Ville分布交叉項的方法。該方法在有效地消除WVD分布交叉項的同時保留了WVD分布的各種優(yōu)良特性。隨后結(jié)合汽輪機(jī)動靜件摩擦故障的實例，用二次EEMD消除Wigner-Ville分布交叉項的方法對故障的振動信號進(jìn)行分析，獲得了故障信號真實的頻率成分并對故障進(jìn)行診斷，所得結(jié)果與HG8904C多通道數(shù)據(jù)采集故障診斷系統(tǒng)所得到的分析結(jié)果一致。與利用EMD消除WVD交叉項的方法對比二次EEMD消除WVD交叉項的方法可獲得更好的效果，從而證明了該方法在實際應(yīng)用中的可行性。

［1］褚福磊，彭志科，馮志鵬，等.機(jī)械故障診斷中的現(xiàn)代信號處理方法［M］.北京:科學(xué)出版社，2009.

［2］曹沖鋒，楊世錫，楊將新.大型旋轉(zhuǎn)機(jī)械非平穩(wěn)振動信號的EEMD降噪方法［J］.振動與沖擊，2009，28(9):33－38.

［3］Huang N E，Shen Z，Long S R，et al.The empiricalmode decomposition and the hilbert spectrum for nonlinear and nonstartionary time series analysis［J］. Proc. Roy. Soc.London，1998，1971(454):903 －995.

［4］Wu Z H，Huang N E.A study of the characteris－tics of white noise using the empiricalmode decomposition method［J］.Proc.Roy.Soc.London，2005，l460A:1579 －1611.

［ 5］Wu Z H，Huang N E.Ensemble empirical mode decomposition:a noise-assisted data analysismethod［J］.Advances in Adaptive Data Analysis，2008，11:1 －4.

［6］韓中合，朱霄珣，李文華，等.基于EMD消除Wigner－vill分布交叉項的研究［J］.汽輪機(jī)技術(shù)，2010，52(3):211－214.

［7］沈國際，陶利民，陳仲生.多頻信號經(jīng)驗?zāi)B(tài)分解的理論研究及應(yīng)用［J］.振動工程學(xué)報，2005，18(1):91－94.

［8］全學(xué)海，丁宣浩，蔣英春.基于經(jīng)驗?zāi)B(tài)分解的篩選條件研究［J］.南陽師范學(xué)院學(xué)報，2010，9(3):22－26.

［9］林 麗，余 輪.基于相關(guān)系數(shù)的EMD改進(jìn)算法［J］.計算機(jī)與數(shù)字工程，2008，36(12):28－29.

［10］盛兆順，尹琦嶺.設(shè)備狀態(tài)監(jiān)測與故障診斷技術(shù)及應(yīng)用［M］.北京:化學(xué)工業(yè)出版社，2009.

［11］齊 天，裘 焱，吳亞鋒.利用聚合經(jīng)驗?zāi)B(tài)分解抑制振動信號中的模態(tài)混疊［J］.噪聲與振動控制，2010，2:103－106.

［12］張 超，陳建軍.EEMD方法和EMD方法抗模態(tài)混疊對比研究［J］.振動與沖擊，2010，29(S):87－90.