孫斌 李宏坤 劉艾強 張孔亮

（大連理工大學(xué)機械工程學(xué)院）

0 引言

行星齒輪箱在工業(yè)應(yīng)用廣泛，如風(fēng)力發(fā)電，航天航空等[1]。行星齒輪箱作為旋轉(zhuǎn)機械的主要傳動部件，長期工作在沖擊、交變載荷和工況大幅度變化的情況下，其關(guān)鍵零部件如齒輪，極易發(fā)生故障，從而引起巨大的經(jīng)濟損失和惡性事故[2-3]。因此，對行星齒輪箱的齒輪故障診斷具有重大的經(jīng)濟與安全意義。

行星齒輪箱在變轉(zhuǎn)速工況下，所測得的振動信號具有典型的非平穩(wěn)和非線性特點。如果利用傳統(tǒng)的信號處理方法，如快速傅里葉變化，得到的頻譜模糊，無法獲得故障特征頻率，難以識別變轉(zhuǎn)速下的齒輪故障[4]。

時變轉(zhuǎn)速下對行星齒輪箱的故障診斷，階次跟蹤被認為是最有效的工具之一[5]。在階次跟蹤中，原始信號以恒定的角度間隔重新被采樣，原始的時域非平穩(wěn)信號轉(zhuǎn)化為角域的平穩(wěn)信號。階次跟蹤的方法可以消除由于轉(zhuǎn)速的變化帶來的頻譜模糊的現(xiàn)象，在傳統(tǒng)的階次跟蹤中，需要安裝硬件設(shè)備獲得轉(zhuǎn)速信息用于信號的重采樣。但是，由于空間和成本有限，安裝硬件設(shè)備在很多場合都難以實現(xiàn)[6]。

為了克服上述問題，從振動信號的時頻表示中提取到瞬時速度信息是比較常用和節(jié)約成本的方法。但是，時頻聚集性對所要提取的瞬時速度精度影響較大。許多學(xué)者將短時傅里葉變換（STFT）和連續(xù)小波變換（CWT）的傳統(tǒng)時頻分析方法，應(yīng)用于機械故障診斷領(lǐng)域。雖然這些傳統(tǒng)的時頻方法能夠處理非平穩(wěn)信號，但是受海森堡（Heisenberg）不確定性原理的限制，時間和頻率分辨率不能同時達到最佳，同時也對噪聲比較敏感[7]。此外，由于所測得行星齒輪箱的信號是典型的調(diào)制信號，使用傳統(tǒng)的時頻分析會出現(xiàn)能量混疊和時頻脊線模糊，導(dǎo)致難以準(zhǔn)確提取到轉(zhuǎn)頻信息[8]。為了提高時頻聚集性，I Daubechies 等提出了同步壓縮變換（SST）,利用時頻重排的方法，獲得了聚集性良好的時頻表示[9]。Gang Yu等人提出了同步提取變換(SET)，與SST的壓縮方式不同，只保留與時變特征最相關(guān)的時頻信息，大大提高了時頻的能量聚集性[10]。但是上述時頻分析方法存在對噪聲敏感，不能很好的表征非線性特征的多分量信號的問題。本文運用一種新型的時頻分析方法，廣義線性調(diào)頻小波變換（GLCT）,可以克服現(xiàn)有時頻方法的一些缺陷，并且對噪聲不敏感[11]。

從時頻表示中，提取到轉(zhuǎn)頻曲線是無轉(zhuǎn)速階次分析中不可或缺的一部分。Zhao等使用峰值搜索算法從時頻圖中提取瞬時頻率，實現(xiàn)了變轉(zhuǎn)速下軸承的故障診斷[12]。陳劍等人提出一種改進的Crazy Cliber算法，從STFT 時頻脊線提取脊線，并成功實現(xiàn)了變轉(zhuǎn)速下的軸承故障特征提取[13]。但是，峰值搜索算法和Crazy Cliber算法容易陷入局部最優(yōu)，考慮不到整個信號的輪廓，如果選擇不合適的起點，提取到的瞬時頻率將會是錯誤的，而且這種方法在低信噪比中是無效的[14]。

