基于ITD改進算法和關(guān)聯(lián)維數(shù)的轉(zhuǎn)子故障診斷方法

楊 宇，王歡歡，喻鎮(zhèn)濤，袁 輝

(1.湖南大學(xué) 汽車車身先進設(shè)計制造國家重點實驗室，長沙 410082;2.東風(fēng)汽車公司 技術(shù)中心試驗部，襄陽 441004)

轉(zhuǎn)子系統(tǒng)振動信號往往會呈現(xiàn)一定程度的非線性特征，當(dāng)轉(zhuǎn)子發(fā)生故障時，其振動信號會混有能量較大的與轉(zhuǎn)速有關(guān)的背景信號和噪聲，傳統(tǒng)的時域或頻域方法很難從這樣的轉(zhuǎn)子振動信號中有效提取故障信息[1]，而分形理論中的分形維數(shù)作為定量刻畫混沌吸引子的一個重要參數(shù)，在研究非線性系統(tǒng)方面具有獨特的優(yōu)勢。其中關(guān)聯(lián)維數(shù)計算簡單，便于從數(shù)據(jù)序列中直接測定，因而應(yīng)用廣泛，并已應(yīng)用于機械故障診斷之中[2]，故可以用來提取轉(zhuǎn)子系統(tǒng)的故障特征。相關(guān)研究表明，關(guān)聯(lián)維數(shù)在較大噪聲水平下難以用關(guān)聯(lián)維數(shù)來準(zhǔn)確刻畫系統(tǒng)的工作狀態(tài)和故障類型[3]，因此要提高計算精度就必須對信號進行降噪處理;而當(dāng)信號的非線性特性突出時，利用線性濾波等傳統(tǒng)方法進行降噪的效果并不明顯[4]。采用 EMD(Empirical Mode Decomposition) 方 法[5]和 LMD (Local Mean Decomposition)方法[6]都可以實現(xiàn)對振動信號的自適應(yīng)濾波，但這兩種方法在理論上都還存在一些問題，如EMD的過包絡(luò)、欠包絡(luò)和頻率混淆問題，LMD也存在由于解調(diào)而引起的信號突變的問題。目前這些問題仍在研究當(dāng)中[7]。為了降低轉(zhuǎn)子振動信號的噪聲，本文采用內(nèi)稟時間尺度分解(ITD)改進算法對轉(zhuǎn)子振動信號進行分解以達到降噪的目的。

ITD方法是由Frei提出的另一種自適應(yīng)時頻分析方法[8]，該方法自適應(yīng)地將非平穩(wěn)信號分解成為若干瞬時頻率具有物理意義的單分量信號(固有旋轉(zhuǎn)分量)。文獻[8]雖然將ITD方法和EMD方法進行了對比，分析了ITD方法在端點效應(yīng)和計算速度上的優(yōu)勢，但是并沒有對ITD方法及固有旋轉(zhuǎn)分量的物理意義進行闡述。因此本文在闡述ITD方法及其固有旋轉(zhuǎn)分量的物理意義的基礎(chǔ)上，提出了基于ITD的改進算法，并利用ITD改進方法對轉(zhuǎn)子振動信號進行分解，將包含主要故障信息的若干分量進行組合作為降噪后的轉(zhuǎn)子振動信號，并對降噪前后轉(zhuǎn)子系統(tǒng)5種狀態(tài)下的信號分別計算了關(guān)聯(lián)維數(shù)，對比發(fā)現(xiàn)，降噪后的關(guān)聯(lián)維數(shù)區(qū)分性明顯增強，因此可以用該方法識別轉(zhuǎn)子的工作狀態(tài)和故障類型。

1 ITD改進方法

1.1 ITD方法

ITD方法將待分析的非平穩(wěn)信號分解成一系列的固有旋轉(zhuǎn)分量和一個單調(diào)的趨勢信號。設(shè)Xt是待分析的原信號，分解前先定義一個基線提取算子L，使得從原始信號中去掉該基線后剩下的余量信號成為一個固有旋轉(zhuǎn)分量。一次分解的表達式為[8]:

式中假設(shè)在[0，τk]上定義了 Lt和 Ht，而 Xt在[0，τk+2]上有定義。在連續(xù)極值點[τk，τk+2]區(qū)間上定義 Xt的基線提取因子L，即:

式中:a是分解時的增益控制參數(shù)，0＜a＜1;Lt保留了信號在各個極值點的單調(diào)性，Ht提取了各個極值點之間疊加的局部高頻分量信號——固有旋轉(zhuǎn)分量。重復(fù)分解過程就可得到一系列固有旋轉(zhuǎn)分量和單調(diào)趨勢信號。

給出如式(4)所示的仿真信號，它由一個調(diào)幅調(diào)頻信號和一個余弦信號組成:

圖1是該仿真信號的波形，圖2是運用ITD算法對信號進行分解得到的結(jié)果。

圖1 仿真信號x(t)Fig.1 Simulation signal x(t)

圖2 仿真信號x(t)的ITD分解結(jié)果Fig.2 Decomposition result of simulation signal x(t)by ITD

從圖2可以看出，分解得到的第1個分量效果較好，但是第2個分量出現(xiàn)了失真，這是因為ITD方法采用線性變換的方法獲得基線信號，使得波形出現(xiàn)了毛刺而失真。因此，為了克服該缺陷，本文提出了改進的ITD方法，即采用三次樣條插值代替了ITD方法中的線性變換。雖然采用了三次樣條插值，但它又不同于EMD中基于局部極值的包絡(luò)均值，因為ITD改進方法每次分解只需要進行一次三次樣條插值。

1.2 ITD 改進方法

文獻[8]中雖然給出了ITD方法的詳細(xì)分解過程，但是并沒有對該方法及分解得到的結(jié)果的物理意義進行闡述;同時經(jīng)過以上的研究發(fā)現(xiàn)，由于ITD方法采用線性變換的方法實現(xiàn)對信號的分解，從第二個分量開始，出現(xiàn)了較明顯的信號失真。因此，本節(jié)首先闡述了ITD分解方法的物理意義，接著提出了ITD改進算法。

EMD與LMD方法在提供有效的自適應(yīng)時頻分析方法的同時，提供了一條很好的自適應(yīng)時頻分析思路，即首先假設(shè)任意一個復(fù)雜信號由若干個瞬時頻率具有物理意義的單分量信號組成，并給出瞬時頻率具有物理意義的單分量信號所需要滿足的條件，然后據(jù)此條件對原始信號進行自適應(yīng)地分解。

在ITD方法中，對于理想的固有旋轉(zhuǎn)分量，其如式(2)所示的基線信號在所有的區(qū)間[τk，τk+2]上都應(yīng)該等于零，因此如式(3)所示的Lk+1應(yīng)該等于零。理論上ITD方法對信號進行分解得到的固有旋轉(zhuǎn)分量應(yīng)該滿足其基線信號控制點Lk+1等于零。本文在此基礎(chǔ)上定義了瞬時頻率具有物理意義的內(nèi)稟尺度分量(Intrinsic Scale Component，ISC)。所定義的ISC分量滿足以下條件:

(1)在整個數(shù)據(jù)段內(nèi)，任意兩個極值點之間呈單調(diào)性;

圖3 ISC需要滿足的條件Fig.3 The needed condition of ISC

(2) 在整個數(shù)據(jù)段內(nèi)，其極值點為 Xk，k=1，2，…，M，各個極值點相對應(yīng)的時刻為 τk，k=1，2，…，M，由任意兩個極大(小)值點(τk，Xk)、(τk+2，Xk+2)連接形成的線段在其中間極小(大)值點(τk+1，Xk+1)相對應(yīng)時刻 τk+1的函數(shù)值Xk)與該極小(大)值 Xk+1的比值關(guān)系不變，即，即滿足，如圖3 所示。其中，a∈(0，1)為一常量，典型地，a=1/2，例如正弦或余弦信號、調(diào)幅信號、調(diào)頻信號、調(diào)幅－調(diào)頻信號等。

以上兩個條件保證了ISC分量任意兩個相鄰極值點之間具有單一的模態(tài)，而且在局部(極值點與相鄰的零交叉點之間)近似吻合標(biāo)準(zhǔn)正弦曲線，因此瞬時頻率具有物理意義。

改進的ITD算法假設(shè)任何復(fù)雜信號是由不同ISC分量組成，任何兩個ISC分量之間相互獨立，這樣就可以將復(fù)雜信號分解為若干個ISC分量之和，具體步驟如下:

(1)確定 X(t)的所有極值點 Xk，k=1，2，…，M 及其相對應(yīng)的時刻 τk，k=1，2，…，M，并設(shè)置參數(shù) a，計算各基線控制點Lk，即式(3);然后對所有Lk進行三次樣條插值，得到基信號線L1。

(2)將L1從原始信號中分離出來，得到P1。理想地，P1為一個ISC分量，則P1為信號x(t)的第1個分量。理想的ISC分量應(yīng)滿足Lk+1等于零，實際中可以設(shè)定變動量Δ，當(dāng)時迭代結(jié)束。

(3)如P1不滿足ISC的條件，則將P1作為原始信號重復(fù)步驟(1)和(2)，循環(huán)k次，直到得到內(nèi)稟尺度分量Pk，Pk即為信號x(t)的第1個分量ISC1。

(4)將ISC1從x(t)中分離出來，得到一個新的信號r1，將r1作為原始信號重復(fù)步驟(1)、(2)和(3)，得到x(t)的第二個滿足ISC條件的分量ISC2。重復(fù)循環(huán)n次，得到信號x(t)的n個滿足ISC條件的分量，直到rn為一單調(diào)函數(shù)為止。這樣便可以將x(t)分解為n個內(nèi)稟尺度分量ISC和一個單調(diào)函數(shù)rn之和，即:

ISC迭代終止判據(jù)Δ的選擇對信號的分解和迭代次數(shù)都有重要的影響，類似于EMD方法內(nèi)稟模態(tài)函數(shù)(IMF)的判據(jù)。常用的判據(jù)方法有，Huang在文獻[9]中介紹的標(biāo)準(zhǔn)偏差法，F(xiàn)landrin在文獻[10]中介紹的三參數(shù)法等。其中標(biāo)準(zhǔn)偏差法與IMF分量的定義無關(guān)。因此，本文在考慮ISC分量定義的基礎(chǔ)上，結(jié)合三參數(shù)法，來確定變動量Δ的取值。

經(jīng)過仿真信號分解實驗發(fā)現(xiàn)，當(dāng)基線控制點Lk+1與信號的對應(yīng)極值點Xk+1之比，小于一定取值范圍，即Δ取值范圍介于0.005時，則得到的分量可以近似為一個ISC分量，此時迭代次數(shù)一般為4～10次。

為了減少ITD改進算法中邊界效應(yīng)的影響，在處理數(shù)據(jù)之前，需對原始數(shù)據(jù)進行延拓處理，處理方法可采用G.Rilling程序的鏡像對稱延拓方法[11]。

在處理了端點效應(yīng)之后，對式(4)所示的仿真信號采用改進ITD算法進行分解，分解時取a=0.5，迭代終止判據(jù)，分解后得到如圖4所示的兩個分量和一個余量，兩次的迭代次數(shù)分別為5次和4次。采用EMD方法對該仿真信號進行分解，在延拓方法和迭代終止判據(jù)相同的情況下，分解結(jié)果如圖5所示，兩次的迭代次數(shù)分別為17次和105次。顯然，ITD改進算法的迭代次數(shù)要少得多，因而分解速率更快。

圖4 仿真信號的改進ITD分解結(jié)果(已處理端點效應(yīng))Fig.4 Decomposition result of simulation signal by improved ITD(with boundary processing)

圖5 仿真信號的EMD分解結(jié)果(已處理端點效應(yīng))Fig.5 Decomposition result of simulation signal by EMD(with boundary processing)

圖6 是由電渦流位移傳感器測得的不平衡故障狀態(tài)下的轉(zhuǎn)子徑向位移振動信號，對其采用ITD改進算法進行分解，共得到3個基本分量和一個余量，如圖7所示，從圖中可以看出，第1個分量頻率成分最高，第2個分量主要為與轉(zhuǎn)速有關(guān)的信號，因此故障信息主要包含在第1和第2個ISC分量中，可以選擇這兩個分量進行組合，得到降噪后的轉(zhuǎn)子信號。

圖6 不平衡故障狀態(tài)下的轉(zhuǎn)子徑向位移振動信號Fig.6 Rotor radial displacement vibration signal with unbalance fault

圖7 不平衡故障狀態(tài)下的轉(zhuǎn)子徑向位移振動信號的改進ITD分解結(jié)果Fig.7 Decomposition result of rotor vibration signal under state of unbalance by improved ITD

2 實驗分析

