基于高階同步壓縮變換的行星齒輪箱聲音信號共振頻帶特征提取

張鄭武，馮志鵬，陳小旺

北京科技大學(xué)機(jī)械工程學(xué)院，北京 100083

行星齒輪箱故障診斷具有重要意義[1?4]. 目前，行星齒輪箱的故障診斷研究主要基于振動信號的測試分析，對聲音信號鮮有涉及. 聲音信號蘊含著機(jī)械設(shè)備自身的結(jié)構(gòu)和運行狀態(tài)等重要的信息. 行星齒輪箱運行過程中產(chǎn)生機(jī)械振動，并引發(fā)空氣振動向四周傳播，從而產(chǎn)生聲音信號. 其中，齒輪嚙合產(chǎn)生的沖擊和振動是行星齒輪箱的主要聲源. 相對于振動信號，聲音信號可通過非接觸的方式采集，具有安裝簡便、操作簡單等優(yōu)點[5]，尤其當(dāng)處在高溫、高腐蝕環(huán)境下振動傳感器無法使用時，仍能從聲音信號中獲得可靠的機(jī)械設(shè)備重要信息. Lin[6]對機(jī)械系統(tǒng)聲音信號進(jìn)行小波降噪，并成功提取了車輛發(fā)電機(jī)的故障特征；Rezaei等[7]利用多個聲音傳感器對不同轉(zhuǎn)速和負(fù)載下的滾動軸承聲音信號進(jìn)行了實驗分析；Metwalley等[8]將齒輪箱正常與故障狀態(tài)下的聲音信號時域、頻域特征進(jìn)行比較，判斷了齒輪故障的嚴(yán)重程度. Bayda和Ball[9]將聲音信號與振動信號進(jìn)行比較分析，發(fā)現(xiàn)聲音信號在早期故障的診斷方面具有優(yōu)勢.

然而，目前針對行星齒輪箱的聲音信號特征提取理論方法還十分有限. 行星齒輪箱具有特殊的齒輪結(jié)構(gòu)和運動形式，齒輪嚙合點與故障接觸位置具有周期性變化的特征，因此振動信號存在復(fù)雜的幅值?頻率調(diào)制特性[10?11]. 進(jìn)一步考慮聲音傳播過程的影響，實際測得的聲音信號結(jié)構(gòu)更加復(fù)雜，其頻率結(jié)構(gòu)及故障特征分布規(guī)律尚未明晰.

傳統(tǒng)的齒輪箱故障特征建模與提取大都以齒輪嚙合頻率或其倍頻作為載波頻率，集中在嚙合頻帶進(jìn)行故障提取. 然而行星齒輪箱常在時變轉(zhuǎn)速條件下運行，其嚙合頻率及調(diào)制邊帶隨轉(zhuǎn)速變化，容易與共振帶交叉，增加故障特征的提取難度. 齒輪故障引起的沖擊會激發(fā)行星齒輪箱共振，并對共振頻率產(chǎn)生幅值調(diào)制和頻率調(diào)制作用[12?13].同時，作為載波頻率的共振頻率不隨轉(zhuǎn)速變化，在共振頻帶進(jìn)行故障特征提取，能夠有效避免頻率交叉現(xiàn)象，為時變工況下的齒輪故障診斷提供新的思路. 因此，本文首先建立行星齒輪箱聲音信號的調(diào)幅?調(diào)頻模型，并推導(dǎo)其故障特征分布規(guī)律.

另外，行星齒輪箱非平穩(wěn)聲音信號的故障特征提取仍是目前研究難題和熱點之一. 傳統(tǒng)的時頻分析方法包括線性時頻分析方法，如短時Fourier變換和連續(xù)小波變換；以及雙線性時頻分布，如Wigner?Ville分布、Cohen類分布等[14]均可用于分析時變信號. 但短時Fourier變換和連續(xù)小波變換受Heisenberg不確定性原則的限制，無法同時獲得良好的時間分辨率和頻率分辨率[15?16]，Wigner?Ville分布雖然改善了時頻分布的分辨率，但存在著交叉項干擾而不適于多分量信號分析[17]. Auger和Flandrin[18]提出了時頻重分配的方法，有效提高了時頻可讀性，但存在計算量較大問題；Daubechies等[19]提出了同步壓縮變換方法，進(jìn)一步提高了時頻分析的可讀性并保留了可逆特性，但難以抑制時間方向上的模糊現(xiàn)象[20?23]；Pham和Meignen[24]以及Oberlin等[25]在同步壓縮變換的基礎(chǔ)上進(jìn)一步提出了高階同步壓縮變換方法，相比傳統(tǒng)同步壓縮變換而言具有更精確的瞬時頻率估計模型，適合于分析快速變化的頻率成分.

