周 智，朱永生，張優(yōu)云，袁 幸，師浩浩

(西安交通大學(xué) 現(xiàn)代設(shè)計及轉(zhuǎn)子軸承系統(tǒng)教育部重點實驗室，西安 710049)

滾動軸承故障信號常為非線性、非平穩(wěn)信號。傳統(tǒng)的頻域消噪方法難以達到理想效果。小波分析作為時頻域分析方法一種，具有多分辨率等優(yōu)點，小波閾值消噪在信號消噪中成為常用方法［1-3］。但該方法與頻域分析法一樣，均基于基函數(shù)，消噪效果仍依賴于基函數(shù)的選取。而基函數(shù)選取主要依賴設(shè)計者的主觀經(jīng)驗，易造成不確定性。經(jīng)驗?zāi)B(tài)分解(Empirical Mode Decomposition，EMD)方法由 Huang 等［4］提出，完全基于信號的自適應(yīng)分解，繼承了多分辨率優(yōu)勢，不存在基函數(shù)選取問題。Flandrin等［5-7］利用EMD對白噪聲特性進行研究，用統(tǒng)計分析進行信號消噪。但該方法僅當(dāng)某些基本模式分量(Intrinsic Mode Function，IMF)全部為噪聲成分，不含有用信號成分時有效，因此不能徹底抑制噪聲。文獻［8-10］受小波閾值消噪思想啟發(fā)，在對原始信號進行EMD分解后直接對每個IMF進行閾值消噪。但此方法未考慮IMF為零均值A(chǔ)M/FM的分量，損失了IMF的大部分有用信息，且仍沿用小波消噪中基于風(fēng)險估計思想的傳統(tǒng)閾值方法［3］。Kopsinis［11-12］提出 EMD 間隔閾值 (EMD Interval-Thresholding，EMD-IT)方法，從每個IMF中任意相鄰兩零點采樣間隔內(nèi)出發(fā)，對間隔內(nèi)唯一極值點進行判斷，再對該間隔內(nèi)序列進行閾值操作，亦沿用傳統(tǒng)的閾值策略，對滾動軸承非平穩(wěn)沖擊信號特征提取效果不佳?；跇O大似然估計(Maximum Likelihood Estimation，MLE)的閾值策略在小波消噪中得到較好應(yīng)用［13］，林京等［14-16］將其用于滾動軸承與齒輪箱中非平穩(wěn)沖擊信號特征提取，效果良好。因此本文將EMD-IT與MLE結(jié)合，用于滾動軸承故障診斷中，并取得良好效果。

1 EMD閾值消噪

EMD可將信號分解成一系列具有不同時間尺度的IMF分量，該分量需同時滿足:① 整個信號數(shù)據(jù)范圍內(nèi)，過零點數(shù)與極值點數(shù)相同或最多相差一個;② 局部最大值與局部最小值包絡(luò)線平均值為0。EMD算法步驟［4］:

(1)找出信號x(t)的所有局部極值點，將所有局部極大值點與局部極小值點分別三次樣條曲線連，得x(t)的上、下包絡(luò)線;

(2)記上、下包絡(luò)線的平均值為m1，x(t)與m1的差為h1，即h1=x(t)-m1;

(3)若h1非基本模式分量，令其為待處理信號，繼續(xù)以上步驟，得到h11=h1-m11;

(4)重復(fù)上述過程k次，直到h1k成為基本模式分量，即h1k=h1(k-1)-m1k;

(5)得第一個基本模式分量h(1)(t)=h1k;

(6)從信號x(t)中分離出基本模式分量h(1)(t)，即r1=x(t)-h(1)(t);

(7)將余量r1作為原始信號，重復(fù)以上步驟，依次得到其它基本模式分量h(1)(t)，h(2)(t)，…，h(L)(t)及余項rn(t)，即:

式中:L為基本模式分量總數(shù)。

滾動軸承故障信號可表示為:

