秦嗣峰，馮志鵬，LIANG Ming

（1.北京科技大學機械工程學院，北京100083；2.渥太華大學機械工程系，安大略渥太華K1N 6N5，加拿大）

Vold-Kalman濾波和高階能量分離在時變工況行星齒輪箱故障診斷中的應用研究

秦嗣峰1，馮志鵬1，LIANG Ming2

（1.北京科技大學機械工程學院，北京100083；2.渥太華大學機械工程系，安大略渥太華K1N 6N5，加拿大）

時變工況下的行星齒輪箱振動信號具有明顯的時變調(diào)制特點，常規(guī)的頻譜分析方法難以識別齒輪故障特征頻率。提出了基于Vold-Kalman濾波和能量分離的時頻分析方法，識別行星齒輪箱的時變特征頻率，診斷齒輪故障。與傳統(tǒng)的時頻分析方法相比，基于Vold-Kalman濾波和能量分離的時頻表示具有良好的時頻分辨率，而且沒有交叉項干擾，能夠有效提取非平穩(wěn)信號中的時變頻率成分。分析了行星齒輪箱時變工況下的實驗信號，準確地診斷了齒輪故障，驗證了該方法的有效性。

故障診斷；行星齒輪箱；高階能量分離；Vold-Kalman濾波；時頻分析

引 言

行星齒輪箱在風力發(fā)電機、直升機和車輛等裝備中應用廣泛，實際運行中，其工況經(jīng)常發(fā)生變化，產(chǎn)生非平穩(wěn)的振動響應。目前的研究主要集中在恒定工況方面［1-2］，而針對時變工況方面的研究并不多見，其中，Williams和Zalubas［3］考慮直升機傳動系統(tǒng)振動信號的非平穩(wěn)特點，應用Wigner-Ville分布檢測行星齒輪箱故障。Bartelmus和Zimroz［4-5］發(fā)現(xiàn)行星齒輪箱故障對載荷變化敏感，應用齒輪嚙合振動幅值與工況之間的線性依賴關系監(jiān)測時變工況下的行星齒輪箱運行狀態(tài)。

行星齒輪箱故障診斷的關鍵是根據(jù)故障特征頻率及其幅值的變化識別故障原因。行星齒輪箱嚙合頻率和齒輪故障特征頻率與轉速直接相關。工況（轉速和負荷）的變化將引起轉速波動，嚙合頻率和齒輪故障特征頻率也會隨之發(fā)生變化。因此，如何從非平穩(wěn)信號中有效地識別齒輪故障特征頻率及其時變特征，是時變工況下行星齒輪箱故障診斷的關鍵問題。

對于時變工況下的行星齒輪箱的非平穩(wěn)振動信號，時頻分析是一種有效的分析手段。但是，常見的時頻分析方法各有優(yōu)缺點［6-7］。例如：線性時頻表示（如短時Fourier變換和小波變換）雖然不存在交叉項干擾問題，但是其時頻分辨率受Heisenberg不確定性原理的限制，不能同時達到最佳；Wigner-Ville分布雖然時頻分辨率高，卻存在固有的交叉項干擾問題，不適合分析復雜多分量信號；以Wigner-Ville分布為基礎的雙線性時頻分布（包括Cohen類分布和仿射類分布）通過各種核函數(shù)進行平滑處理，雖然抑制了交叉項干擾，但是會造成信號自項畸變，降低時頻分辨率；Hilbert-Huang變換雖然對信號的時變具有自適應性，而且具有良好的局部時頻聚集能力，但是對信號中的奇異點敏感，容易產(chǎn)生模式混淆，得到虛假的本質(zhì)模式函數(shù)，影響分析結果。

行星齒輪箱振動信號具有明顯的調(diào)制特點，在時變工況下，將表現(xiàn)出時變調(diào)制特征。非線性信號處理中的能量算子方法為分析這種時變調(diào)制信號提供了有效的途徑。能量算子通過信號的瞬時幅值及其微分的非線性組合估計信號的“能量”，在其基礎上提出的能量分離算法可以計算任意調(diào)制信號的幅值包絡和瞬時頻率［8-12］。能 量算 子 和 能 量分離在信號及其導數(shù)的基礎上進行計算，算法簡單、時間分辨率高，對信號的瞬態(tài)變化具有良好的適應性。但是，該方法僅適用于單分量信號。對于成分復雜的行星齒輪箱振動信號，需要將其分解為單分量再進行分析。

