仝 蕊， 康建設(shè)， 李寶晨， 張星輝

(1. 陸軍工程大學(xué)石家莊校區(qū)裝備指揮與管理系， 河北 石家莊 050003；2. 陸軍工程大學(xué)科研學(xué)術(shù)處， 江蘇 南京 210007)

作為機械傳動系統(tǒng)的核心部件，齒輪的損傷或失效會影響整個系統(tǒng)的穩(wěn)定運行,而齒輪發(fā)生故障時的振動信號呈現(xiàn)非平穩(wěn)、非線性等特征[1]。為提高故障診斷的準(zhǔn)確程度，需要尋求一種更為有效的非線性信號處理方法。

典型的非線性信號處理方法有小波分析、經(jīng)驗?zāi)B(tài)分解(Empirical Mode Decomposition， EMD)、局域均值分解(Local Mean Decomposition，LMD)和形態(tài)學(xué)分析(Morphological Analysis，MA)等。小波分析有自身缺陷，在重構(gòu)采樣中會遺漏部分信息[2]。EMD雖具有自適應(yīng)性，可根據(jù)信號特點進(jìn)行分解，但存在模態(tài)混疊現(xiàn)象[3]。LMD與EMD類似，是一種能自適應(yīng)處理非平穩(wěn)信號的方法，它在一定程度上改善了EMD易模態(tài)混疊的缺陷，但存在端點效應(yīng)抑制的問題。MA通過結(jié)構(gòu)元素探針在信號中的移動來提取信號的沖擊特征，具有抑制脈沖干擾能力較強的優(yōu)點，但存在盲目選擇結(jié)構(gòu)元素和過于依賴相關(guān)先驗知識的缺陷。GOUTSIAS等[4]針對以上問題提出了形態(tài)小波(Morphological Wavelet，MW)分解，將形態(tài)濾波方法引入小波多分辨率分解中，在信號處理上具有良好的降噪性能，但信號分解時“隔二抽取”方式會造成逐層信息減半的問題。ZHANG等[5]提出了一種形態(tài)非抽樣小波(Morphological Un-Decimated Wavelet，MUDW)分解，它省去分解過程下抽樣和重構(gòu)過程上抽樣，解決了MW在信號分解時逐層信息減半的問題；但MUDW分解以信號最高分解層的近似信號作為預(yù)處理結(jié)果，部分故障信息仍會丟失[6]。

基于此，筆者提出一種MUDW分解和峭度融合處理齒輪故障振動信號的方法。首先，應(yīng)用MUDW對信號進(jìn)行分解，在此基礎(chǔ)上對MUDW初始參數(shù)進(jìn)行分析，并利用網(wǎng)格搜索法對參數(shù)組合進(jìn)行優(yōu)選，以避免主觀經(jīng)驗對信號預(yù)處理效果的影響；然后，以峭度作為沖擊信號的評估指標(biāo)來確定包絡(luò)分析的最優(yōu)頻帶，并對各分解層近似信號進(jìn)行加權(quán)融合，以解決傳統(tǒng)MUDW分解方法存在的部分信息遺漏問題，更好地提高特征信息比重；最后，利用仿真信號和齒輪故障的實測振動信號對該方法進(jìn)行驗證。

1 基于MUDW的信號預(yù)處理

1.1 MUDW分解的基礎(chǔ)運算

假設(shè)集合Vi為第i層信號空間，Wi為第i層細(xì)節(jié)空間，T(·)為形態(tài)算子，ZHANG等[5]給出了傳統(tǒng)MUDW分解方法的一般框架：

(1)

(2)

(3)

MUDW是基于非抽樣算法和形態(tài)濾波算子構(gòu)造的，其關(guān)鍵的運算在于形態(tài)算子的選擇。T(·)通常為腐蝕、膨脹、開運算、閉運算4種形態(tài)運算的任意組合形式，常用算子有形態(tài)梯度算子、形態(tài)差值算子和混合算子[10-11]。因齒輪故障信號中包含著不同尺度的形態(tài)特征信息，僅用單尺度結(jié)構(gòu)元素對信號進(jìn)行MUDW分解不能有效提取多尺度特征信息，需進(jìn)行多尺度形態(tài)差值濾波算子的MUDW變換，同時提取信號中的正負(fù)脈沖。所以，MUDW的基礎(chǔ)運算可描述為[12]

