基于經(jīng)驗(yàn)?zāi)Ｊ椒纸夂团帕徐氐妮S承故障特征提取

王 濤，胡定玉，丁亞琦，廖愛華，師 蔚

（1.上海工程技術(shù)大學(xué) 城市軌道交通學(xué)院，上海201620；2.上海地鐵維護(hù)保障有限公司 車輛分公司，上海200235）

滾動(dòng)軸承作為減少摩擦損失的精密元件，被廣泛用于旋轉(zhuǎn)機(jī)械中，但高速重載的環(huán)境導(dǎo)致軸承極易發(fā)生故障，影響列車的安全出行，因此對(duì)軸承的運(yùn)行狀況檢測十分必要。當(dāng)滾動(dòng)軸承出現(xiàn)故障時(shí)，內(nèi)外滾道與滾動(dòng)部件之間發(fā)生碰撞，會(huì)產(chǎn)生一個(gè)強(qiáng)迫短脈沖信號(hào)，通過對(duì)這些信號(hào)進(jìn)行分析可以有效檢測故障從而減少事故的發(fā)生。與傳統(tǒng)的振動(dòng)分析相比，采用聲信號(hào)監(jiān)測軸承的健康狀況更加方便快捷，無需粘貼傳感器，與檢測對(duì)象相互獨(dú)立，該方法對(duì)軸承早期故障的檢測也有較好的效果。

傳統(tǒng)的時(shí)域與頻域分析方法難以有效分析軸承故障信號(hào)的時(shí)變特性，無法準(zhǔn)確體現(xiàn)其故障特征。針對(duì)該類信號(hào)，通常采用小波分解[1]和經(jīng)驗(yàn)?zāi)Ｊ椒纸鈁2]等時(shí)頻分析方法。時(shí)頻分析將時(shí)頻兩域結(jié)合起來，可以同時(shí)考慮到兩個(gè)域的性能；在時(shí)頻相平面上可以精確地定位出某一時(shí)刻出現(xiàn)了哪些頻率分量。小波分解在分析信號(hào)之前設(shè)置一些參數(shù)（小波基函數(shù)），將信號(hào)分解成不同頻率范圍的子信號(hào)[3-4]，不同的小波基會(huì)產(chǎn)生不同的分解結(jié)果，最佳小波基的選擇成為了較大的難題。經(jīng)驗(yàn)?zāi)B(tài)分解（EMD）[5]可以自適應(yīng)地將信號(hào)分解為固有模態(tài)函數(shù)和殘差，可以對(duì)不同的分量進(jìn)行獨(dú)立分析[6]。

目前EMD 已在滾動(dòng)軸承故障診斷領(lǐng)域廣泛應(yīng)用，楊世錫等[7]將經(jīng)驗(yàn)?zāi)Ｊ椒纸夂虷ilbert 變換引入旋轉(zhuǎn)機(jī)械的故障診斷中，并與小波時(shí)頻分析相比較，研究結(jié)果表明基于EMD 的Hilbert 變換具有更好的效果，但出現(xiàn)了一些虛假分量；張志剛等[8]提出一種改進(jìn)經(jīng)驗(yàn)?zāi)B(tài)分解和譜峭度相結(jié)合的算法，剔除了分解過程中出現(xiàn)的虛假分量，但缺點(diǎn)是其受噪聲影響較大。

本文提出一種基于經(jīng)驗(yàn)?zāi)Ｊ椒纸夂团帕徐氐妮S承故障分析的信號(hào)處理方法，為EMD分解后IMF分量的選擇提供了新的方法。首先介紹相關(guān)基礎(chǔ)理論，然后通過數(shù)值仿真與傳統(tǒng)的包絡(luò)分析進(jìn)行比較，最后通過實(shí)驗(yàn)數(shù)據(jù)的分析驗(yàn)證該方法的有效性。

1 基礎(chǔ)理論

1.1 經(jīng)驗(yàn)?zāi)Ｊ椒纸馑惴?/h3>

EMD是一種自適應(yīng)時(shí)頻信號(hào)處理方法，通過對(duì)復(fù)雜信號(hào)進(jìn)行分解，獲得一系列用于表征原始信號(hào)特性的IMF 分量，篩選一些合適的本征模態(tài)分量實(shí)現(xiàn)對(duì)原始信號(hào)所包含的信息特征進(jìn)行分析。分解出的IMF 必須滿足以下兩個(gè)條件[9-11]：(1)在所需分解的數(shù)據(jù)段中，極點(diǎn)數(shù)和過零點(diǎn)數(shù)的差值不能超過一個(gè)；(2)信號(hào)的任意時(shí)刻，由極大值點(diǎn)插值所構(gòu)造的上包絡(luò)線和極小值點(diǎn)構(gòu)成的下包絡(luò)線均值為零。