近年來提出的Vold-Kalman濾波方法［13-15］在最小化結構和數(shù)據(jù)方程誤差的基礎上，可以有效分解復雜多分量信號的頻率成分，為解決頻率交錯重疊的行星齒輪箱振動信號單分量分解問題提供了一種有效手段。

針對時變工況下行星齒輪箱故障診斷面臨的復雜非平穩(wěn)振動信號分析問題，考慮Vold-Kalman濾波方法在單分量分解以及能量分離方法在瞬時頻率計算方面的獨特優(yōu)勢，提出了基于二者的時頻分析方法。

1 Vold-Kalman濾波

Vold-Kalman濾波可以將復雜多分量信號分解為單分量。與傳統(tǒng)的濾波方法相比，該方法直接從時域中提取感興趣的信號分量，避免了由時域至頻域變換帶來的相位偏差；與經(jīng)驗模式分解方法相比，Vold-Kalman濾波的中心頻率可以根據(jù)瞬時頻率自適應調(diào)節(jié)，能夠有效分離在時頻域內(nèi)鄰近甚至交叉的信號分量［14-15］，避免 了 由經(jīng)驗 模 式分 解 方法帶 來的模式混淆現(xiàn)象。

任意調(diào)制信號可以表示為

式中k為階次，Ak（t）為第k個分量的幅值包絡，Θk（t）為載波信號。

式中ω（τ）為瞬時角頻率（τ）dτ為瞬時相位。

與載波信號相比，幅值包絡以較低頻率緩慢變化，可以用低階多項式表示。對于離散信號

式中 ▽為差分算子，s為差分階次，εk（n）為非齊次項。設多項式階次為2，根據(jù)式（3）得

當n=1時，Ak（0）-2Ak（1）+Ak（2）=εk（1）。由于實際信號中，Ak（0）不存在，因此將Ak（0）視為0，可得-2Ak（1）+Ak（2）=εk（1）。同樣，當n=N時，Ak（N+1）不存在，將Ak（N+1）視為0，可得Ak（N-1）-2Ak（N）=εk（N），式中N為采樣點數(shù)。

將式（4）展開成矩陣形式

不同階次的多項式具有相同的矩陣形式

式中M為N×N矩陣。

實際測試信號y（n）由各階分量與噪聲或誤差δ（n）組成，則

式中y（n）是實測振動信號，δ（n）為噪聲或誤差項。

為了獲得某一感興趣的信號分量xk（n），由式（7）得

其矩陣形式為

根據(jù)常規(guī)的時頻分析（如短時Fourier變換）的時頻脊線，可以估計感興趣信號分量的瞬時頻率ωk（n），從而獲得其載波矩陣C。

根據(jù)式（6）和（9），應用最小二乘法使非齊次項和噪聲或誤差項的平方和最小［16］，即

式中C*為C共軛轉置；A*為A共軛轉置；r為加權因子，決定了濾波器的頻率分辨力。r越大，濾波器的分辨力越高，但幅值包絡的收斂速度越慢；相反，r越小，濾波器的分辨力越低，但幅值包絡的收斂速度越快。為了保證幅值包絡的收斂性及計算的可行性，r的取值不能過大，其最大值如表1所示［17］。

表1 加權因子最大值Tab.1 The maximum value of the weighting factor

由式（10）求得各信號分量的幅值包絡矩陣A為

由載波矩陣C和幅值包絡矩陣A可以重構感興趣的各信號分量

2 高階能量分離算法

能量算子和能量分離方法完全由信號數(shù)據(jù)驅(qū)動，無需構造任何基函數(shù)，能夠適應信號的瞬時變化，時間分辨率高，適合檢測信號幅值和頻率的非平穩(wěn)變化特征，結合時頻分析方法，能夠分析處理任意時變調(diào)制信號。高階能量算子不僅具有普通二階能量算子的各種優(yōu)點，而且對噪聲干擾的魯棒性更好［11，12］，因此 本 文應用高 階 能量 算 子進行 能 量分離，計算信號的幅值包絡和瞬時頻率。

對于任意連續(xù)信號x（t），高階能量算子定義 為［10］

式中k為高階能量算子的階次，為任意整數(shù)；x（k）為x（t）的第k階導數(shù)

當k為負數(shù)時，x（k）（t）由微分變?yōu)榉e分，將削弱高階能量算子檢測瞬態(tài)成分的能力，因此本文只關注k為正數(shù)時的高階能量算子。通常，高階能量算子可由低階能量算子通過遞歸推導得到

由式（13），（14）可知，當k=2時，γ2（x）=（˙x）2-x¨x，即為二階Teager能量算子［8］；零階能量算子為γ0（x）=x-x2；一階微分能量算子對于任何信號都為零。在高階能量算子中，三階和四階能量算子比較實用，分別定義為