本文建立了行星齒輪箱近場聲源聲音信號在共振頻帶內(nèi)的齒輪故障調(diào)幅?調(diào)頻模型，推導(dǎo)時變工況下的聲音信號的時變頻率結(jié)構(gòu)；并利用高階同步壓縮變換分析方法對聲音信號進(jìn)行故障特征提取.

1 聲音信號共振頻帶調(diào)制模型

根據(jù)聲源與傳感器距離的不同，聲場模型可分為近場模型和遠(yuǎn)場模型. 本文重點研究行星齒輪箱聲音信號的近場模型. 根據(jù)聲學(xué)理論，聲壓方程可表示為：

式中，p為聲壓，r為聲源輻射半徑，ω 為聲源角頻率，c0為聲速，C0為常數(shù)，t為聲音傳播時間，j為虛數(shù)單位，J0(·)為第0階貝塞爾函數(shù)，可表示為：

式（1）描述了聲音傳播過程中，聲壓與傳播距離及聲源頻率的關(guān)系，聲音信號中聲源頻率成分保持不變，但幅值在傳播過程中逐漸衰減，衰減形式如式（2）所示. 幅值衰減因子C0J0(ωr/c0)可用來描述聲音信號的幅值衰減. 當(dāng)行星齒輪箱齒輪出現(xiàn)故障時，故障齒輪在運行過程中與其他齒輪嚙合產(chǎn)生沖擊，激發(fā)齒輪箱共振. 在系統(tǒng)阻尼作用下，共振迅速衰減，又在后續(xù)沖擊激勵下重復(fù)出現(xiàn)，對共振頻率產(chǎn)生幅值調(diào)制和頻率調(diào)制作用. 這種重復(fù)出現(xiàn)的共振產(chǎn)生聲音輻射，因此聲音信號在共振頻帶包含齒輪故障特征信息. 通常，近場聲源的聲音信號采集時，聲音傳感器的位置相對于聲源距離較近且固定不變，加之共振頻率僅受系統(tǒng)動力學(xué)參數(shù)決定. 因此，聲音信號的幅值衰減因子C0J0(ωr/c0)可近似為常數(shù). 結(jié)合行星齒輪箱振動信號的調(diào)制形式[12]，由振動產(chǎn)生的對外輻射的聲音信號可表示為：

其中，C為幅值因子，A代表共振信號強(qiáng)度，a1、a2為幅值調(diào)制系數(shù)，a3為頻率調(diào)制系數(shù)，?1、?2和?3為對應(yīng)的初始相位. 對（3）進(jìn)行Fourier變換展開，得到聲音信號內(nèi)含有的特征頻率結(jié)構(gòu)為：

式中，m代表Bessel函數(shù)展開的階數(shù)，fn為共振頻率，ft(t)為太陽輪旋轉(zhuǎn)頻率或行星架旋轉(zhuǎn)頻率（太陽輪或行星輪出現(xiàn)故障時由于傳遞路徑的時變現(xiàn)象產(chǎn)生的調(diào)制頻率，因齒圈固定，當(dāng)齒圈故障時此項可略去），ff(t)為故障特征頻率，Jm(·)為第一類Bessel函數(shù)，δ (·)為Dirac函數(shù).