(4)

f(xi)°(i+1)g]，

(5)

(6)

式中：f(xi)為原始信號；“°”和“·”分別表示開運算和閉運算；(i+1)g表示對結(jié)構(gòu)元素g進(jìn)行i次膨脹操作。

1.2 MUDW參數(shù)優(yōu)選

與線性小波分解類似，MUDW的分解層數(shù)N和分解元素長度L也是需要注意的問題，通常是根據(jù)待提取的形態(tài)特征來選擇，但缺乏科學(xué)方法指導(dǎo)，可考慮刻畫一種指標(biāo)來衡量能達(dá)到信號預(yù)處理良好效果的最優(yōu)參數(shù)組合。引用小波能量熵的概念：將小波分析和信息熵結(jié)合，反映信號的能量分布信息。將小波能譜[13]表征為E，即信號函數(shù)x(t)在尺度為a時的能量值。x(t)經(jīng)分解后在N層上的小波能量為

E=(E1,E2,…,EN)

，

定義相應(yīng)的小波能量熵

(7)

定義相對小波能量熵Hw為衡量對特征信息利用程度的指標(biāo)：

(8)

Hw越小，信號成分越單純，特征信息越明顯，信號處理效果越好[14]。

在此基礎(chǔ)上采用網(wǎng)格搜索法[15]對(N，L)參數(shù)組合進(jìn)行優(yōu)選，具體步驟如下：

1) 初始化網(wǎng)格搜索中的搜索范圍和搜索步長。設(shè)置N的搜索范圍，搜索步長為1；L的最大值不超過信號的脈沖周期，結(jié)合參數(shù)N可確定L的搜索范圍，搜索步長為1；在N和L的坐標(biāo)系上構(gòu)造一個二維網(wǎng)格。

2) 計算所有網(wǎng)格點的目標(biāo)函數(shù)值。即相對小波能量熵Hw，根據(jù)最小Hw評價當(dāng)前參數(shù)。

3) 將當(dāng)前參數(shù)存放于記憶器中，若滿足Hw最小，則搜索結(jié)束，將Hw最小的那組(N，L)作為最優(yōu)解。

2 基于峭度的故障信號融合

2.1 MUDW融合指標(biāo)構(gòu)建

通過分析MUDW分解算法可知：分解后的每一層近似信號都不同程度地包含著沖擊信息。為了充分利用這些信息，以確保融合過程中能有效提高故障信息含量并減少干擾噪聲，需要選擇能敏感反映沖擊信息的融合指標(biāo)。

反映沖擊程度的常用指標(biāo)有均方、脈沖、峭度和特征頻率幅值。其中，峭度是時域統(tǒng)計指標(biāo)中無量綱的特征分析值，在機械設(shè)備盲信號處理過程中是非高斯性的自然度量指標(biāo)，它反映了信號概率密度函數(shù)分布與高斯分布的偏離程度。峭度對沖擊信息比較敏感，能夠有效判斷故障振動信號的沖擊程度，故障越嚴(yán)重(沖擊信號越明顯)，偏離程度越高，峭度值越大。因此，為了突出信號中的沖擊特征，有效篩選出對沖擊貢獻(xiàn)大的信號頻段，筆者采用峭度[16]作為評價指標(biāo)對各頻段進(jìn)行加權(quán)融合。

計算所有頻帶包絡(luò)信號的峭度值，最大峭度值對應(yīng)的頻帶則被確定為包含沖擊信號強度最大的頻帶。峭度較好地反映了振動信號中沖擊能量的大小，是歸一化的4階中心矩，其計算公式為[16-17]

(9)

分解后不同頻段信號峭度值反映了本段近似信號的沖擊程度，峭度值越大，信號中沖擊脈沖信號所占的比重越大[18]，以此衡量含故障信息的貢獻(xiàn)度，計算MUDW融合權(quán)值。

2.2 MUDW融合權(quán)值的計算及信號重構(gòu)