其中

對于簡諧信號x（t）=Acos（ωt+θ），可以得到高階能量算子輸出的遞歸方程

初始條件為E0=-A2，E1=0。由式（18）遞推，可得

對于簡諧信號x（t）=Acos（ωt+θ），由四階和二階能量算子可以估計信號的頻率和幅值

對于調(diào)幅-調(diào)頻信號

若幅值包絡a（t）和瞬時頻率ω（t）變化的速率和幅度相對于載波頻率來說不大，則依然可以由四階和二階能量算子估計信號的時變頻率和包絡

3 分析過程

（1）應用短時Fourier變換等常規(guī)時頻分析方法對信號進行初步時頻分析，判斷感興趣的信號成分。根據(jù)時頻脊線估計感興趣信號成分的瞬時頻率。

（2）根據(jù)已獲得的瞬時頻率曲線，構造其載波信號，從而獲得感興趣信號成分的載波矩陣。

（3）根據(jù)感興趣信號成分的載波矩陣，應用Vold-Kalman濾波將感興趣的信號成分分解為單分量，滿足能量分離方法的要求。

（4）應用高階能量分離方法計算各單分量信號成分xi（t）的四階“能量”Ei（t）和瞬時頻率fi（t）。

（5）根據(jù)各單分量信號成分的四階“能量”和瞬時頻率，構造信號的高階時頻“能量”分布

式中Ei（t）為第i個單分量在t時刻的四階“能量”，fi（t）為第i個單分量在t時刻的瞬時頻率，δ（·）為Dirac函數(shù)。

4 實驗信號分析

4.1 實驗說明

本實驗在Ottawa大學機器監(jiān)控實驗室完成。圖1為該實驗室的風電機組行星齒輪傳動實驗臺，其中行星齒輪箱和定軸齒輪箱均為兩級結構，齒輪參數(shù)分別見表2和3所示。該實驗臺的具體傳動路線為：電機→定軸齒輪箱→第1級行星齒輪箱→第2級行星齒輪箱。為了模擬行星齒輪箱的齒輪故障，在第1級太陽輪的某個輪齒上加工了剝落損傷，如圖2所示。在行星齒輪箱箱體頂部安裝了加速度傳感器，采集行星齒輪箱的振動信號。在自然停機過程中采集了振動信號，采樣頻率為20 k Hz，行星齒輪箱第2級輸出軸承受的負荷為16.284 N·m。為了計算的可行性，僅取60至40 Hz降速過程內(nèi)的振動信號進行分析。

圖1 風電機組行星齒輪傳動實驗系統(tǒng)Fig.1 Wind turbine planetary gearbox experimental system

表2 行星齒輪箱參數(shù)Tab.2 Planetarygearboxconfigurationparameter

表3 定軸齒輪箱參數(shù)Tab.3 Fixed-shaftgearboxconfigurationparameter

圖2 太陽輪損傷Fig.2 Sun gear damage

根據(jù)行星齒輪箱和定軸齒輪箱的參數(shù)及電機轉速，可計算出第1級行星齒輪箱齒輪特征頻率與電機轉速之間的關系：

嚙合頻率

太陽輪旋轉頻率

行星架旋轉頻率

太陽輪局部故障特征頻率

行星輪局部故障特征頻率

齒圈局部故障特征頻率

4.2 實驗信號分析

4.2.1 正常信號

圖3為正常行星齒輪箱降速過程中驅(qū)動電機的轉速變化曲線。在降速過程中，驅(qū)動電機的最高轉速為60 Hz，相應地，第1級行星齒輪箱的最高嚙合頻率為222.222 Hz。由于高階倍頻處特征頻率的變化規(guī)律與基頻處相同，不失一般性，以下分析的頻率范圍為0～400 Hz，覆蓋第1級行星齒輪箱最高嚙合頻率的3/2倍，振動信號中包含該齒輪箱的健康狀態(tài)信息。

圖3 正常狀態(tài)電機轉速Fig.3 Motor speed of normal gearbox

圖4（a）為正常信號的短時Fourier變換譜圖。從圖中可以看出信號成分主要包括嚙合頻率fm、嚙合頻率與太陽輪故障特征頻率之差fm-fs、以及嚙合頻率、太陽輪故障特征頻率的2倍頻和太陽輪旋轉頻率的組合fm±2fs-fsr，但是時頻分辨率較低。

圖4 正常信號Fig.4 Normal gearbox signal

