基于LCD與雙譜分析的齒輪故障特征提取

仝 蕊， 康建設(shè)， 李寶晨， 鐘 文

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

振動(dòng)信號(hào)處理和特征提取是齒輪故障診斷及預(yù)測(cè)的核心問題，有效的特征提取方法對(duì)提高齒輪箱故障預(yù)測(cè)準(zhǔn)確程度具有重要意義[1]。齒輪振動(dòng)信號(hào)通常表現(xiàn)出復(fù)雜的時(shí)變特性，尤其是故障信號(hào)具有非平穩(wěn)和非高斯分布的特征。譜分析方法是處理非高斯信號(hào)的有效方法，能夠有效提取振動(dòng)信號(hào)特征，被廣泛應(yīng)用于機(jī)械設(shè)備的信號(hào)處理中。

目前，常用的譜分析方法主要有功率譜分析、倒譜分析、奇異譜分析和高階譜分析等方法。其中：功率譜分析方法采用時(shí)域統(tǒng)計(jì)指標(biāo)、功率譜密度分析信號(hào)的頻域結(jié)構(gòu)，其對(duì)非線性較強(qiáng)的信號(hào)處理效果不佳[2]；倒譜分析方法是在對(duì)功率譜取對(duì)數(shù)的基礎(chǔ)上進(jìn)行譜分析，在放大特征分量的同時(shí)也放大了噪聲和其他分量信息，因此影響了特征提取效果[3]；奇異譜分析方法是基于二階統(tǒng)計(jì)的協(xié)方差矩陣進(jìn)行譜分析，但分析非線性信號(hào)效果不理想[4]；雙譜分析方法是一種典型的高階譜分析方法，它既可以消除非線性信號(hào)中高斯干擾成分，同時(shí)又保留了其中的非高斯成分，適于處理具有周期性和調(diào)制性的齒輪振動(dòng)信號(hào)[5]。但雙譜分析方法從理論上僅能抑制高斯噪聲，而對(duì)非高斯噪聲無能為力，且直接對(duì)信號(hào)整體進(jìn)行分析，無法剔除諧波及其他不必要分量信號(hào)的干擾。

通過局部特征尺度分解(Local Characteristic-scale Decomposition，LCD)[6]可提取雙譜分析所需的包含敏感信息的信號(hào)分量，其在迭代速度和分解能力等方面具有明顯優(yōu)勢(shì)，但對(duì)內(nèi)稟尺度分量(Intrinsic Scale Component,ISC)的篩選并無一套系統(tǒng)科學(xué)的方法。為彌補(bǔ)上述缺陷，筆者提出一種利用LCD分解和雙譜分析相結(jié)合的故障特征提取方法,其先對(duì)LCD算法進(jìn)行優(yōu)化，篩選有用分量，再結(jié)合雙譜分析計(jì)算出雙譜熵(Bispectral Entropy,BIE)作為故障特征量，最后通過齒輪預(yù)置故障試驗(yàn)對(duì)該方法進(jìn)行驗(yàn)證。

1 局部特征尺度分解的優(yōu)化

1.1 LCD算法

利用LCD對(duì)任意一個(gè)復(fù)雜的振動(dòng)信號(hào)x(t)進(jìn)行分解[7]，分解為有限個(gè)不同ISC分量之和，且任意2個(gè)ISC分量相互獨(dú)立，其計(jì)算步驟如下：

1)選擇合適參數(shù)構(gòu)造高頻諧波信號(hào)H(t)，進(jìn)行并疊加：

H(t)=Ucos(2πft)，

(1)

x′(t)=x(t)+H(t)。

(2)

式中：U為幅值，常取U=max(x(t))；f為頻率，常取f=fs/2[7],其中fs為采樣頻率。

2) 確定x′(t)的極值點(diǎn)Xk′(k′=1,2,…,M)及對(duì)應(yīng)時(shí)刻τk′，其中M為極值點(diǎn)數(shù)。

3) 計(jì)算Xk′的函數(shù)值

(3)

在任意2個(gè)極值點(diǎn)之間對(duì)x′(t)進(jìn)行線性變換,可得

Lk′=aAk′+(1-a)Xk′,k′=2,3,…,M-1。

(4)

式中：a=0.5[6]。通過延拓獲得左右端點(diǎn)值，具體方法為[7]：在左端增加一個(gè)極值點(diǎn)進(jìn)行左延拓，即令X0=X2，τ0=2τ1-τ2，由式(4)可得到L1；同理，在右端點(diǎn)增加一個(gè)極值點(diǎn)進(jìn)行右延拓，即令XM+1=XM-1，τM+1=2τM-τM-1，由式(4)可得到LM。利用式(4)計(jì)算A1和AM。