EMD的實(shí)施步驟如下：

(1)找出原始信號(hào)所有的極大值點(diǎn)和極小值點(diǎn)，利用插值的方法分別對(duì)極大值點(diǎn)和極小值點(diǎn)進(jìn)行連接來擬合原始信號(hào)的波形，形成上包絡(luò)線和下包絡(luò)線。

(2)計(jì)算上下包絡(luò)線的平均值m1，求得：

若h1滿足上述條件，則為信號(hào)的一個(gè)IMF分量，記為IMF1，進(jìn)入下一個(gè)循環(huán)。

(3)若h1不滿足條件，則用h1替代初始數(shù)據(jù)x(t)，重復(fù)步驟(a)、(b)，得到h2，并判斷h2是否滿足(b)的條件，若還是不滿足條件，則繼續(xù)循環(huán)迭代直到獲得滿足條件的IMF1分量。

(4)將IMF1分量從信號(hào)x(t)中分離出來，得到新的循環(huán)信號(hào)：

(5)將R1作為原始信號(hào)，重復(fù)步驟(1)～(5)，當(dāng)Rn的極值點(diǎn)數(shù)小于2，則循環(huán)結(jié)束，得到n個(gè)IMF 分量和殘余分量Rn，信號(hào)x(t)可表示為

在EMD的分解過程中，由于需要人為設(shè)置迭代次數(shù)和終止條件，因此會(huì)產(chǎn)生一些無關(guān)信息成分，即虛假分量。針對(duì)該問題，可以利用相關(guān)系數(shù)來識(shí)別虛假分量[12]。

1.2 排列熵算法

排列熵作為描述復(fù)雜系統(tǒng)的算法，對(duì)突變的信號(hào)較為敏感。當(dāng)滾動(dòng)軸承系統(tǒng)發(fā)生故障時(shí)，聲信號(hào)包含的沖擊成分，會(huì)對(duì)排列熵值造成影響，因此可以利用排列熵對(duì)軸承進(jìn)行早期故障診斷。

給定一個(gè)序列{X(K),K=1,2,…,N}，對(duì)其進(jìn)行相空間重構(gòu)得：

其中：d為嵌入維數(shù)，τ為延遲時(shí)間，k為重構(gòu)分量。對(duì)矩陣Z中的每一行序列進(jìn)行升序排列：

式中：t表示索引的列數(shù)，R1,R2,…,Rd為X(K)中各元素所處的位置。定義xdi為重構(gòu)后的任意一組序列，0 ≤Ri≤d!，對(duì)于d維相空間映射有d!種排列的可能性，計(jì)算每一種序列出現(xiàn)的可能性P1,P2,…,Pi，對(duì)于時(shí)間序列X(K)的i種排列熵Hp(d)：

當(dāng)Pi=m!時(shí)，Hp(d)就達(dá)到最大值In(m!)，為了方便熵值的比較，往往會(huì)對(duì)排列熵值進(jìn)行歸一化處理：

Hp的值越小則表示該序列越有序，反之該序列的復(fù)雜程度越高。由于發(fā)生故障的軸承元件在運(yùn)行過程中與其他元件發(fā)生周期性的沖擊，相較于正常軸承信號(hào)，故障信號(hào)更加有序，排列熵值較小。因此可以用排列熵值的突變來判斷軸承是否發(fā)生故障。

在計(jì)算排列熵時(shí)，插入維數(shù)d過小會(huì)導(dǎo)致重構(gòu)信號(hào)狀態(tài)量缺失，對(duì)沖擊信號(hào)的敏感性將大大下降，算法失去有效性；插入維數(shù)過大會(huì)使得信號(hào)在相空間的重構(gòu)過程中均勻化時(shí)間序列，因而一般維數(shù)選擇為3～7，本文取5，時(shí)間延遲τ=1[13]。

EMD將信號(hào)分解為多個(gè)IMF分量，高頻IMF分量中包含較多沖擊成分，但包含大量噪聲，低頻IMF分量噪聲成分較少，但可能是虛假分量。為篩選本征模態(tài)分量中包含信息量較大的真實(shí)分量進(jìn)行信號(hào)的重構(gòu)，可以利用排列熵值和相關(guān)系數(shù)來作為篩選指標(biāo)。排列熵值Hp與原始信號(hào)的有效成分呈反相關(guān)，相關(guān)系數(shù)r與原始信號(hào)的線性相關(guān)程度呈正相關(guān)。基于此本文提出聯(lián)合系數(shù)Co作為IMF 分量的篩選指標(biāo)，取聯(lián)合系數(shù)最大的前兩個(gè)本征模態(tài)分量進(jìn)行信號(hào)的重構(gòu)：

1.3 快速譜峭度算法