假設(shè)MUDW各分解層近似信號為xi(i=1,2,…,N)，其對應(yīng)的峭度為Ki，則其融合權(quán)值為

ki=Ki∑Ki。

(10)

由Ki的特性分析可知：Ki值越大，其對應(yīng)的近似信號對特征的貢獻(xiàn)度也相對越大，對應(yīng)的ki也越大。因此，加權(quán)融合后重構(gòu)的信號

(11)

因重構(gòu)信號包含了各分解層特征信息，因此較融合前信號內(nèi)包含的特征信息量得到了有效改善，提高了信噪比。

MUDW分解后，利用峭度指標(biāo)進(jìn)行加權(quán)融合的振動信號預(yù)處理流程如圖2所示。

3 仿真信號分析

為驗證方法的有效性，采用仿真信號模擬齒輪局部異常時的振動信號。設(shè)采樣頻率fs=2 048 Hz，采樣時間t=0.5 s，采樣長度為1 024的仿真信號

x(t)=x1(t)+x2(t)+n(t)，

(12)

式中：x1(t)模擬故障產(chǎn)生的轉(zhuǎn)頻沖擊信號，其沖擊頻率f0=32 Hz；x2(t)模擬正常的嚙合振動信號及其他傳遞到測點的信號，其齒數(shù)為20，將正常的嚙合振動信號設(shè)為sin(2π×640t)，其他部件的振動信號設(shè)為0.5(sin(2π×15t)+sin(2π×100t))；n(t)模擬白噪聲。

仿真信號x(t)的時域、頻域圖分別如圖3、4所示。由圖4可見：信號故障特征頻率32 Hz基本淹沒在噪聲中，故障特征頻率和噪聲混雜在一起。因此，需利用本文所提出的方法對故障信號進(jìn)行預(yù)處理。

3.1 基于網(wǎng)絡(luò)搜索法的參數(shù)(N，L)優(yōu)選

首先,確定N的搜索范圍。因由MUDW分解得到的信號最小帶寬只有大于特征頻率的3倍[18]，才能確保其包含的沖擊信號足夠長，且分解重構(gòu)后信號與原信號長度保持一致，由此確定其最大分解層數(shù)不會過大。由于MUDW的離散性，其分解層數(shù)N與信號長度相關(guān)，即信號采樣長度為2N。分解初始層數(shù)由2開始，在分解到最大階次時伸縮的小波母函數(shù)長度不應(yīng)大于待分析信號的長度，因仿真信號采樣長度為1 024，所以N最大值不能超過10，由此設(shè)置N的搜索范圍為[2，10]。

然后，確定L的搜索范圍。每個參數(shù)N對應(yīng)L的一個搜索范圍，根據(jù)前面分析并參考文獻(xiàn)[6]，L搜索范圍為[Nmin+1,?(fs/f0+N-1)/N」]，其中?·」表示向下取整運算，且fs=2 048 Hz，f0=32 Hz，在坐標(biāo)系上構(gòu)造(N，L)參數(shù)組合的二維網(wǎng)格。

最后，在(N，L)二維網(wǎng)格內(nèi)，根據(jù)式(8)利用最小Hw評價當(dāng)前參數(shù)，將當(dāng)前最優(yōu)參數(shù)存放于記憶器中。

在(N，L)數(shù)值組合的搜索范圍中，N可取2～10的任意值，且都能對應(yīng)得到L的搜索范圍，具體變化情況如表1所示。

表1 N和L搜索范圍

圖5為不同(N，L)參數(shù)組合下Hw的數(shù)值變化情況。由式(8)計算出每個參數(shù)組合下所對應(yīng)的Hw，N=2，3，…，10時對應(yīng)的最小Hw分別為0.051 8、0.041 5、0.041 8、0.035 3、0.036 8、0.040 5、0.040 6、0.042 9、0.043 1。計算得到當(dāng)N=5，L=7時，minHw=0.035 3，說明此時MUDW分解處理效果最好，信號構(gòu)成成分簡單，特征信息成分突出。

3.2 基于MUDW和峭度的振動信號預(yù)處理