針對上述存在時頻聚集性不高和提取到的瞬時頻率準(zhǔn)確性低的問題，本文提出了一種基于GLCT和FPO相結(jié)合的方法，用于診斷在時變轉(zhuǎn)速下行星齒輪箱齒輪的故障。采用GLCT對振動信號進行時頻表示，該時頻方法比傳統(tǒng)的時頻分析方法如STFT，時頻聚集性高，并且對噪聲不敏感；再采用FPO 從時頻矩陣中提取瞬時速度信息，該提取算法克服了峰值搜索算法的噪聲干擾和局部最優(yōu)[15]，然后利用所提取到的瞬時速度信息進行信號重采樣，得到角域平穩(wěn)信號，再進行包絡(luò)階次分析，得到故障特征階次，達到故障診斷的目的。仿真和實驗的振動信號證明了該方法在變轉(zhuǎn)速運行下行星齒輪箱齒輪故障診斷的有效性。

1 基本原理介紹

1.1 廣義線性調(diào)頻小波變換(GLCT)

標(biāo)準(zhǔn)的短時傅里葉變換（STFT）如式(1)：

將引進一個解調(diào)因子，去消除調(diào)制成分的影響?？紤]到調(diào)制成分是隨著時間變化的，所以引入的解調(diào)因子也必須是時變的。引入時變的解調(diào)因子s(t) 的STFT表達式為：

s(t)的STFT幅值用瞬時頻率φ(t')可以表達為：

從上式可以看出，如果解調(diào)因子和分析信號的調(diào)制成分相一致，則瞬時頻率φ(t')的STFT 幅值達到最大，然后時頻能量在瞬時頻率處具有很高的聚集性。但是在大多數(shù)情況下，一個信號的IF 特征不能預(yù)先得知，則解調(diào)因子e-ic(t')(u-t')2/2難以準(zhǔn)確確定，特別是含有復(fù)雜的多成分信號中，更加困難。為了解決上述問題，用離散的解調(diào)因子無限去逼近最優(yōu)的解調(diào)因子e-ic(t')(u-t')2/2，則用離散化的解調(diào)因子表示STFT 的公式如式(4)：

如果離散的解調(diào)因子與信號的調(diào)制成分相近時，則在時頻表示中，瞬時頻率的時頻點(t',ω)有很高的能量聚集，則|S(t',ω,c) |幅值達到最大值。從而能夠得到最佳的參數(shù)c為：

則時頻表達式為：

頻譜的定義式為：

為了得到最優(yōu)的解調(diào)因子e-ic(t')(u-t')2/2，引入?yún)?shù)α，通過旋轉(zhuǎn)arctan(-c)，對時頻結(jié)果產(chǎn)生影響，則時頻平面旋轉(zhuǎn)為：

則上述的時頻表達為：

因為參數(shù)α有N個值，則時頻平面能均分N+1 個部分：

由上述可知，GLCT比STFT多引入了一個參數(shù)N,如果N=1時，GLCT即為STFT。

1.2 快速路徑最優(yōu)算法(FPO)