本文采用南京東大測振儀器廠生產(chǎn)的ZT－3型轉(zhuǎn)子振動模擬實驗臺進行轉(zhuǎn)子系統(tǒng)故障實驗，實驗裝置如圖8所示。實驗時鍵相傳感器采用電渦流傳感器，可以提供相位和轉(zhuǎn)速信號;轉(zhuǎn)子徑向位移振動信號由垂直和水平安裝的電渦流傳感器拾取，信號經(jīng)過信號調(diào)理箱處理后，送入數(shù)據(jù)采集系統(tǒng)。實驗?zāi)M了不平衡、不對中、碰摩和油膜渦動四種故障狀態(tài)[12]，信號采集時對應(yīng)轉(zhuǎn)速為3 000 r/min，采樣率為4 096 Hz。

圖9～12分別是轉(zhuǎn)子在正常、不對中、碰摩、油膜渦動狀態(tài)下的徑向位移振動信號(不平衡狀態(tài)下的振動信號見圖6)。

圖8 轉(zhuǎn)子試驗臺裝置Fig.8 Rotor test rig

圖9 正常狀態(tài)下的轉(zhuǎn)子徑向位移振動信號時域波形Fig.9 Rotor radial displacement vibration signal under normal state

圖10 不對中故障的轉(zhuǎn)子徑向位移振動信號時域波形Fig.10 Rotor radial displacement vibration signal with misalignment fault

圖11 碰摩故障的轉(zhuǎn)子徑向位移振動信號時域波形Fig.11 Rotor radial displacement vibration signal with rub-impact fault

圖12 油膜渦動故障的轉(zhuǎn)子徑向位移振動信號時域波形Fig.12 Rotor radial displacement vibration signal with oil whirl fault

首先，原始振動信號不進行降噪處理，直接進行分形維數(shù)的計算，計算結(jié)果如圖13所示。考慮到不同嵌入維數(shù)下計算得到的關(guān)聯(lián)維數(shù)也不同，因此選取了一定范圍的嵌入維數(shù)，從圖中可以看出，各狀態(tài)對應(yīng)的分形維數(shù)相互交叉，因此對原始信號直接計算關(guān)聯(lián)維數(shù)并不能清楚地分辨出轉(zhuǎn)子的工作狀態(tài)和故障類型。

圖13 不同狀態(tài)轉(zhuǎn)子原始信號在不同嵌入維數(shù)下的關(guān)聯(lián)維數(shù)Fig.13 The correlation dimension of vibration signals under different embedding dimension

然后，采用ITD改進算法對轉(zhuǎn)子各個狀態(tài)的原始信號進行降噪處理。即對各信號進行分解，得到若干分量。研究發(fā)現(xiàn)，各狀態(tài)下轉(zhuǎn)子振動信號的故障信息主要集中在前兩個分量中，選取了不同狀態(tài)下轉(zhuǎn)子振動信號的前兩個分量重構(gòu)信號。進一步對這些重構(gòu)信號計算不同嵌入維數(shù)下的關(guān)聯(lián)維數(shù)，結(jié)果如圖14所示。從圖中可以看出，經(jīng)過ITD改進算法降噪后的各狀態(tài)轉(zhuǎn)子信號在一定嵌入維數(shù)范圍內(nèi)的關(guān)聯(lián)維數(shù)明確分離，可以明顯地區(qū)別各工作狀態(tài)和故障類型，可見轉(zhuǎn)子信號的改進ITD降噪結(jié)合關(guān)聯(lián)維數(shù)可以準(zhǔn)確識別轉(zhuǎn)子的工作狀態(tài)和故障類型。嵌入維數(shù)的選取恰當(dāng)與否對相空間的重構(gòu)質(zhì)量有重要的影響。如果選取過小，將無法容納系統(tǒng)的吸引子，會出現(xiàn)偽鄰近點;反之太大，不僅造成計算工作量大，還減少可使用數(shù)據(jù)長度，使所建相空間中的相點過于稀疏。關(guān)于嵌入維數(shù)的選取，實際工程應(yīng)用中，一般情況，缺少對系統(tǒng)分形維數(shù)的先驗知識。對于實際系統(tǒng)嵌入維數(shù)m的選擇比較困難，文獻[13]中給出了經(jīng)驗公式m≥D，本文中嵌入維數(shù)的范圍在5～35時各狀態(tài)關(guān)聯(lián)維數(shù)的區(qū)分效果明顯。