式中:x(t)為實測信號，y(t)為純故障信號，n(t)為隨機噪聲，N為采樣點數(shù)。對x(t)進行EMD分解，得各基本模式分量h(i)(t)(i=1，2，…，L)，利用EMD直接閾 值 消 噪 (EMD Direct-Thresholding，EMD -DT)［8-10］，其硬、軟閾值形式為:

式中:Ti為第i層IMF中所用閾值。

圖1EMD-DT與EMD-ITFig.1 EMD-DT and EMD-IT

本文將EMD-IT閾值消噪用于滾動軸承故障診斷。其算法流程為［12］:

(1)對觀察到的原始信號x(t)進行EMD分解，得基本模式分量h(i)(t)，(i=1，2，…L);

2 MLE 閾值策略［13－16］

滾動軸承故障信號多為超高斯分布，或即稀疏分布。與信號成分相比，噪聲則非稀疏分布。因此，滾動軸承故障信號處理問題可簡化為稀疏分布形式信號的消噪或特征提取問題。稀疏性可采用熵測度量化，熵值越小，分布越稀疏，均勻性越差。MLE閾值策略即建立在“稀疏碼”基礎(chǔ)上。

假設(shè)m為反映故障的振動信號，v為服從均值0、方差σ2的高斯噪聲，則傳感器所測的信號s為:

軸承發(fā)生故障時，振動信號m常呈現(xiàn)出非高斯特性，Hyvarinen對m的能量密度函數(shù)p(m)進行估計:

式中:d為s的標(biāo)準(zhǔn)差，α為控制能量密度函數(shù)稀疏程度參數(shù)。為提取故障沖擊響應(yīng)信號，Lin等［14］將α設(shè)為1。由式(9)利用MLE方法得到m的估計為:

3 基于EMD-IT與MLE的故障診斷方法

在環(huán)境惡劣、白噪聲干擾強度大情況下，傳統(tǒng)的小波閾值消噪存在小波基函數(shù)難以選擇及硬、軟閾值消噪效果差等缺點。為凸顯故障信息，提高信噪比，本文結(jié)合EMD-IT與MLE策略，對滾動軸承故障信號進行特征信號提取，即在EMD-IT閾值消噪算法步驟(7)的閾值操作中采用MLE閾值方法:

4 數(shù)字仿真驗證

為驗證本文所提方法的正確性與有效性，采用數(shù)字仿真信號x(t)，其中y(t)模擬滾動軸承故障產(chǎn)生的沖擊響應(yīng)信號，調(diào)制頻率fr為10 Hz，載波頻率為500 Hz，并混入信噪比為-14 dB的白噪聲n(t)。

采樣頻率fs設(shè)為2 000 Hz，采樣點數(shù)N設(shè)為2 048，y(t)及x(t)時域波形分別如圖2、圖3所示，對x(t)的FFT分析見圖4。由圖3、圖4，x(t)時域波形中無法看出沖擊響應(yīng);而在頻譜圖中，500 Hz及由10 Hz調(diào)制產(chǎn)生的邊頻帶有用信息已被白噪聲淹沒，較難看出故障頻率。因此對原始信號用小波消噪處理，基函數(shù)選“db6”，結(jié)果見圖5。由圖5知小波消噪結(jié)果與y(t)相差較大，效果不佳。為克服該缺點，對原始信號分別采用基于EMD-IT的硬、軟閾值與極大似然估計消噪，結(jié)果如圖6～圖8所示。比較分析可見，采用EMD-IT的硬、軟閾值消噪效果較好，提取到部分故障沖擊響應(yīng)信號，但仍消噪有限，而EMD-IT-MLE閾值消噪方法更具優(yōu)勢，能更好提取出沖擊信號。分析比較圖8、圖2可知，本文方法提取的故障特征基本與構(gòu)造的故障特征相同，但在沖擊

圖2 y(t)時域波形圖Fig.2 Time-domain waveform of y(t)