當(dāng)N=5，L=7時，設(shè)結(jié)構(gòu)元素g0={0，0，0，0，0，0，0}，由式(4)-(6)對信號進(jìn)行分解。由式(9)計算各層近似信號的峭度分別為：K1=13.301,K2=10.667,K3=11.393,K4=13.805,K5=10.336。根據(jù)式(10)計算相應(yīng)的融合權(quán)值分別為：k1=0.223 5，k2=0.179 3，k3=0.191 5，k4=0.232，k5=0.173 7。根據(jù)式(11)進(jìn)行加權(quán)融合，可得重構(gòu)信號

xFinal= 0.223 5x1+0.179 3x2+0.191 5x3+

0.232x4+0.173 7x5。

圖6為MUDW融合峭度指標(biāo)對信號預(yù)處理后的頻譜圖?？梢钥闯觯号c圖4原始信號相比，能夠清晰地看到特征頻率f0=32 Hz及其2倍頻、3倍頻，說明噪聲及諧波干擾得到了有效抑制。

4 試驗驗證

為驗證本文所提方法的實用性，將該方法應(yīng)用于齒輪裂紋故障振動信號預(yù)處理中。將齒輪齒根裂紋故障加工在輸出軸大齒輪齒根上，進(jìn)行齒輪齒根裂紋故障試驗。齒根裂紋采用線切割以角度α進(jìn)行加工，以深度5 mm進(jìn)行試驗。齒根裂紋加工位置及寬度如圖7所示。

在機械振動及故障模擬平臺上進(jìn)行試驗，如圖8所示，主要設(shè)備有二級平行軸變速箱、4 kW三相電磁調(diào)速電動機和風(fēng)冷磁粉制動器。電機輸出軸與變速箱輸入軸通過聯(lián)軸器聯(lián)接，在聯(lián)軸器的右側(cè)安裝了一個轉(zhuǎn)速轉(zhuǎn)矩傳感器，軸旋轉(zhuǎn)一轉(zhuǎn)產(chǎn)生60個脈沖。將4個3056B4型壓電加速度傳感器安裝在變速箱上，利用數(shù)據(jù)采集板及LabVIEW軟件將采集到的振動信號存入電腦。

該變速箱齒輪的齒數(shù)分別為：低速軸齒輪(故障齒輪)齒數(shù)Z1=81，中間軸大齒輪齒數(shù)Z2=64，中間軸小齒輪齒數(shù)Z3=19，高速軸齒輪齒數(shù)Z4=35。輸入軸(軸3)轉(zhuǎn)速為800 r/min，傳感器位置表示為①、②、③、④，傳感器位置及變速箱內(nèi)部結(jié)構(gòu)如圖9所示。試驗負(fù)載為10 N·m，采樣頻率為20 kHz，信號采樣長度為1 024，計算可得軸3轉(zhuǎn)頻f3=800/60=13.33 Hz。設(shè)輸出軸(軸2)轉(zhuǎn)頻為f2，平行軸齒輪傳動比i1=Z4/Z2=f2/f3，則計算可得f2=7.29 Hz，一級嚙合頻率fm1=f1×Z4=f2×Z2=466.6 Hz，軸2轉(zhuǎn)速為f2×60=437 r/min。設(shè)故障齒輪軸(軸1)轉(zhuǎn)頻為f1，齒輪傳動比i2=Z3/Z1=f1/f2，則計算可得f1=1.71 Hz，二級嚙合頻率fm2=f2×Z3=f1×Z1=138.5 Hz。一般來說，具有局部異常故障的齒輪(如裂紋故障)將以轉(zhuǎn)頻為主要頻率特征，主要分析其嚙合頻率及其邊頻帶。齒輪箱含裂齒時，其振動能量將大幅增加，嚙合頻率及其諧波周圍會產(chǎn)生邊頻帶，邊頻帶寬且高。

以傳感器③ 6 s采集的振動數(shù)據(jù)作為原始樣本，其所采集的故障振動信號的時域和頻域圖(功率譜)分別如圖10、11所示?？梢姡汗收闲盘柎嬖诿黠@的調(diào)制現(xiàn)象，故障特征頻率淹沒在噪聲干擾中，且信號特征不明顯。