圖14 轉(zhuǎn)子信號經(jīng)改進ITD分解剔除噪聲和背景信號后在不同嵌入維數(shù)下的關(guān)聯(lián)維數(shù)Fig.14 The correlation dimension of vibration signals denoised by improved ITD under different embedding dimension

3 結(jié)論

本文將改進ITD算法和關(guān)聯(lián)維數(shù)相結(jié)合應(yīng)用于轉(zhuǎn)子的故障診斷中，先通過改進算法對轉(zhuǎn)子信號進行分解降噪，再計算關(guān)聯(lián)維數(shù)來提取故障特征，以識別轉(zhuǎn)子的不同工作狀態(tài)和故障類型。由于關(guān)聯(lián)維數(shù)受噪聲的影響較大，而轉(zhuǎn)子原始信號中往往含有噪聲信號和與轉(zhuǎn)速有關(guān)的背景信號，在求取關(guān)聯(lián)維數(shù)之前，通過ITD改進算法，可將轉(zhuǎn)子原始信號的故障主要信息提取出來，以剔除噪聲和背景信號，再計算關(guān)聯(lián)維數(shù)，可獲得較好的效果，從而為轉(zhuǎn)子系統(tǒng)故障診斷提供了一種新的方法。值得提出的是，ITD方法剛剛提出，還有許多理論問題需要研究和完善，如迭代終止條件、端點效應(yīng)等問題，因此對ITD方法的研究還不是很深入。盡管如此，本文所采用的先定義瞬時頻率具有物理意義的單分量信號，然后再據(jù)此進行自適應(yīng)信號分解的方法，開拓了時頻分析的新思路，對于信號分析和故障診斷都有較大的借鑒意義。

[1]McFadden P D.Window functions for calculating of the time domain averages of the vibration of the individual planet gears and sun gear in an epicyclical gearbox[J].Journal of Vibration and Acoustics，1994，116(1):179 －187.

[2]Alberto R，Maria M.A method for the correlation dimension estimation for on-line condition monitoring of large rotating machinery[J].Mechanical Systems and Signal Processing.2005，19:939 －954

[3]Argyris J.The influence of noise on the correlation dimension of chaotic attractors[J].Chaos Solitions ＆ Fractals，1998，9(3):343－361.

[4]徐美芳，韓 焱，董劍龍.數(shù)字射線圖像降噪中的線性濾波方法[J].無損檢測，2005，27(5):259 －263.XU Mei-fang，HAN Yan，DONG Jian-long.Linear filter methods for image denoising I digital radiography[J].Nondestructive Inspection，2005，27(5):259 －263.

[5]Huang N E，Shen Z，Long S R.A new view of nonlinear water waves:the Hilbert spectrum[J].Annu.Rev.Fluid Mech，1999，31:417－457.

[6]Smith J S.The local mean decomposition and its application to EEG perception data[J].Journal of the Royal Society Interface，2005，2(5):443－454.

[7]Lei Y G，He Z J，Zi Y Y.Application of the EEMD method to rotor fault diagnosis of rotating machinery[J].Mechanical Systems and Signal Processing，2009，23:1327－1338.

[8]Frei M G，Osorio I.Intrinsic time-scale decomposition:timefreqeucny-energy analysis and real-time filtering of nonstationary signals[J].Proceedings of the Royal Society，A，2007，463:321－342.

[9]Huang N E，Shen Z，Long S R.The empirical mode decomposition and the Hilbert spectrum for nonlinear nonstationary time Series Analysis[J]. Proc. Roy. Soc.London，1998，454:903 －995.

[10]Rilling G，F(xiàn)landrin P.On empirical mode decomposition and its algorithms[J].IEEE Sig.Pro.Lett，2003:326 －329.

[11]Rilling G，F(xiàn)landrin P，Goncalves P.On empirical mode decomposition and its algorithms[J].In IEEE-EURASIP Workshop on Nonlinear Signaland Image Processing，Grado，Italy，June 8－11，2003.

[12]趙利華.旋轉(zhuǎn)機械故障機理與軸心軌跡識別方法研究[D].大連:大連理工大學(xué)，2010.

[13]Maestrello L.Chaotic response of panel vibrations forced by turbulent boundary layer and sound[J].AIAA Journal，1999，37(3):289 －294.