4) 利用三次樣條函數(shù)[7]對(duì)L1，L2，…，LM進(jìn)行擬合，對(duì)得到的基線B1(t)分離

h1(t)=x(t)-B1(t),

(5)

判斷h1(t)是否滿足ISC分量條件[8]：若滿足，則輸出ISC1=h1(t)；否則，以h1(t)作為原始數(shù)據(jù)，重復(fù)步驟1)-4)，直到得到ISC分量。

5) 從原始數(shù)據(jù)中分離ISC1分量，將剩余信號(hào)作為新的原始數(shù)據(jù)u1(t)，重復(fù)步驟2)-4)，直到提取ISC分量n′次后滿足un′(t)為單調(diào)函數(shù)或?yàn)槌：瘮?shù)為止，最終得到ISC1，ISC2，…，ISCn′及趨勢(shì)項(xiàng)un′(t)。

1.2 ISC分量的篩選

分析振動(dòng)信號(hào)的LCD分解過程發(fā)現(xiàn)：原始信號(hào)所包含的故障信息被細(xì)化到各個(gè)分量信號(hào)中，產(chǎn)生了包含關(guān)鍵故障信息的ISC分量。分量信號(hào)不僅含有較多的故障特征信息，還含有較多的干擾信息。為了減少干擾分量的影響，對(duì)LCD算法進(jìn)行優(yōu)化，并選擇合適的方法對(duì)分解出的ISC分量進(jìn)行篩選，具體方法如下：

1) ISC分量個(gè)數(shù)m的估計(jì)。將LCD分解所得的ISC分量和趨勢(shì)項(xiàng)un′(t)與分解前的原始振動(dòng)信號(hào)x(t)進(jìn)行重組，構(gòu)成新信號(hào)xob(t)={x(t),ISC1(t),…,ISCn′(t),un′(t)} 。構(gòu)建xob(t)的相關(guān)矩陣R，并進(jìn)行奇異值分解[8]：

(6)

式中：u′為xob(t)的維數(shù)；Λd為包含d個(gè)主特征值的主特征值矩陣；Vd為奇異值矩陣；Λu′-d為包含u′-d個(gè)噪聲特征值的噪聲特征值矩陣。Λu′-d的噪聲特征值理論上與噪聲功率相等，且取值非常小，故可根據(jù)矩陣R中最小特征值的個(gè)數(shù)直接確定主特征值數(shù)，估算出m值。但實(shí)際應(yīng)用中，噪聲大小并不能控制和估算，尤其是在強(qiáng)噪聲背景下無法有效確定出主特征值與噪聲特征值間的閾值。因此，本文以R的奇異值分解結(jié)果中的非零特征值為參考，結(jié)合貝葉斯信息準(zhǔn)則(Bayesian Information Criterion,BIC)[7,9]

(7)

2) 有用ISC分量的確定。峭度可有效體現(xiàn)信號(hào)中的沖擊特征，為了提取有用的ISC分量，可將計(jì)算的ISC分量和原始信號(hào)的峭度時(shí)間序列互相關(guān)系數(shù)作為篩選有用ISC分量的評(píng)判指標(biāo)。對(duì)采樣長(zhǎng)度為N的原始信號(hào)x(t)，因需要足夠多的樣本計(jì)算才有意義，故選擇計(jì)算峭度時(shí)的長(zhǎng)度a′=30。以xi(i=1,2,…,N-a′+1)為起點(diǎn)，從原始振動(dòng)信號(hào)中依次向后截取長(zhǎng)度為a′的數(shù)據(jù)序列[10]，計(jì)算其峭度值，從而得到長(zhǎng)度為N-a′+1的峭度時(shí)間序列K(i)：

(8)

,

(9)

式中：ρ為互相關(guān)函數(shù)取得最大值的間隔數(shù)。

互相關(guān)系數(shù)的大小反映了不同ISC分量與原始信號(hào)之間的相關(guān)性。信號(hào)經(jīng)過LCD分解后，各有用分量與原始信號(hào)的相關(guān)性約等于該分量的自相關(guān)；而非有用分量與原始信號(hào)的相關(guān)性很小[12]?；ハ嚓P(guān)系數(shù)越大，包含的敏感信息越多；反之，干擾成分越多。因此，結(jié)合前面估計(jì)的m個(gè)ISC分量，選取互相關(guān)系數(shù)最大的前m個(gè)ISC分量作為有用分量，即為后續(xù)雙譜分析融合的信號(hào)分量。