譜峭度為能量歸一化的4 階譜累積量，利用計(jì)算每根譜線的峭度值來反映瞬時(shí)沖擊的大小，快速譜峭度(Fast-kurtogram)則以STFT作為窗口，根據(jù)快速譜峭度圖選擇最佳濾波參數(shù)。

在信號(hào)x(t)的激勵(lì)下系統(tǒng)響應(yīng)y(t)Wold-cramer分解的頻域表達(dá)式為

其中：H(t,f)為y(t)在頻率f處的復(fù)包絡(luò)信號(hào)。y(t)的四階譜累積量為

式中：S2ny(f)為2n階譜瞬時(shí)距，用以表示復(fù)包絡(luò)信號(hào)的能量大小，定義為

式中，df表示對(duì)頻率f求導(dǎo)。譜峭度定義為

為減少運(yùn)算時(shí)間，可以應(yīng)用于實(shí)際的工程測量，采用塔式快速譜峭度對(duì)信號(hào)進(jìn)行2倍頻采樣濾波迭代，將原信號(hào)進(jìn)行L層濾波分解，每層濾波后的信號(hào)頻段為2L個(gè)，一般情況下L為整數(shù)。為了提高信號(hào)分解的精度，在L+1 層和L+2 層中插入3×2L個(gè)濾波器，即L+1.6層，如當(dāng)L=0時(shí)，第1層和第2層之間分解層數(shù)為1.6層時(shí)，表明該層將原始信號(hào)分解為3個(gè)信號(hào)頻段。根據(jù)快速譜峭度圖上峭度值的大小，即顏色最深的頻段選取最佳濾波中心和濾波帶寬，進(jìn)行平方包絡(luò)分析和故障提取。

2 EMD-PE故障診斷方法流程

應(yīng)用快速譜峭度算法診斷軸承故障雖然取得了一定的效果，但對(duì)信號(hào)進(jìn)行預(yù)處理可以進(jìn)一步提高診斷效果，本文將EMD 和排列熵值相結(jié)合，提出一種基于聲信號(hào)滾動(dòng)軸承故障診斷的EMD-PE 方法，可以有效降低噪聲的干擾并確定合適的濾波器參數(shù)。該方法流程圖如圖1所示。具體過程如下。

(1)對(duì)原始信號(hào)進(jìn)行EMD分解，把非平穩(wěn)、非線性的聲信號(hào)分解為8 個(gè)本征模態(tài)分量IMF1，IMF2，···，IMF8。

圖1 EMD-PE診斷故障軸承流程圖

(2)計(jì)算各個(gè)本征模態(tài)分量的排列熵值Hp1，Hp2，··· ，Hp8(0 ＜Hpi＜1)和相關(guān)系數(shù)r1，r2，··· ，r8(0 ＜ri＜1)；

(3)計(jì)算聯(lián)合系數(shù)，根據(jù)聯(lián)合系數(shù)最大原則，篩選系數(shù)值最大的前兩個(gè)本征模態(tài)分量進(jìn)行信號(hào)的重構(gòu)。

(4)對(duì)重構(gòu)信號(hào)進(jìn)行快速譜峭度分析，根據(jù)峭度最大化原則選擇最佳濾波中心和濾波帶寬。

(5)將從中提取到的軸承特征頻率與理論故障特征頻率對(duì)比，進(jìn)行故障診斷。

3 軸承故障信號(hào)仿真分析

為了證明該方法在滾動(dòng)軸承故障診斷中的有效性，現(xiàn)對(duì)模擬的軸承故障仿真信號(hào)進(jìn)行分析，設(shè)軸承發(fā)生單點(diǎn)局部故障，傳聲器接收到信號(hào)的信號(hào)數(shù)學(xué)模型為

其中：x(t)為軸承故障仿真信號(hào)，A為沖擊信號(hào)幅值，s(t)為脈沖沖擊函數(shù)，Ti為特征頻率出現(xiàn)的周期，fs為軸承座、傳感器系統(tǒng)的高頻固有頻率，ξ為幅值衰減系數(shù)，n(t)為高斯白噪聲信號(hào)。

設(shè)仿真參數(shù)A=1，特征頻率fi=1/Ti=130 Hz，固有頻率fs=5 000 Hz，衰減系數(shù)ξ=1000，采樣頻率Fs=51200 Hz，采樣點(diǎn)數(shù)N=10 240。仿真脈沖信號(hào)如圖2(a)所示，對(duì)仿真信號(hào)添加高斯白噪聲，SNR=10，如圖2(b)所示。

