仝 蕊， 康建設， 李寶晨， 張星輝

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

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

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

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

1 基于MUDW的信號預處理

1.1 MUDW分解的基礎運算

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

(1)

(2)

(3)

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

(4)

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

(5)

(6)

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

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

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

E=(E1,E2,…,EN)

，

定義相應的小波能量熵

(7)

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

(8)

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

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

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

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

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

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

2.1 MUDW融合指標構建

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

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

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

(9)

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

2.2 MUDW融合權值的計算及信號重構

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

ki=Ki∑Ki。

(10)

由Ki的特性分析可知：Ki值越大，其對應的近似信號對特征的貢獻度也相對越大，對應的ki也越大。因此，加權融合后重構的信號

(11)

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

MUDW分解后，利用峭度指標進行加權融合的振動信號預處理流程如圖2所示。

3 仿真信號分析

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

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

(12)

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

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

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

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

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

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

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

表1 N和L搜索范圍

圖5為不同(N，L)參數(shù)組合下Hw的數(shù)值變化情況。由式(8)計算出每個參數(shù)組合下所對應的Hw，N=2，3，…，10時對應的最小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。計算得到當N=5，L=7時，minHw=0.035 3，說明此時MUDW分解處理效果最好，信號構成成分簡單，特征信息成分突出。

3.2 基于MUDW和峭度的振動信號預處理

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

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

0.232x4+0.173 7x5。

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

4 試驗驗證

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

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

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

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

采用本文所提方法對齒輪故障信號進行預處理，具體如下：

1) 采用網(wǎng)格搜索法計算分解層數(shù)N和分解元素長度L，根據(jù)采樣長度可得N的搜索范圍為[2，10]，N=2，3，…，10時對應的最小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。計算得到當N=5，L=4時，minHw=0.034 8，表明此時MUDW預處理效果最佳，信號構成成分簡單，特征信息成分突出。

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

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

0.206x4+0.164 8x5。

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

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

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

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

5 結論

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

參考文獻：

[1] 蘇勛文,王少萍.齒輪故障演化的非線性動力學機理與實驗研究[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的液壓泵振動信號預處理方法[J].振動與沖擊,2015,34(21):93-99.

[7] 章立軍,陽建宏,徐金梧,等.形態(tài)非抽樣小波及其在沖擊信號特征提取中的應用[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].農業(yè)機械學報,2010,41(2):204-207.

[10] 沈路,周曉軍,張杰.形態(tài)非抽樣小波與灰色關聯(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] 申弢,黃樹紅,楊叔子.旋轉機械振動信號的信息熵特征[J].機械工程學報,2001,37(6):94-98.

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

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

[16] 郭慶豐,王成棟,劉佩森.時域指標和峭度分析法在滾動軸承故障診斷中的應用[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.