采用本文所提方法對齒輪故障信號進(jìn)行預(yù)處理，具體如下：

1) 采用網(wǎng)格搜索法計算分解層數(shù)N和分解元素長度L，根據(jù)采樣長度可得N的搜索范圍為[2，10]，N=2，3，…，10時對應(yīng)的最小Hw分別為0.0573、0.050 8、0.041 3、0.034 8、0.041 7、0.043 1、0.040 1、0.042 1、0.048 2。計算得到當(dāng)N=5，L=4時，minHw=0.034 8，表明此時MUDW預(yù)處理效果最佳，信號構(gòu)成成分簡單，特征信息成分突出。

2)設(shè)置結(jié)構(gòu)元素g0={0，0，0，0}，由式(4)-(6)對信號進(jìn)行分解。由式(9)計算各層信號的峭度分別為：K1=6.027 8，K2=4.226 2，K3=5.131 1，K4=5.036 5，K5=4.028 6。根據(jù)式(10)計算相應(yīng)的融合權(quán)值分別為：k1=0.246 5，k2=0.172 8，k3=0.209 9，k4=0.206，k5=0.164 8，根據(jù)式(11)進(jìn)行加權(quán)融合，可得重構(gòu)信號

xFinal= 0.246 5x1+0.172 8x2+0.209 9x3+

0.206x4+0.164 8x5。

基于MUDW和峭度的故障信號預(yù)處理效果如圖12所示。

由圖12可見：通過MUDW分解和采用峭度指標(biāo)加權(quán)融合對故障信號進(jìn)行預(yù)處理后，能較清晰地看到輸入軸轉(zhuǎn)頻f3=13.33 Hz、故障齒輪的嚙合頻率fm2=138 Hz及其倍頻，以及邊頻帶中含故障特征頻率(裂齒轉(zhuǎn)頻)的倍頻，原振動信號所存在的噪聲和諧波干擾現(xiàn)象得到了有效抑制。

由于瞬態(tài)沖擊能量大，激勵起齒輪固有頻率，產(chǎn)生了齒輪轉(zhuǎn)頻調(diào)制現(xiàn)象，出現(xiàn)了以齒輪嚙合頻率為中心頻率、以裂齒所在軸的轉(zhuǎn)頻及其高次諧波為調(diào)制頻率的調(diào)制邊頻帶。對比正常齒輪振動信號并對信號進(jìn)行預(yù)處理，效果如圖13所示?？梢钥闯觯赫｝X輪的振動信號經(jīng)預(yù)處理后可以清晰地觀察到輸入軸轉(zhuǎn)頻，由于齒輪未發(fā)生故障，因此其他轉(zhuǎn)頻及嚙合頻率并沒有明顯的沖擊信號。對比圖12可知：故障齒輪的嚙合頻率沖擊信號較明顯(138 Hz)，可以較好地判斷出齒輪發(fā)生故障的位置。

為了進(jìn)一步驗證加權(quán)重構(gòu)方法的有效性，與傳統(tǒng)MUDW信號預(yù)處理方法進(jìn)行對比，效果如圖14所示。可見：傳統(tǒng)MUDW信號預(yù)處理方法受噪聲干擾，且部分特征頻率不明顯或淹沒在噪聲中，同時遺漏了部分特征信息(如故障軸嚙合頻率2倍頻等)，導(dǎo)致整體含特征信息比重相對較小；而基于MUDW和峭度的信號預(yù)處理方法因采用融合處理及參數(shù)優(yōu)選，有效利用了包含在各分解層的特征信息并利用權(quán)重突出了沖擊信號，比傳統(tǒng)MUDW信號預(yù)處理方法特征信息比重大且故障特征明顯，處理效果更清晰有效。

5 結(jié)論

針對齒輪故障振動信號不平穩(wěn)、非線性強，導(dǎo)致處理效果不理想的問題，筆者在MUDW分解運算的基礎(chǔ)框架下，提出了一種基于MUDW和峭度融合處理齒輪故障振動信號的方法，通過對信號各分解層的加權(quán)融合，有效提高了特征信息的比重，同時對其中的參數(shù)進(jìn)行了優(yōu)選，避免了主觀經(jīng)驗對信號預(yù)處理效果的影響，最后利用仿真信號及齒輪實測信號驗證了方法的有效性和實用性。結(jié)果表明：該方法對齒輪故障振動信號的預(yù)處理效果理想，能有效提取故障信號特征，為后續(xù)的故障診斷研究奠定了基礎(chǔ)。