2 雙譜分析

2.1 雙譜估計(jì)

雙譜定義為三階累積量的二維離散傅里葉變換，用Bx(ω1,ω2)表示，其表達(dá)式可參考文獻(xiàn)[13]。由于高斯平穩(wěn)隨機(jī)信號(hào)的三階累積量為0，其雙譜值也為0，因此理論上雙譜分析能消除信號(hào)中的高斯噪聲干擾，同時(shí)將信號(hào)中的非高斯成分保留。由于所采集的齒輪振動(dòng)信號(hào)為一定長(zhǎng)度的離散序列，無法準(zhǔn)確計(jì)算雙譜值，因此可采用基于傅里葉變換的非參數(shù)化估計(jì)方法進(jìn)行雙譜估計(jì)。具體方法如下[14]：

1) 將采樣長(zhǎng)度為N的輸入信號(hào)x(t)分成X段，每段W個(gè)數(shù)據(jù)點(diǎn)，并將每段數(shù)據(jù)去均值處理。

2) 計(jì)算第j(j=1,2,…,X)段數(shù)據(jù)xj(t)(t=0,1,…,W-1)的離散傅里葉變換系數(shù)

(10)

式中：η=0,1,…,W/2，為離散傅里葉變換后各數(shù)據(jù)點(diǎn)對(duì)應(yīng)的編號(hào)。

3) 計(jì)算各段數(shù)據(jù)的三重相關(guān)

Y(j)(η2+X2)Y(j)(-η1-X1-η2-X2),

(11)

式中：Δ0=fs/N0，為所要求的頻率樣本間隔，其中N0為樣本數(shù)；η1+η2≤fs/2，0≤η2≤η1,其中η1、η2分別為任2段信號(hào)數(shù)據(jù)X1、X2離散傅里葉變換后對(duì)應(yīng)編號(hào)；W=(2L1+1)N0。

4) 計(jì)算雙譜估計(jì)值

(12)

2.2 雙譜熵的計(jì)算

2) 計(jì)算定義域頻域內(nèi)所有點(diǎn)的雙譜幅值概率

(13)

3) 計(jì)算雙譜熵值

(14)

由于BIE包含了對(duì)有用ISC分量信息的融合處理，因此其對(duì)頻帶能量的變化更為敏感：故障程度越輕時(shí)，其能量在各頻帶分布較均勻，對(duì)應(yīng)BIE值越大；故障程度越重時(shí)，其能量主要集中在少數(shù)一些特征頻帶上，對(duì)應(yīng)BIE值越小。在此，將經(jīng)LCD分解獲得ISC分量，然后利用雙譜分析得到的BIE作為故障特征量，稱為L(zhǎng)CD-BIE,上述特征提取過程稱為L(zhǎng)CD-BIE特征提取方法。

3 試驗(yàn)分析驗(yàn)證

3.1 齒輪故障特征的提取

為驗(yàn)證LCD-BIE特征提取方法的實(shí)用性，將該方法應(yīng)用于齒輪預(yù)置故障試驗(yàn)的實(shí)測(cè)振動(dòng)信號(hào)中。以齒輪正常狀態(tài)、裂紋故障、斷齒故障為例。試驗(yàn)臺(tái)包括二級(jí)平行軸變速箱、4 kW三相電磁調(diào)速電動(dòng)機(jī)、風(fēng)冷磁粉制動(dòng)器各1個(gè)，還包括1套數(shù)據(jù)采集系統(tǒng)、1臺(tái)計(jì)算機(jī)、4個(gè)3056B4型壓電加速度傳感器、PXI-1031機(jī)箱、數(shù)據(jù)采集卡及LabVIEW軟件。該變速箱齒輪的齒數(shù)分別為：高速軸齒輪35個(gè)齒；中間軸大齒輪64個(gè)齒；中間軸小齒輪19個(gè)齒；低速軸齒輪81個(gè)齒。在齒輪箱殼體安裝傳感器，采集齒輪箱振動(dòng)信號(hào)。試驗(yàn)負(fù)載為10 N·m，轉(zhuǎn)速為800 r/min，采樣頻率為20 kHz，采樣時(shí)間為6 s。

在輸出軸大齒輪齒根上以線切割角度α進(jìn)行裂紋故障加工，深度分別為1、2、5、8 mm。在中間軸大齒輪上進(jìn)行斷齒故障加工，斷齒寬度分別為2、5、10 mm。將采集到的齒輪箱振動(dòng)信號(hào)按正常、裂紋深度由淺到深、斷齒寬度由窄到寬依次標(biāo)記為xF1,xF2,…,xF8，每個(gè)信號(hào)均包含120 000個(gè)數(shù)據(jù)點(diǎn)。以裂紋深度為5 mm的齒輪箱振動(dòng)信號(hào)xF4為例進(jìn)行分析，其步驟如下：