根據(jù)（4）可總結(jié)出時變轉(zhuǎn)速工況下行星齒輪箱齒輪故障時聲音信號的時變故障頻率特征：在太陽輪故障條件下，聲音信號的共振頻帶內(nèi)存在分布于共振頻率fn兩側(cè)的、間隔為太陽輪故障特征頻率ff(t)和太陽輪旋轉(zhuǎn)頻率ft(t)的組合的一系列時變邊帶成分，即fn±m(xù)ff(t)±nft(t)（其中m,n=0,1,2···代表公式4中的階數(shù)，下同）；在行星輪故障條件下，聲音信號共振頻帶內(nèi)存在與共振頻率fn間隔為行星輪故障特征頻率ff(t)和行星架旋轉(zhuǎn)頻率ft(t)的組合的時變邊帶，即fn±m(xù)ff(t)±nft(t)；在齒圈故障條件下，聲音信號的共振頻帶內(nèi)存在與共振頻率fn間隔為齒圈故障特征頻率ff(t)的時變邊帶fn±m(xù)ff(t). 其中，故障特征頻率對稱分布在共振頻率fn兩側(cè)，可根據(jù)聲音信號的時變邊帶對稱中心定位共振頻率，并可通過識別時變邊帶信息，有效檢測及定位齒輪故障.

2 高階同步壓縮變換

短時Fourier變換（STFT）是最為常用的時頻分析方法之一，其本質(zhì)是加窗的Fourier變換. 對于一個有限能量信號x(t)，其STFT的表達(dá)式為：

其中，τ為時間單位，g(t?τ)為窗函數(shù). 由于窗函數(shù)尺寸固定，其時頻分辨能力固定，在Heisenberg不確定性的制約下，難以同時得到復(fù)雜多變的瞬時頻率細(xì)節(jié)的最佳表達(dá).

基于短時Fourier變換的同步壓縮變換方法（FSST）的關(guān)鍵思想是通過對瞬時頻率進(jìn)行精確估計，從而在頻率方向上對時頻能量進(jìn)行重分配，F(xiàn)SST方法中，瞬時頻率可根據(jù)短時Fourier變換的時頻分布 S TFTx(t,f)計算為：

其中，? {·}表示實數(shù)部分，基于此得到STFT頻率重分配后的FSST結(jié)果為：

其中，g(0)為窗函數(shù)，由（7）可知，通過將頻率方向上能量模糊壓縮到瞬時頻率脊線附近，可有效提高時頻聚集性，有利于準(zhǔn)確識別變化的瞬時頻率特征.

然而，雖然FSST提高了時頻分布的頻域聚集性，但由于僅在頻率方向上進(jìn)行重分配，時域的能量溢出現(xiàn)象未得到有效解決，依然存在時頻模糊問題. 為了克服這個不足，基于更精確的瞬時頻率估計的二階同步壓縮變換（FSST2）被提出[25]. FSST2首先定義一個二階局部調(diào)制系數(shù)：

該調(diào)制系數(shù)為FSST瞬時頻率估計（6）對時間t的一階導(dǎo)數(shù). 則二階瞬時頻率估計可表達(dá)為：

因而，信號x(t)的二階同步壓縮變換如下：

FSST2被證明能夠提高具有高斯幅值調(diào)制信號的時頻分布的效果. Pham和Meignen[24]針對其他更為復(fù)雜的快速變化頻率信號，進(jìn)一步提出了具有更為通用形式的高階同步壓縮變換方法（FSSTN）. FSSTN基于信號的幅值與相位的泰勒展開定義了新的瞬時頻率估計：

類似于FSST2，F(xiàn)SSTN的表達(dá)式為：

利用FSSTN方法，可以獲得更加精確的瞬時頻率估計，從而有效抑制時頻模糊現(xiàn)象，使時頻分布的可讀性進(jìn)一步提高. 根據(jù)實際應(yīng)用需求，本文采用四階同步壓縮變換（FSST4），即取N=4.

3 仿真信號分析

在本節(jié)中，通過對時變轉(zhuǎn)速工況下行星齒輪箱聲音信號進(jìn)行仿真分析，驗證前文總結(jié)的時變故障特征規(guī)律形式，并驗證FSST4方法提取共振頻帶內(nèi)故障特征的有效性. 以太陽輪故障為例，由式（3）可得聲音信號模型如下：

式中，fn為 行星齒輪箱共振頻率，為太陽輪旋轉(zhuǎn)頻率，fs(t)為太陽輪故障特征頻率. 取初始相位?=θ=0，共振頻率fn=2000 Hz，幅值調(diào)制系數(shù)a=1.2，頻率調(diào)制系數(shù)B=0.1，幅值衰減因子C=0.98，加入信噪比為0 dB的Gauss白噪聲 ξ (t). 令太陽輪旋轉(zhuǎn)頻率太陽輪故障特征頻率為信號的采樣頻率為6000 Hz，時間長度為3 s.