經(jīng)過上述GLCT 得到的時頻結(jié)果為GS(t',ω)，計算時頻表示中每個時刻幅值的局部極大值，然后得到該時刻所有極大值所對應(yīng)的頻率值為：

其中，tn(n=1,2,...,N)表示時頻表示的時刻點；Np(tn)表示時刻tn的局部極大值的個數(shù)；fm(tn)是時刻tn處第m個極值所對應(yīng)的頻率值，其幅值的大小記為SPm(tn)。

采用FPO 確定每個時刻點處應(yīng)該被提取出來的幅值極大值，極大值所對應(yīng)的頻率連接起來即為目標(biāo)脊線，具體算法如下：

其中，mc(tn)為時刻tn應(yīng)該被提取的幅值極大值；[F]表示優(yōu)化方程，該方程的解就是目標(biāo)脊線。優(yōu)化方程為：

其中，Δfd是fd的導(dǎo)數(shù)；percp[f(t)]表示f(t) 的第p分位數(shù)。則優(yōu)化方程的解為：

其中，n=1,2,...,N,m=1,2,...,Np(tn)，q(m,tn)確定每個時刻幅值的極大值與前一時刻哪一個極值點構(gòu)成優(yōu)化路徑，規(guī)劃所有路徑的可能性。U(m,tn)表示優(yōu)化結(jié)果的中間向量，選擇該時刻U(m,tn)最大值的幅值極大值點對應(yīng)的頻率值作為目標(biāo)脊線。

2 提取方法流程

為解決時變工況下無轉(zhuǎn)速信號和傳統(tǒng)時頻方法時頻聚集性不高導(dǎo)致提取的瞬時轉(zhuǎn)速誤差較大的問題，提出了基于GLCT 和FPO 的時變轉(zhuǎn)速行星齒輪箱故障診斷方法。所提方法的流程圖如圖1，具體實施步驟如下：

圖1 所提方法流程圖Fig.1 The flowchart of proposed method

步驟1：采集時變轉(zhuǎn)速下行星齒輪箱的振動信號。

步驟2：對獲取的振動信號進行GLCT 時頻表示GS(t',ω)。

步驟3：設(shè)置最大搜索頻率，在GS(t',ω)時頻域中采用FPO算法提取瞬時轉(zhuǎn)頻趨勢線。

步驟4：采用最小二乘法對提取的趨勢線進行擬合，得到目標(biāo)脊線；進而利用提取的轉(zhuǎn)速曲線對原始振動信號進行等角度重采樣，將非平穩(wěn)信號轉(zhuǎn)換為角域平穩(wěn)信號。

步驟5：對角域信號進行包絡(luò)階次分析，從包絡(luò)階次譜中提取故障特征階次并判斷故障類型。

3 仿真驗證

為了證明上述方法的有效性，仿真一個如下時變轉(zhuǎn)速下的太陽輪故障信號：

其中，xs(t)為軸的振動信號；xm(t)為太陽輪局部故障振動信號；n(t)為高斯白噪聲，具體數(shù)學(xué)模型如式(17)和(18)：

式中，(t)為軸旋轉(zhuǎn)頻率；fs(t)為太陽輪故障特征頻率；fm(t)為齒輪嚙合頻率。

考慮到時變轉(zhuǎn)速，則上述表達式為：

在這個仿真信號中，采樣頻率Fs=1024Hz,采樣時間為10 秒，n(t)=-5dB。而仿真信號中，太陽輪故障階次Os為：

上述仿真信號的時域波形及其FFT 頻譜和包絡(luò)譜如圖2所示。由圖2可知，故障沖擊被完全淹沒在噪聲中，此外受時變轉(zhuǎn)速影響，故障特征頻率以及調(diào)制成分在頻譜圖中出現(xiàn)混疊和模糊現(xiàn)象，從而不能提取到有用的故障特征信息。

圖2 仿真信號時域圖及其相關(guān)譜圖Fig.2 Time-domain diagram of simulated signal and its related spectrum diagram

圖3 為本文所提方法的分析結(jié)果，圖3(a)為仿真信號GLCT時頻表示，瞬時頻率最低的時頻脊線為理論轉(zhuǎn)頻曲線。因此，以瞬時頻率最低的時頻脊線為基頻O即為所提取的目標(biāo)脊線，進行角度重采樣。用FPO從圖3(a)的時頻表示中提取到的目標(biāo)脊線如圖3(b)的黑色曲線所示，因所提取的目標(biāo)脊線不夠光滑，本文采用最小二乘法進行曲線擬合，圖3(b)的藍色線為擬合曲線，紅色線為理論曲線。從圖3(b)可以明顯看出擬合后的曲線和理論曲線基本重合。根據(jù)上述提取出來的曲線進行角度重采樣包絡(luò)階次譜，從圖中可以明顯得到O,Os-O,Os,Os+O,2Os。從上述的分析結(jié)果可以得出本文所提出的方法能夠準(zhǔn)確提取到變轉(zhuǎn)速下仿真信號的太陽輪故障特征。

圖3 本文所提方法分析結(jié)果Fig.3 The analysis result of the method proposed in this article

為了驗證本文所提方法的優(yōu)越性，分別采用STFT與FPO 和GLCT 與峰值搜索算法進行對比。仿真信號的STFT時頻表示如圖4(a)所示，用FPO從圖4(a)中提取到的目標(biāo)脊線如圖4(b)黑色曲線所示，紅色線為最小二乘擬合曲線，藍色線為理論曲線，使用擬合后的曲線得到的包絡(luò)階次譜如圖4(c)所示?？梢钥闯觯捎谒捎玫臅r頻分析方法時頻聚集性不高，并且對噪聲比較敏感，導(dǎo)致提取出來的目標(biāo)脊線與理論曲線誤差較大，進而最后的診斷效果不理想。

圖4 STFT與FPO分析結(jié)果Fig.4 STFT and FPO analysis results

從圖5可以看出，由于峰值搜索算法考慮不到時頻矩陣的全局信息，容易出現(xiàn)局部最優(yōu)的現(xiàn)象，導(dǎo)致所提取到的目標(biāo)曲線與理論曲線誤差較大，最后的診斷效果不理想。

圖5 GLCT與峰值搜索算法結(jié)果包絡(luò)階次譜Fig.5 GLCT and peak search analysis results envelope order spectrum

綜上，本文所提出的方法，能克服傳統(tǒng)時頻分析方法的時頻聚集性不高，而帶來的目標(biāo)脊線提取精度不高以及傳統(tǒng)的時頻脊線提取算法由于考慮不到時頻矩陣的全局信息，導(dǎo)致提取到的目標(biāo)脊線局部最優(yōu)的問題，最后能夠準(zhǔn)確提取時變轉(zhuǎn)速下太陽輪故障特征，從而驗證了本文所提的方法在時變轉(zhuǎn)速下對行星齒輪箱齒輪故障特征提取的有效性。

4 實驗驗證

4.1 實驗說明

設(shè)計了行星齒輪箱太陽輪和行星輪故障在時變轉(zhuǎn)速下的模擬實驗，實驗系統(tǒng)如圖6 所示，采集實驗是在圖6所示的實驗臺上進行的，主要由變頻器、驅(qū)動電機、行星齒輪箱、磁粉制動器、轉(zhuǎn)速傳感器、加速度傳感器（型號為美國DYTRAN 公司生產(chǎn)的3035B 型傳感器）和安裝NI9234采集卡的工控機組成。

圖6 時變轉(zhuǎn)速下行星齒輪箱實驗系統(tǒng)Fig.6 Planetary gearbox experimental system under time-varying speed

在本實驗中，行星齒輪箱的齒圈固定，太陽輪為輸入端，行星架為輸出端，電機的轉(zhuǎn)速在20s內(nèi)從0～1800r/min線性增加，為了驗證提取到的目標(biāo)脊線的準(zhǔn)確度，以轉(zhuǎn)速傳感器測得的瞬時轉(zhuǎn)速為理論轉(zhuǎn)速。采樣頻率為12800Hz，為了增加分析的可靠性，本文取5～15秒的數(shù)據(jù)。磁粉制動器施加的扭矩為0，行星齒輪箱的參數(shù)見表1。為了模擬行星齒輪箱齒輪的局部故障，通過線切割技術(shù)在太陽輪的某個輪齒上切割一部分和切掉行星輪一個輪齒，這樣形成太陽輪斷齒和行星輪缺齒局部故障，故障齒輪如圖7所示。行星齒輪箱的故障特征階次見表2。

表1 行星齒輪箱的參數(shù)Tab.1 Parameters of planetary gearbox

圖7 齒輪故障件Fig.7 Gear failure parts

表2 行星齒輪箱各部分故障特征階次Tab.2 Characteristic order of failure of each part of planetary gearbox

4.2 實驗驗證1：太陽輪故障信號分析

太陽輪缺齒故障振動信號的時域波形及其對應(yīng)的傅里葉頻譜和包絡(luò)譜如圖8所示?？梢钥吹?，隨著轉(zhuǎn)速的升高，沖擊逐漸變大；同時，受時變轉(zhuǎn)速影響，譜圖中會出現(xiàn)頻譜模糊現(xiàn)象，因此根據(jù)傳統(tǒng)的時頻分析很難提取到有價值的故障特征信息。

圖8 太陽輪斷齒時域圖及其相關(guān)譜圖Fig.8 Time domain diagram of sun gear broken tooth and its related spectrum diagram

為提升運算效率，本文對原始信號進行降采樣預(yù)處理，降采樣的頻率設(shè)置為512Hz。圖9(a)為降采樣信號GLCT 時頻表示，從圖中可以明顯的看出來，瞬時頻率最低的時頻脊線為太陽輪輸入端的瞬時速度曲線，即目標(biāo)提取脊線。用FPO 從圖9(a)提取到的目標(biāo)瞬時頻率脊線如圖9(b)的黑色曲線所示，本文采用最小二乘法擬合所提取出來的目標(biāo)曲線為圖9(b)的紅色曲線，藍色線為理論曲線，從圖中可以看出擬合后的曲線和理論轉(zhuǎn)頻曲線基本重合。根據(jù)上述所擬合后的曲線進行角度重采樣后的包絡(luò)階次譜如圖9(d)所示，可以清晰地觀察到O,Os-O,Os,Os+O,2Os,2O的故障特征信息。上述分析結(jié)果表明太陽輪出現(xiàn)故障，從而驗證了所提算法的有效性。

圖9 本文所提方法分析結(jié)果Fig.9 The analysis result of the method proposed in this article

為了驗證本文所提方法的優(yōu)越性，分別采用STFT 與FPO 和GLCT 與峰值搜索算法進行對比。圖10(a)為太陽輪斷齒信號的STFT 時頻表示，用FPO 從STFT 時頻表示中提取到的目標(biāo)脊線如圖10(b)的黑色曲線所示，紅色線為擬合曲線，藍色線為理論曲線，使用擬合曲線得到的包絡(luò)階次譜如圖10(c)所示。從最后的診斷結(jié)果可以看出，所采用的時頻分析方法時頻聚集性不高，并且對噪聲敏感，導(dǎo)致提取出來的目標(biāo)曲線與理論的曲線誤差較大，進而最后的診斷效果不理想。

圖10 STFT與FPO分析結(jié)果Fig.10 STFT and FPO analysis results

GLCT 和峰值搜索的包絡(luò)階次譜如圖11 所示?？梢钥闯?，從最后診斷結(jié)果可以看出，因為峰值搜索算法一旦出現(xiàn)選錯頻率起始點，將導(dǎo)致所提取到的目標(biāo)曲線錯誤，從而使該方法失效。

圖11 GLCT與峰值搜索算法結(jié)果包絡(luò)階次譜Fig.11 GLCT and peak search analysis results envelope order spectrum

4.3 實驗驗證2：行星輪故障信號分析

行星輪缺齒故障振動信號的時域波形及其對應(yīng)的傅里葉頻譜和包絡(luò)譜如圖12 所示?？梢钥吹?，隨著轉(zhuǎn)速的升高，沖擊逐漸變大；同時，受時變轉(zhuǎn)速影響，譜圖中會出現(xiàn)頻譜模糊現(xiàn)象，因此根據(jù)傳統(tǒng)的時頻分析很難提取到有價值的故障特征信息。

圖12 行星輪缺齒時域圖及其相關(guān)譜圖Fig.12 Time domain diagram of planetary gear missing teeth and its related spectrum diagram

圖13(a)為降采樣信號GLCT時頻表示，從圖13(a)可以明顯的看出來，瞬時頻率最低的時頻脊線為目標(biāo)提取脊線。用FPO從圖13(a)提取到的目標(biāo)脊線如圖13(b)黑色曲線所示，紅色線為擬合后的曲線，藍色線為理論曲線，從圖中可以看出擬合后的曲線和理論轉(zhuǎn)頻曲線基本重合。根據(jù)上述所擬合后的曲線進行角度重采樣，最后的包絡(luò)階次譜如圖13(d)所示，可以明顯得到O,Op-Oc,Op,Op+Oc,2Op,2O,4Op的故障信息。上述分析結(jié)果表明行星齒輪出現(xiàn)故障，從而驗證了所提算法的有效性。

圖13 本文所提方法分析結(jié)果Fig.13 The analysis result of the method proposed in this article

為了驗證本文所提方法的優(yōu)越性，分別采用STFT與FPO 和GLCT 與峰值搜索算法進行對比。行星輪缺齒信號的STFT時頻表示如圖14(a)所示，用FPO從圖14(a)中提取到的目標(biāo)曲線如圖14(b)黑色曲線所示，紅色線為擬合曲線，藍色線為理論曲線,使用擬合曲線得到最后的包絡(luò)階次譜如圖14(c)所示。從最后的診斷結(jié)果可以看出，因為所采用的時頻分析方法時頻聚集性不高，并且對噪聲敏感，導(dǎo)致提取出來的時頻曲線與理論的時頻曲線誤差較大，導(dǎo)致最后的診斷效果不理想。

圖14 STFT與FPO分析結(jié)果Fig.14 STFT and FPO analysis results

GLCT和峰值搜索的包絡(luò)階次譜如圖15所示。從最后診斷結(jié)果可以看出，峰值搜索算法一旦選錯頻率起始點，將導(dǎo)致所提取到的目標(biāo)曲線錯誤，從而使該方法失效。

圖15 GLCT與峰值搜索算法結(jié)果包絡(luò)階次譜Fig.15 GLCT and peak search analysis results envelope order spectrum

綜上，本文提出的方法，能在強噪聲的背景下和時變轉(zhuǎn)速工況下，準(zhǔn)確提取到行星輪缺齒故障信息。提出的GLCT方法能夠生成時頻聚集性較高的時頻圖，克服了傳統(tǒng)的時頻方法在強噪聲背景下時頻聚集性不高的問題。提出的FPO脊線提取算法能夠考慮時頻矩陣的全局信息，準(zhǔn)確提取到目標(biāo)脊線，克服了傳統(tǒng)脊線提取算法容易陷入局部最優(yōu)的問題。

5 結(jié)論

針對無轉(zhuǎn)速計時變工況下行星齒輪箱故障特征難以提取的問題，提出了基于GLCT與FPO的故障特征提取算法。通過構(gòu)造仿真信號和實際應(yīng)用對其有效性和優(yōu)越性進行分析，可以得到如下結(jié)論：

1）在噪聲干擾比較嚴重的情況下，GLCT的時頻聚集性和噪聲魯棒性優(yōu)于傳統(tǒng)的STFT 方法，從而可以清晰地展現(xiàn)齒輪箱瞬時頻率曲線，為后續(xù)目標(biāo)脊線提取提供基礎(chǔ)。

2）FPO時頻脊線搜索算法能夠考慮到時頻矩陣的全局信息，提取到的時變轉(zhuǎn)速曲線與理論曲線能夠完全重合，克服了峰值搜索算法容易陷入局部最優(yōu)和提取目標(biāo)脊線不準(zhǔn)確的問題。

3）該算法能在不安裝轉(zhuǎn)速計等硬件的條件下提取轉(zhuǎn)速曲線，隨后對原始振動信號進行角域重采樣，在其包絡(luò)階次譜中能夠準(zhǔn)確識別到行星齒輪箱故障特征階次，為變轉(zhuǎn)速工況下的行星齒輪箱故障診斷提供了一種新思路。