1) 利用帶有高頻諧波的LCD算法對(duì)xF4進(jìn)行分解，得到10個(gè)ISC分量以及1個(gè)趨勢(shì)項(xiàng)u(t)。

2) 對(duì)重組信號(hào)篩選。根據(jù)式(6)、(7)估算融合所需的ISC分量個(gè)數(shù)為2，并利用式(8)、(9)分別計(jì)算每個(gè)ISC分量的峭度時(shí)間序列與原始信號(hào)峭度時(shí)間序列的互相關(guān)系數(shù)，結(jié)果依次為4.27×10-3、4.21×10-2、1.89×10-2、5.87×10-2、2.15×10-2、4.78×10-3、6.03×10-3、4.99×10-2、1.31×10-1、9.58×10-2。

以互相關(guān)系數(shù)絕對(duì)值最大為篩選準(zhǔn)則，所選取的分量為ISC9和ISC10，其頻域如圖1所示?？梢钥闯觯航?jīng)LCD分解和ISC分量篩選后，信號(hào)分量能夠有效地捕捉到故障頻率132 Hz及其倍頻、周邊故障信息豐富的邊頻帶，以及頻率13.33 Hz等這些敏感特征信息，這就可以剔除噪聲分量等干擾，抓取包含斷齒故障頻率的敏感特征信息。

3) 取X=10，對(duì)ISC9、ISC10兩個(gè)分量分段去均值處理并進(jìn)行雙譜估計(jì)。齒輪正常、裂紋故障和斷齒故障的雙譜等高線分布圖分別如圖2-4所示。

由圖3可以看出：1)雙譜幅值的非零成分為信號(hào)中的非高斯噪聲成分，這表明振動(dòng)信號(hào)具有明顯的非高斯特性；2)圖譜聚集性隨裂紋深度增加而不斷增大，表現(xiàn)出了復(fù)雜的耦合特性。這是因?yàn)椋汗收习l(fā)生后，齒輪內(nèi)部運(yùn)行狀態(tài)受到影響，固有頻率因發(fā)生改變而影響周圍零部件的耦合程度，使得振動(dòng)信號(hào)中的非高斯成分大幅增加，其雙譜幅值明顯增大。

由圖4可以看出：1)隨著斷齒寬度的增大，圖譜逐漸向中心聚集，故障頻率位于圖譜峰值處，而邊緣的高頻諧波逐漸減少；2)從幅值上看，振動(dòng)信號(hào)中的非高斯成分隨斷齒寬度的增大而增多，其特征頻率處的雙譜幅值也隨之增大。

由圖2-4可以看出：與故障狀態(tài)相比，正常狀態(tài)下的雙譜幅值明顯較小。這是因?yàn)椋涸谡顟B(tài)下，齒輪箱振動(dòng)信號(hào)的高斯成分占主要部分，而高斯噪聲的雙譜幅值為0。

結(jié)合圖3、4可知：雙譜分布能在一定程度上反映故障變化過程；同類型故障的雙譜分布具有一定相似性[16]，且與不同類型故障相比存在明顯差異，由此可大致區(qū)分不同的故障類型。但僅在平面圖上分析振動(dòng)信號(hào)雙譜幅值的變化難以對(duì)故障程度量化描述。因此，為更好研究故障變化規(guī)律并為剩余壽命預(yù)測(cè)奠定基礎(chǔ)，需要在雙譜分析基礎(chǔ)上進(jìn)行量化處理。

4) 通過式(13)、(14)可得到xF4的LCD-BIE特征量值為6.838。按上述方法分別計(jì)算齒輪箱振動(dòng)信號(hào)的LCD-BIE及其相應(yīng)的特征量值。

3.2 與其他特征提取方法的對(duì)比

為說明本文特征提取方法的有效性，分別計(jì)算xF1，xF2，…，xF8經(jīng)LCD后各階段的能譜熵(LCD-spectrum Entropy,LE)[17]，及經(jīng)LCD分解但未進(jìn)行分量篩選與雙譜分析處理后的雙譜熵BIE1，結(jié)果見表1?？梢钥闯觯弘S齒輪故想障程度的加劇，各特征量熵值大體呈下降趨勢(shì)；LCD-BIE中敏感信息的能量較為集中，其熵值要小于LE和BIE1，說明其特征表征能力較好。

表1 不同故障狀態(tài)的特征量值