圖1 仿真信號. （a）波形；（b）Fourier頻譜；（c）STFT 時頻分布；（d）Wigner?Ville 分布；（e）FSST；（f）FSST4Fig.1 Simulation signal:(a) waveform; (b) Fourier spectrum; (c) time-frequency representation(TFR) by STFT; (d) Wigner-Ville distribution; (e) timefrequency representation by FSST; (f) time-frequency representation by FSST4

仿真信號的分析結(jié)果如圖1所示. 圖1（a）和（b）分別為信號的時域波形和Fourier頻譜分析結(jié)果. 圖 1（c）、（d）、（e）和（f）為時頻分布結(jié)果，右側(cè)顏色條表征信號的幅值大小，顏色越深幅值越大（下同）. 由圖 1（c）可見，STFT 雖然正確地顯示出信號的瞬時頻率變化趨勢，但由于邊帶成分密集，且時頻分辨率不足，時變故障特征提取困難；圖1（d）Wigner?Ville時頻分布雖然具有高時頻分辨率特點，但存在嚴(yán)重交叉項干擾，無法準(zhǔn)確辨識故障特征頻率；圖1（e）FSST時頻分布盡管得到了較為清晰的瞬時頻率變化趨勢，但瞬時頻率變化較快時，仍存在時間方向的模糊現(xiàn)象，如圖1（e）中紅色標(biāo)記區(qū)域，因此難以精準(zhǔn)識別密集分布的時變邊帶；圖1（f）的FSST4時頻分析方法相較于FSST采用了更為精確的瞬時頻率估計公式，能夠更有效地捕捉快速變化的頻率信息，顯著提高了時頻分辨率.

4 實驗信號分析

4.1 實驗說明

本節(jié)進(jìn)一步通過實驗信號分析驗證聲音信號共振頻帶的時變故障特征，以及FSST4時頻分析方法在實際行星齒輪箱故障診斷中的效果. 該實驗在加拿大Ottawa大學(xué)機(jī)械工程系實驗室完成.行星齒輪箱實驗裝置如圖2，實驗臺由電動機(jī)、兩級行星齒輪箱、定軸齒輪箱、轉(zhuǎn)速計、加速度傳感器、聲壓傳感器和磁粉制動器構(gòu)成. 各級齒輪箱齒輪齒數(shù)見表1. 根據(jù)各齒輪箱的齒數(shù)參數(shù)以及輸入轉(zhuǎn)頻fd(t)，可計算得各齒輪故障特征頻率，見表1.其中分別為一級行星齒輪箱的太陽輪旋轉(zhuǎn)頻率、太陽輪故障特征頻率、行星輪故障特征頻率和齒圈故障特征頻率.

實驗包含健康和故障狀態(tài)下兩組. 在健康狀態(tài)下，各齒輪理論上均無故障；在故障狀態(tài)下，人為在第一級行星齒輪箱的太陽輪上加工輪齒磨損，如圖2所示. 兩組實驗的運行工況相同，電動機(jī)轉(zhuǎn)頻在30 Hz上下正弦波動，其中最大轉(zhuǎn)頻為40 Hz，最小轉(zhuǎn)頻 20 Hz，電機(jī)轉(zhuǎn)頻表達(dá)式近似為fd(t)=30+10sin(0.06πt). 聲壓傳感器的采樣頻率為 20 kHz，采樣時間為 25 s.

4.2 實驗信號分析

4.2.1 正常狀態(tài)信號分析