應用Vold-Kalman濾波分離上述4個主要信號分量，然后應用高階能量分離算法計算各分量的4階能量和瞬時頻率，最后應用各分量的瞬時頻率和4階能量構造高階時頻能量分布，結果如圖4（b）所示。該方法得到的時頻分布的分辨率高，而且沒有交叉項干擾?？梢郧逦刈R別出嚙合頻率fm、嚙合頻率與太陽輪故障特征頻率之差fm-fs、以及嚙合頻率、太陽輪故障特征頻率的2倍頻和太陽輪旋轉頻率的組合fm±2fs-fsr，但它們的幅值均較小。齒輪制造和安裝誤差以及微小的損傷在所難免，導致振動信號中出現(xiàn)上述頻率成分，這些現(xiàn)象符合正常信號的理論預期結果。

4.2.2 太陽輪故障信號

圖5為太陽輪故障狀態(tài)下驅(qū)動電機的轉速變化曲線。在降速過程中，驅(qū)動電機的最高轉速仍為60 Hz，因此，以下分析頻率范圍仍為0～400 Hz。

圖5 太陽輪故障狀態(tài)電機轉速Fig.5 Motor speed of faulty sun gearbox

圖6（a）為太陽輪故障信號的短時Fourier變換譜圖。從圖中可以看出信號成分主要包括嚙合頻率fm、嚙合頻率與太陽輪故障特征頻率及其倍頻的組合fm±nfs、以及嚙合頻率、太陽輪故障特征頻率及其倍頻和太陽輪旋轉頻率及其倍頻的組合fm±nfs-nfsr等7個頻率成分，但是時頻分辨率較低。

應用Vold-Kalman濾波和高階能量分離方法構造的高階時頻能量分布如圖6（b）所示。嚙合頻率fm、嚙合頻率與太陽輪故障特征頻率之差fmfs、以及嚙合頻率、太陽輪故障特征頻率的2倍頻和太陽輪旋轉頻率的組合fm±2fs-fsr依然存在，但在1～2 s階段，嚙合頻率與太陽輪故障特征頻率之差fm-fs的能量明顯增強。此外，還出現(xiàn)了嚙合頻率與太陽輪故障特征頻率及其3倍頻之和fm+fs，fm+3fs以及嚙合頻率、太陽輪故障特征頻率和太陽輪旋轉頻率2倍頻的組合fm+fs-2fsr等3個頻率成分。與正常信號相比，該信號的高階時頻能量分布中出現(xiàn)了更多的頻率成分，且這些頻率均和太陽輪故障直接相關，說明太陽輪出現(xiàn)了故障，符合實驗中的實際情況。

圖6 太陽輪故障信號Fig.6 Faulty sun gearbox signal

5 結 論

風力發(fā)電機組、直升機和車輛等裝備中的行星齒輪箱的運行工況經(jīng)常變化，振動信號具有非平穩(wěn)性。針對時變工況下行星齒輪箱振動信號的時變調(diào)制特點，發(fā)揮Vold-Kalman濾波單分量分解和高階能量分離時間分辨率高的優(yōu)勢，提出了一種具有良好的時頻聚集能力、無交叉項干擾的時頻分析方法。應用該方法分析了行星齒輪箱時變工況下的實驗信號，準確識別信號中的時變頻率成分，正確診斷了齒輪故障，驗證了該方法的有效性。

［1］ 雷亞國，何正嘉，林京，等.行星齒輪箱故障診斷技術的研究進展［J］.機械工程學報，2011，47（19）：59—67. Lei Yaguo，He Zhengjia，Lin Jing，et al.Research advances of fault diagnosis technique for planetary gearboxes［J］.Journal of Mechanical Engineering，2011，47（19）：59—67.

［2］ Lei Yaguo，Lin Jing，Zuo M J，et al.Condition monitoring and fault diagnosis of planetary gearboxes：a review［J］.Measurement，2014，48：292—305.

［3］ Williams W J，Zalubas E J.Helicopter transmission fault detection via time-frequency，scale and spectral methods［J］.Mechanical Systems and Signal Processing，2000，14（4）：545—559.

［4］ Bartelmus W，Zimroz R.Vibration condition monitoring of planetary gearbox under varying external load［J］.Mechanical Systems and Signal Processing，2009，23：246—257.

［5］ Bartelmus W，Zimroz R.A new feature for monitoring the condition of gearboxes in non-stationary operation conditions［J］.Mechanical Systems and Signal Processing，2009，23：1 528—1 534.