對(duì)表1中各特征量進(jìn)行歸一化處理，繪制成曲線，如圖5所示。圖5描述了歸一化后不同故障下各特征量的變化過程。一般認(rèn)為一個(gè)好的特征曲線應(yīng)是隨著故障程度的加深單調(diào)遞增或單調(diào)遞減。

通過對(duì)比圖5中的3條曲線可知：

1) LE曲線對(duì)裂紋故障由1 mm加深至2 mm過程反映較敏感，但對(duì)其他故障狀態(tài)表征效果不理想。這是因?yàn)椋何春Y選過的ISC分量中仍存在其他分量的干擾；能譜熵不具備雙譜熵消除高斯噪聲的優(yōu)勢(shì)。

2) BIE1曲線對(duì)整個(gè)故障變化過程具有相對(duì)較好的表征能力，但由于分量融合前未做篩選處理，導(dǎo)致所提取故障特征包含成分復(fù)雜，整體對(duì)故障的敏感度不夠高，尤其體現(xiàn)在特征曲線單調(diào)遞減趨緩。

3) LCD-BIE曲線避免了干擾成分和消除了高斯噪聲影響，其熵值最小，且能提取更完整的特征信息，對(duì)故障狀態(tài)具有良好的表征能力。

3.3 故障特征值的檢驗(yàn)

故障狀態(tài)識(shí)別是評(píng)價(jià)故障特征表征能力的重要方面。為了檢驗(yàn)本文故障特征提取方法的有效性，采用訓(xùn)練好的支持向量機(jī)(Support Vector Machine,SVM)定量評(píng)價(jià)所提取的特征量對(duì)不同故障的識(shí)別能力[18]。將試驗(yàn)中所采集的F2，F(xiàn)3，…，F(xiàn)8狀態(tài)下各信號(hào)的2 400個(gè)數(shù)據(jù)點(diǎn)作為一組樣本，各狀態(tài)下各取50組樣本，共350組數(shù)據(jù)樣本，計(jì)算LE、BIE1、LCD-BIE特征量值。所有的特征量不同特征量對(duì)歸一化到[0，1]，利用訓(xùn)練好的SVM進(jìn)行分類識(shí)別，計(jì)算2種故障狀態(tài)的識(shí)別率，結(jié)果如表2、3所示。

表2 不同特征量對(duì)裂紋故障的識(shí)別率

表3 不同特征量對(duì)斷齒故障的識(shí)別率

由表2、3可以看出：LE識(shí)別率較低；BIE1對(duì)斷齒故障的識(shí)別率高，但對(duì)裂紋故障的識(shí)別能力較差；LCD-BIE整體識(shí)別率最高，僅對(duì)個(gè)別故障狀態(tài)的識(shí)別存在一定局限性。因此，LCD-BIE是有效的特征量，能夠相對(duì)準(zhǔn)確地識(shí)別設(shè)備退化狀態(tài)。分析其原因?yàn)椋篖CD-BIE經(jīng)過信號(hào)篩選可以更有效提取關(guān)鍵故障信息；而LE包含較多噪聲成分,穩(wěn)定性相對(duì)較差；BIE1比LE敏感度較高，但未篩選的信號(hào)成分依然復(fù)雜，因此穩(wěn)定性差于LCD-BIE。

4 結(jié)論

筆者首先利用LCD分解去除信號(hào)高頻噪聲，再采用相關(guān)方法篩選出有用ISC分量，最后對(duì)融合后信號(hào)進(jìn)行雙譜分析，能夠有效地提取出故障特征率。主要結(jié)論如下：

1) 齒輪故障信號(hào)具有非高斯特性，而噪聲和其他干擾分量往往是高斯信號(hào)，利用雙譜分析能夠抑制振動(dòng)信號(hào)中的高斯噪聲，準(zhǔn)確地提取信號(hào)中的故障特征信息。

2) 采用BIC準(zhǔn)則和峭度時(shí)間序列互相關(guān)系數(shù)相結(jié)合的方法對(duì)ISC分量進(jìn)行選取，能夠有效地去除干擾分量，抓取有用信息，提高特征的敏感度，從而進(jìn)一步改善對(duì)故障狀態(tài)的表征能力。

3) LCD-BIE故障特征提取方法能有效提取正常狀態(tài)下的齒輪故障特征，提高了對(duì)故障的識(shí)別能力。但在變工況條件下，振動(dòng)信號(hào)將更加復(fù)雜，故障信號(hào)與噪聲之間的關(guān)系有待于研究，算法的適用性還要進(jìn)一步驗(yàn)證。