對(duì)含有高斯白噪聲的仿真信號(hào)進(jìn)行EMD分解，對(duì)分解后的各個(gè)IMF分量進(jìn)行排列熵值與相關(guān)系數(shù)的計(jì)算，并根據(jù)式(8)求得聯(lián)合系數(shù)，如表1所示。根據(jù)聯(lián)合系數(shù)最大化原則，篩選出兩個(gè)最佳的IMF 分量(IMF6，IMF7)進(jìn)行信號(hào)的重組，圖3顯示的是重組信號(hào)進(jìn)行快速譜峭度分析得到快速譜峭度圖。由圖3可以看出，分解層數(shù)為1.6，帶寬為8 333.33 Hz，濾波中心為23 333.33 Hz時(shí)，譜峭度值達(dá)到最大。在濾波后的平方包絡(luò)圖中，準(zhǔn)確檢測出了130 Hz 的特征頻率及其倍頻程如圖4所示，這與預(yù)設(shè)的故障特征頻率相符合，于是實(shí)現(xiàn)了故障的診斷，說明基于EMD 和排列熵的快速譜峭度的軸承故障診斷是有效的。

表1 仿真數(shù)據(jù)

圖2 仿真信號(hào)時(shí)域圖

圖3 仿真信號(hào)快速譜峭度圖

圖4 基于EMD-PE算法的軸承故障仿真試驗(yàn)結(jié)果

4 實(shí)驗(yàn)分析

實(shí)驗(yàn)采用SKF 公司生產(chǎn)的滾動(dòng)軸承故障實(shí)驗(yàn)臺(tái)，故障軸承型號(hào)為SKF6016深溝球軸承,內(nèi)圈存在一條裂紋，深約1 mm，其他部位保持完好，表2為該軸承的參數(shù)。采樣頻率是51 200 Hz，采樣時(shí)間為10s，電機(jī)轉(zhuǎn)速為404 r/min，轉(zhuǎn)頻為7 Hz，負(fù)載為5 N。

表2 滾動(dòng)軸承參數(shù)

對(duì)于外圈固定，內(nèi)圈轉(zhuǎn)動(dòng)的滾動(dòng)軸承，其內(nèi)圈故障的特征頻率為

將軸承參數(shù)代入式(14)中可得軸承內(nèi)圈故障特征頻率為55.69 Hz。

聲信號(hào)的希爾伯特包絡(luò)變換的頻譜圖如圖5所示，Hilbert頻譜圖中雖然可以得到故障頻率，但卻不夠明顯。

對(duì)信號(hào)進(jìn)行EMD 分解，得到12 個(gè)本征模態(tài)分量，計(jì)算各個(gè)分量的排列熵值、相關(guān)系數(shù)以及聯(lián)合系數(shù)，如表3所示，由于篇幅有限，這里只顯示前8 個(gè)IMF的計(jì)算結(jié)果，其中原始信號(hào)的排列熵值為0.77。

由表3可知IMF4和IMF5的聯(lián)合系數(shù)Co較大，原始信號(hào)中的突變特性成分較多，對(duì)IMF4和IMF5分量進(jìn)行重組，得到的重組信號(hào)的排列熵值為0.29，相比于原始故障信號(hào)有用信息比例有所提高。

圖5 轉(zhuǎn)速為404r min時(shí)聲信號(hào)的希爾伯特包絡(luò)變換頻譜圖

表3 各IMF分量排列熵值、相關(guān)系數(shù)及聯(lián)合系數(shù)

對(duì)合成信號(hào)進(jìn)行快速譜峭度分析得到快速譜峭度圖和平方包絡(luò)圖,如圖6(a)、圖6(b)所示。選取分解層數(shù)為4，帶通濾波器的濾波中心18 000 Hz，而濾波帶寬為1 200 Hz,在此濾波范圍內(nèi)峭度值達(dá)到最大，對(duì)濾波后的信號(hào)進(jìn)行平方包絡(luò)分析可以得到較為明顯的滾動(dòng)軸承內(nèi)圈故障頻率53.93 Hz 及其倍頻，證明了本文方法的有效性。

圖6 重構(gòu)信號(hào)的快速峭度譜圖和平方包絡(luò)譜圖

5 結(jié)語

本文提出了一種基于經(jīng)驗(yàn)?zāi)Ｊ椒纸夂团帕徐氐妮S承故障分析的信號(hào)處理方法，利用聯(lián)合系數(shù)代替了傳統(tǒng)的峭度值，對(duì)IMF 分量進(jìn)行篩選重構(gòu)。此外針對(duì)濾波參數(shù)較難確定的問題，譜峭度對(duì)聲信號(hào)中瞬態(tài)沖擊較為敏感，可以自動(dòng)確定帶通濾波器的參數(shù)。經(jīng)驗(yàn)?zāi)Ｊ椒纸夂团帕徐叵嘟Y(jié)合的滾動(dòng)軸承故障診斷方法，能夠有效提取出故障特征頻率及其倍頻。仿真及實(shí)驗(yàn)結(jié)果驗(yàn)證了該方法的有效性。