［6］ 馮志鵬，范寅夕，Liang Ming，等.行星齒輪箱故障診斷的非平穩(wěn)振動信號分析方法［J］.中國電機工程學報，2013，33（17）：105—110.Feng Zhipeng，F(xiàn)an Yinxi，Liang Ming，et al.A nonstationary vibration signal analysis method for fault diagnosis of planetary gearboxes［J］.Proceedings of the CSEE，2013，33（17）：105—110.

［7］ Feng Zhipeng，Liang Ming，Chu Fulei.Recent advances in time-frequency analysis methods for machinery fault diagnosis：A review with application examples［J］.Mechanical Systems and Signal Processing，2013，38：165—205.

［8］ Kaiser J F.On a simple algorithm to calculate the ‘energy’of a signal［A］.Proceedings of IEEE International Conference on Acoustics，Speech and Signal Processing［C］.Albuquerque，USA.1990，1：381—384.

［9］ Maragos P，Kaiser J F，Quatieri T F.Energy separation in signal modulations with application to speech analysis［J］.IEEE Transactions on Signal Processing，1993，41（10）：3 024—3 051.

［10］Maragos P，Potamianos A.Higher order differential energy operators［J］.IEEE Signal Processing，1995，2 （8）：152—154.

［11］Diop E H S，Boudraa A O，Salzenstein F.A joint 2D AM-FM estimation based on higher order Teager-Kaiser energy operators［J］.Signal Image and Video Processing，2011，5：61—68.

［12］Salzenstein F，Boudraa A O，Chonavel T.A new class of multi-dimensional Teager-Kaiser and higher order operators based on directional derivatives［J］.Multidimensional Systems and Signal Processing，2013，24 （3）：543—572.

［13］Vold H，Herlufsen H.Order tracking at extreme slew rates，using Kalman tracking filters［J］.SAE Paper Number 931288，1993.

［14］Vold H，Mains M，Blough J.Theoretical foundations for high performance order tracking with the Vold-Kalman tracking filter［J］.SAE Paper Number 972007，1997.

［15］Vold H，Herlufsen H，Mains M.Multi axle order tracking with the Vold-Kalman tracking filter［J］. Journal of Sound and Vibration，1997，13（5）：30—34.

［16］Feldbauer Ch，Holdrich R.Realization of a Vold-Kalman tracking filter—A least squares problem［A］. Proceedings of the COST G-6 Conference on Digital Audio Effects（DAFX-000）［C］.Verona Italy，December 7-9，2000：1—4.

［17］Tuma J.Setting the pass bandwidth in the Vold-Kalman order tracking filter［A］.Twelfth International Congress on Sound and Vibration［C］.Lisbon.2005，Paper 719.

Application of Vold-Kalman filter and higher order energy separation to fault diagnosis of planetary gearbox under time-varying conditions

QIN Si-feng1，F(xiàn)ENG Zhi-peng1，LIANG Ming2
（1.School of Mechanical Engineering，University of Science and Technology Beijing，Beijing 100083，China；2.Department of Mechanical Engineering，University of Ottawa，Ottawa，ON K1N 6N5，Canada）

The vibration signals of planetary gearboxes under nonstationary running conditions have significant time-varying modulation feature.Conventional spectral analysis methods are unable to identify the gear fault characteristic frequencies from such nonstationary signals.A method based on Vold-Kalman filter and higher order energy separation is proposed to analyze the vibration signals of planetary gearboxes under nonstationary conditions in time-frequency domain，thus to identify the timevarying characteristic frequencies and diagnose the gear faults.The time-frequency representation derived from Vold-Kalman filter and higher order energy separation provides nice time-frequency resolution and is free from cross term interference，and thus it can effectively pinpoint the time-varying constituent frequency components of nonstationary signals.The proposed method is validated with analysis of lab experimental signals of a planetary gearbox under time-varying running conditions.

fault diagnosis；planetary gearbox；higher order energy separation；Vold-Kalman filter；time-frequency analysis

TH165+.3；TH132.425

A

1004-4523（2015）05-0839-07

10.16385/j.cnki.issn.1004-4523.2015.05.20

秦嗣峰（1989—），男，碩士研究生。電話：18866226821；E-mail：qsf0312@163.com

馮志鵬（1973—），男，博士，教授，博士生導師。電話：13520978607；E-mail：fengzp@ustb.edu.cn

2014-05-13<； class="emphasis_bold">；修訂日期：2；

2014-11-05

國家自然科學基金資助項目（11272047，51475038）和教育部新世紀優(yōu)秀人才支持計劃項目（NCET-12-0775）