圖3 x(t)時域波形圖Fig.3 Time-domain waveform of x(t)

圖4 的頻譜圖Fig.4 Spectrum of x(t)

圖5 小波消噪結(jié)果Fig.5 Results by using wavelet denoising

圖6 EMD-IT-硬閾值消噪Fig.6 Signals by using EMD-IT hard threshold denoising

圖7 EMD-IT-軟閾值消噪Fig.7 Signals by using EMD-IT soft threshold denoising

圖8 EMD-IT-MLE閾值消噪Fig.8 Signals by using EMD-IT-MLE

圖9 本文方法消噪后頻譜圖Fig.9 Spectrum of denoising signal by using the proposed method

圖10 本文方法消噪后包絡(luò)譜圖Fig.10 Envelope spectrum of denoising signal by using the proposed method

響應(yīng)衰減至幅值較小處，由于其值小于閾值而被置零，但不影響本文方法提取故障特征的優(yōu)勢。對本文方法消噪后信號進行頻譜分析(圖9)，比較消噪前情況(圖4)可知，消噪后500 Hz及因故障頻率調(diào)制產(chǎn)生的邊頻帶信息凸顯，有效抑制了噪聲干擾。包絡(luò)譜分析結(jié)果(圖10)顯示，故障頻率10 Hz及倍頻成分清晰可見。

5 實驗驗證

為進一步驗證本文方法的有效性，采用實驗數(shù)據(jù)進行驗證。實驗用ABLT-1A型軸承加速疲勞壽命試驗機基本原理見圖11。由HK9141加速度傳感器、PCI-1715U數(shù)據(jù)采集卡、工控機構(gòu)成數(shù)據(jù)采集硬件平臺，利用NI-Labvivew軟件實現(xiàn)數(shù)據(jù)采集。實驗工況為:轉(zhuǎn)速3 000 r/min，徑向載荷16.17 kN，無軸向加載。四個測試軸承均為深溝球軸承6309，參數(shù)為:滾動體數(shù)8個，直徑 17.5 mm;軸承節(jié)徑 72.5 mm;接觸角 0°。采樣頻率48 kHz，采樣點數(shù)為32 768。由滾動軸承故障特征頻率理論計算公式得到外圈故障特征頻率fw=151.3 Hz，內(nèi)圈故障特征頻率fn=247.9 Hz，滾動體故障特征頻率fg=97.5 Hz，保持架故障特征頻率fb=18.9 Hz。

原始信號x(t)的時域波形見圖12。為提高信噪比，對原始信號用“db6”小波基函數(shù)進行小波消噪，見圖13。

圖11 實驗測試平臺Fig.11 Experiment platform

圖12 實驗數(shù)據(jù)原始時域波形Fig.12 Original time-domain waveform of data

圖13 實驗數(shù)據(jù)小波消噪效果圖Fig.13 Experiment result by using wavelet denoising

圖14 實驗數(shù)據(jù)本文方法消噪效果圖Fig.14 Result by using the proposed method

用本文方法對原始信號進行處理，結(jié)果見圖14。對比可知，基于EMD-IT與MLE閾值消噪能有效提取滾動軸承故障沖擊信號。對消噪后信號包絡(luò)譜分析結(jié)果見圖15。

圖15 消噪后信號的包絡(luò)譜圖Fig.15 Envelope spectrum of denoising signal

由圖15看出，外圈故障特征頻率152.3 Hz≈fw=151.3 Hz及倍頻 2fw≈306.2 Hz、3fw≈458.5 Hz、4fw≈612.3 Hz、5fw≈764.6 Hz、6fw≈918.5 Hz 清晰可見，其中49.8 Hz為轉(zhuǎn)頻fz;判斷為外圈故障。實物檢查發(fā)現(xiàn)外圈有疲勞剝落凹坑(圖16)，凹坑寬6.25 mm、長3.59 mm。由此表明，本文所用故障診斷方法有效。