參考文獻(xiàn)：

[1] 蘇勛文,王少萍.齒輪故障演化的非線性動力學(xué)機理與實驗研究[J].機械傳動,2015,39(2):19-24.

[2] SARAVANAN N,RAMACHANDRA K I.Incipient gearbox fault diagnosis using Discrete Wavelet Transform (DWT) for feature extraction and classification using Artificial Neural Network (ANN) [J].Expert systems with applications,2010,37(6):4168-4181.

[3] FREI M G,OSORIO I.Intrinsic time-scale decomposition:time-frequency-energy analysis and real-time filtering of non-stationary signals[C]∥Proceedings of the Royal Society A:Mathematical,Physical & Engineering Sciences.England:The Royal Society,2007,463(2078):321-342.

[4] GOUTSIAS J,HEIJMANS H J A M.Nonlinear multi resolution signal decomposition schemes,part 1:morphological pyramids [J].IEEE transaction on image processing,2000,9(11):1862-1876.

[5] ZHANG J F,SMITH J S,WU Q H.Morphological un-decimated wavelet decomposition for fault location on power transmission lines [J].IEEE transactions on circuits and systems,2006,53(6):1395-1402.

[6] 孫健,李洪儒,王衛(wèi)國,等.一種基于WMUWD的液壓泵振動信號預(yù)處理方法[J].振動與沖擊,2015,34(21):93-99.

[7] 章立軍,陽建宏,徐金梧,等.形態(tài)非抽樣小波及其在沖擊信號特征提取中的應(yīng)用[J].振動與沖擊,2007,26(10):56-59.

[8] COUSTY J,NAJMAN L,DIAS F,et al.Morphological filtering on graphs[J].Computer vision and image understanding,2013,117(4):370-385.

[9] 黃兵峰,沈路,周曉軍,等.基于形態(tài)非抽樣小波分解的滾動軸承故障特征提取[J].農(nóng)業(yè)機械學(xué)報,2010,41(2):204-207.

[10] 沈路,周曉軍,張杰.形態(tài)非抽樣小波與灰色關(guān)聯(lián)度的軸承故障診斷[J].電子機械工程,2012,28(5):60-64.

[11] HE W,JIANG Z N,QIN Q.A joint adaptive wavelet filter and morphological signal processing method for weak mechanical impulse extraction[J].Journal of mechanical science and technology,2011,24(8):1709-1716.

[12] 孫健,李洪儒,王衛(wèi)國,等.基于形態(tài)非抽樣融合與DCT高階奇異熵的液壓泵退化特征提取[J].振動與沖擊,2015,34(22):54-61.

[13] 申弢,黃樹紅,楊叔子.旋轉(zhuǎn)機械振動信號的信息熵特征[J].機械工程學(xué)報,2001,37(6):94-98.

[14] 艾延廷,付琪,田晶,等.基于融合信息熵距的轉(zhuǎn)子裂紋-碰摩耦合故障診斷方法[J].航空動力學(xué)報,2013,28(10):2161-2166.

[15] 呂中亮,湯寶平,周憶,等.基于網(wǎng)格搜索法優(yōu)化最大相關(guān)峭度反卷積的滾動軸承早期故障診斷方法[J].振動與沖擊,2016,35(15):29-34.

[16] 郭慶豐,王成棟,劉佩森.時域指標(biāo)和峭度分析法在滾動軸承故障診斷中的應(yīng)用[J].機械傳動,2016,40(11):172-175.

[17] WANG D,TSE P W,TSUI K L.An enhanced kurtogram method for fault diagnosis of rolling element bearings[J].Mechanical systems and signal processing,2013,35:176-199.

[18] ZHANG X H,KANG J S.Rolling element bearings fault diagnosis based on correlated kurtosis kurtogram [J].Journal of vibroengineering,2015,17(6):3023-3034.