正常狀態(tài)下，傳統(tǒng)的STFT時頻分析方法和FSST4時頻分析方法被用于該聲音信號分析，所得到的時頻分布如圖 3（a）和圖 3（b）. 由圖可知，F(xiàn)SST4時頻分析方法獲得了更好的時頻分析效果，可讀性更高. 根據(jù)行星齒輪箱共振頻率不隨轉(zhuǎn)速變化的性質(zhì)，以及在共振頻率兩側(cè)呈現(xiàn)對稱分布的調(diào)制頻率邊帶特點，識別出共振頻率為5003 Hz. 進(jìn)一步驗證共振頻帶內(nèi)的時變頻率邊帶，可得調(diào)制頻率為電機(jī)轉(zhuǎn)頻，即存在時變頻率成分如fn±fd(t)、fn±3fd(t)和fn±5fd(t). 由于制造及安裝誤差，齒輪旋轉(zhuǎn)頻率對振動及聲音信號的調(diào)制作用難以避免，由于未檢測到明顯的故障頻率調(diào)制邊帶，表明齒輪未出現(xiàn)故障，與實際情況相符.

圖2 實驗裝置. 1—電動機(jī)；2—定軸齒輪箱；3,5—行星齒輪箱；4—聲壓傳感器；6—磁粉制動器；7—太陽輪故障Fig.2 Test rig:1—motor; 2—fixed-shaft gearbox; 3, 5—planetary gearbox; 4—microphone; 6—magnetic powder brake; 7—sun gear fault

表1 齒輪箱主要參數(shù)Table 1 Main parameters of gearboxes

4.2.2 太陽輪故障信號分析

引入太陽輪故障后，在相同工況下采集的聲音信號，利用STFT和FSST4時頻分析方法分別對太陽輪故障狀態(tài)信號進(jìn)行分析，結(jié)果如圖4（a）和（b）所示. 同樣，F(xiàn)SST4獲得更好的時頻分析結(jié)果，依據(jù)共振頻率的時不變特性及共振邊帶對稱分布特點，可識別出共振頻率為5003 Hz，與正常狀態(tài)下的共振頻率一致. 利用FSST4方法的高時頻分辨率優(yōu)勢，可識別出一系列共振調(diào)制邊帶. 相較于圖3（b）中健康狀態(tài)下的信號時頻分布結(jié)果，故障狀態(tài)下還存在額外的系列邊帶，包括太陽輪故障頻率fs(t)在電機(jī)轉(zhuǎn)頻的調(diào)制邊帶fn±fd(t)、fn±3fd(t)成 分 兩 側(cè) 的 額 外 邊 帶fn+fd(t)±fs(t)，fn?fd(t)±fs(t)，以及fn+3fd(t)±fs(t)，fn?3fd(t)±fs(t).通過與健康狀態(tài)的分析結(jié)果對比，這些故障相關(guān)的時變邊帶表明行星齒輪箱的太陽輪出現(xiàn)了故障，與實際情況完全吻合.

圖3 正常狀態(tài)聲音信號分析. （a）STFT 時頻分布；（b）FSST4 時頻分布Fig.3 Acoustic signal analysis under normal conditions:(a) TFR by STFT; (d) TFR by FSST4

圖4 太陽輪故障狀態(tài)聲音信號分析；（a）STFT 時頻分布；（b）FSST4 時頻分布Fig.4 Acoustic signal analysis under sun gear fault conditions:(a) TFR by STFT; (b) TFR by FSST4

5 結(jié)論

聲音信號中蘊含著機(jī)械裝備運行過程中的健康狀態(tài)信息. 本文首先建立了行星齒輪箱近場聲源聲音信號共振頻帶調(diào)制模型，并推導(dǎo)了非平穩(wěn)工況下齒輪故障特征時變邊帶在共振頻率兩側(cè)的分布規(guī)律. 進(jìn)一步利用共振頻率不隨轉(zhuǎn)速變化的特性，在行星齒輪箱聲音信號中準(zhǔn)確識別了共振頻率. 最后，在行星齒輪箱聲音信號共振頻帶故障特征分布規(guī)律的指導(dǎo)下，應(yīng)用高階同步壓縮變換方法，克服了傳統(tǒng)時頻分析方法時頻分辨率低及干擾項嚴(yán)重的不足，實現(xiàn)了密集時變邊帶的準(zhǔn)確識別，準(zhǔn)確提取了時變的行星齒輪箱太陽輪故障特征. 通過數(shù)值仿真和實驗信號分析，驗證了聲音信號共振頻帶調(diào)制模型和故障特征規(guī)律的準(zhǔn)確性，以及共振頻率識別方法和高階同步壓縮變換提取時變故障特征的有效性，成功診斷了行星齒輪箱的齒輪故障.