圖16 外圈故障滾動軸承Fig.16 Rolling bearing with outer race fault

6 結(jié)論

(1)本文將EMD間隔閾值消噪與極大似然估計相結(jié)合，并應(yīng)用于滾動軸承故障診斷中。與小波消噪不同，該方法完全基于信號的自適應(yīng)分解，同時結(jié)合極大似然估計提取故障沖擊特征優(yōu)點，克服了小波基函數(shù)選取困難及傳統(tǒng)硬、軟閾值消噪效果欠佳的缺點;

(2)通過對仿真信號及實驗數(shù)據(jù)分析發(fā)現(xiàn)，在滾動軸承故障診斷中，本文方法明顯優(yōu)于小波消噪法，能消除噪聲干擾，凸顯故障沖擊特征。

［1］Mallat S.A wavelet tour of signal processing［M］.New York:Academic，1999.

［2］Antoniadis A，Bigot J.Wavelet estimators in nonparametric regression:a comparative simulation study［J］.Journal of Statistical Software，2001，6(6):1-83.

［3］Donoho D L，Johnstone I M.Ideal spatial adaptation by wavelet shrinkage［J］.Biometrika，1994，81(3):425-455.

［4］ 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］.Proc.R.Soc.Lond.A，1998，454(1971):903-995.

［5］ Flandrin P，RillingG， GoncalvesP. Empiricalmode decomposition as a?lter bank［J］.IEEE Signal Process Letters，2004，11(2):112-114.

［6］ Wu Z，Huang N E.A study of the characteristics of white noise using the empirical mode decomposition method［J］.Proc.Roy.Soc.London A，2004，460(2046):1597-1611.

［7］ Huang N E，Wu Z.Hilbert-Huang trans-form and its applications［M］.Singapore:World Scientific，2005.

［8］Boudraa A O，Cexus J C.Denoising via empirical mode decomposition［M］.Proc.Isccsp，2006.

［9］ Mao Y，Que P.Noise suppression and flaw detection of ultrasonic signals via empirical mode decomposition［J］.Russian Journal of Nondestructive Testing，2007，43(3):196-203.

［10］鄒 清，湯井田，唐 艷.Hilbert-Huang變換應(yīng)用于心電信號消噪［J］.中國醫(yī)學(xué)物理學(xué)雜志，2007，24(4):309-312.

ZOU Qing，TANG Jing-tian，TANG Yan. Hilbert-Huang transform for ECG de-noising［J］.Chinese Journal of Medical Physics，2007，24(4):309-312.

［11］ Kopsinis Y，McLaughlin S.Empirical mode decomposition based soft-thresholding［C］. Proceedings ofthe 16th European Signal Processing Conference，2008.

［12］ KopsinisY， McLaughlin S. DevelopmentofEMD-based denoising methods inspired by wavelet thresholding［J］.IEEE Transactions on Signal Processing，2009，57(4):1351-1362.

［13］ Hyvarinen A.Sparse code shrinkage:denoising of nongaussian data by maximum likelihood estimation ［J］.Neural Computation，1999，11(7):1739-1768.

［14］ Lin J，Zuo M J，F(xiàn)yfe K R.Mechanical fault detection based on the wavelet de-noising technique［J］.Journal of Vibration and Acoustics，Transactions of the ASME，2004，126(1):9-16.

［15］林 京.基于最大似然估計的小波閾值消噪技術(shù)及信號特征提?。跩］.儀器儀表學(xué)報，2005，26(9):923-927.

LIN Jing.Wavelet de-noising based on maximun likelihood estimation and its application for feature extraction［J］.Chinese Joumal of Scientifie Instrument，2005，26(9):923-927.

［16］He W，Jiang Z N，F(xiàn)eng K.Bearing fault detection based on optimal wavelet filter and sparse［J］.Measurement，2009，42(7):1092-1102.