基于聽覺模型和極值點(diǎn)概率密度的斷齒故障特征提取方法研究

吳文壽， 李允公， 王 波， 李國萌， 石悅紅

(東北大學(xué) 機(jī)械工程與自動化學(xué)院，沈陽 110819)

?

基于聽覺模型和極值點(diǎn)概率密度的斷齒故障特征提取方法研究

吳文壽， 李允公， 王 波， 李國萌， 石悅紅

(東北大學(xué) 機(jī)械工程與自動化學(xué)院，沈陽 110819)

齒輪斷齒故障的重要特征是嚙合過程中在斷齒處產(chǎn)生碰撞與沖擊?？紤]到人耳聽覺系統(tǒng)對于突發(fā)的瞬態(tài)聲信號具有本能的反應(yīng)，為提取斷齒故障誘發(fā)的瞬態(tài)沖擊響應(yīng)成分，提出一種基于聽覺模型和信號極值點(diǎn)概率密度的特征提取方法。該方法首先對信號進(jìn)行GT帶通濾波、相位調(diào)整及極值點(diǎn)提取，然后計(jì)算各極值點(diǎn)的幅值概率密度，通過對其求導(dǎo)判斷各濾波通道中是否存在瞬態(tài)沖擊成分，繼而提取與之相關(guān)的極值點(diǎn)。同時(shí)，由于系統(tǒng)振動時(shí)會產(chǎn)生與斷齒沖擊無關(guān)的極值點(diǎn)，為準(zhǔn)確提取斷齒沖擊，根據(jù)瞬態(tài)信號頻帶連續(xù)性和多頻段分布特點(diǎn)，設(shè)計(jì)了相應(yīng)的提取方法。經(jīng)實(shí)測信號驗(yàn)證表明，所提方法能準(zhǔn)確刻畫及提取斷齒故障特征，可以在含有多種類型的瞬態(tài)沖擊響應(yīng)成分中提取出只由斷齒故障所誘發(fā)的沖擊成分，且提取結(jié)果精確度較高。

斷齒故障；聽覺模型；概率密度；瞬態(tài)信號；故障診斷；特征提取

齒輪傳動是一種最為常見的傳動形式。據(jù)統(tǒng)計(jì)，在機(jī)械設(shè)備故障或狀態(tài)不良的原因中，約60%[1]與齒輪失效有關(guān)，而輪齒斷裂在齒輪失效形式中所占比例又最高，約為41%[2]。斷齒使齒輪副不能正常嚙合，并產(chǎn)生周期性沖擊響應(yīng)成分。斷齒故障的動力學(xué)特性十分豐富[3]，但最本質(zhì)、最能區(qū)別于其他齒輪故障特征的在于斷齒會誘發(fā)出周期性瞬態(tài)沖擊響應(yīng)成分[4]，且斷齒越嚴(yán)重，沖擊越明顯。因此，如果能提取斷齒故障的沖擊響應(yīng)成分，并且能夠識別出沖擊響應(yīng)出現(xiàn)的周期，則可對齒輪斷齒故障診斷提供有效依據(jù)。在瞬態(tài)沖擊成分提取方面，栗茂林等[5]提出一種對連續(xù)小波系數(shù)進(jìn)行非線性化，實(shí)現(xiàn)沖擊成分特征提取。WANG等[6]提出一種基于Levenberg-Marquardt方法的瞬態(tài)模型和參數(shù)識別迭代提取方法,并最終應(yīng)用于軸承和齒輪故障特征提取。嚴(yán)?？档萚7]提出一種利用形態(tài)提升方法，通過放大微弱沖擊成分來實(shí)現(xiàn)微沖擊特征提取，并成功運(yùn)用到滾動軸承早期微弱脈沖故障檢測中。近年以來，利用譜峭度[8]方法也受到學(xué)者深入研究，主要利用譜峭度對瞬態(tài)沖擊成分較為敏感，依此來發(fā)現(xiàn)沖擊成分所占據(jù)的頻段，以便確定濾波器的帶寬和中心頻率，此方法對于提取瞬態(tài)沖擊信號表現(xiàn)良好。

人耳聽覺系統(tǒng)對突發(fā)的瞬態(tài)聲音信號具有本能的反應(yīng)[9]，基于此，李允公等[10]提出一種基于極值點(diǎn)概率密度和聽覺模型的瞬態(tài)信號提取方法，該方法從瞬態(tài)信號導(dǎo)致的概率密度曲線的長拖尾現(xiàn)象出發(fā)，并結(jié)合聽覺模型提出相應(yīng)的存在瞬態(tài)信號判斷方法。此方法提取的瞬態(tài)信號起始時(shí)間的準(zhǔn)確性較高，且當(dāng)背景信號或干擾噪聲較強(qiáng)時(shí)也可獲得滿意效果，但未考慮在對不同類型瞬態(tài)信號的區(qū)分。

考慮到即使在正常的齒輪系統(tǒng)中也有可能產(chǎn)生瞬態(tài)振動成分(如輕載時(shí)側(cè)隙過大)，為提取只與斷齒有關(guān)的瞬態(tài)沖擊成分，本文基于聽覺模型和信號極值點(diǎn)幅值概率密度，提出了專門針對由斷齒故障所誘發(fā)的沖擊成分提取方法。經(jīng)實(shí)驗(yàn)驗(yàn)證表明該方法能很好的提取僅與斷齒故障相關(guān)的沖擊響應(yīng)成分，并且能夠識別沖擊響應(yīng)出現(xiàn)的周期。

1 方法概述

1.1 基本原理

本文所提方法的實(shí)現(xiàn)過程如圖1所示，首先利用Gammatone濾波器組對信號進(jìn)行帶通濾波，為消除濾波時(shí)出現(xiàn)的相位不一致，對濾波之后的結(jié)果再采用逆序?yàn)V波，然后計(jì)算各濾波通道內(nèi)幅值極大值點(diǎn)的概率密度及其導(dǎo)數(shù)。根據(jù)處理結(jié)果判斷是否存在瞬態(tài)沖擊成分，繼而篩選有效幅值極值點(diǎn)并進(jìn)行區(qū)分提取，最終提取得到只含有斷齒故障的瞬態(tài)信號成分。

圖1 方法實(shí)現(xiàn)過程原理圖Fig.1 The principle diagram of the method implementation process

1.2 Gammatone濾波器

Gammatone濾波器組[11]是通過模擬人耳耳蝸基底膜的濾波特性所提出的一種帶通濾波方法。設(shè)振動信號為x(t)，共有N個(gè)濾波器，則第i個(gè)濾波器的沖擊響應(yīng)表達(dá)式為

h(i,t)=Bntn-1exp(-2πBt)cos(2πfit+φi)

(1)

式中：n為濾波器階數(shù)；φi為相位；fi為第i個(gè)濾波器的中心頻率。

參數(shù)B與fi的關(guān)系為：

B=1.019×ERB(fi)=1.019(24.7+0.108fi)

(2)

式中:ERB(fi)為Gammatone濾波器的等效矩形帶寬。

濾波器的時(shí)域波形如圖2所示。

圖2 濾波器時(shí)域波形Fig.2 The filter waveform in time domain

則帶通濾波為

y(i,t)=x(t)*h(i,t)

(3)

式中:*代表卷積。

為避免同一頻率成分在不同濾波通道中出現(xiàn)相位差異，在基于傅里葉變換進(jìn)行卷積運(yùn)算時(shí)忽略濾波器h(i,t)的相頻特性，實(shí)現(xiàn)逆序?yàn)V波，即

(4)

式中:F-1表示傅里葉逆變換，X(f)和H(i,f)分別為x(t)和h(i,t)的頻域表達(dá)。

隨后對逆序?yàn)V波之后的信號提取所有極大值點(diǎn)，其余幅值均置為零，得到y(tǒng)1(i,t)。

本文選擇200個(gè)Gammatone濾波器來實(shí)現(xiàn)帶通濾波，濾波器中心頻率通常按照對數(shù)均勻分布設(shè)置，圖3為不同中心頻率下濾波器的頻譜圖(圖中只繪出16個(gè)濾波器)。

圖3 濾波器頻譜圖Fig.3 The spectrum of filter

1.3 瞬態(tài)信號識別方法

若機(jī)械系統(tǒng)中存在碰撞與沖擊，一般情況下，其振動信號中會存在幅值明顯的瞬態(tài)沖擊響應(yīng)成分，信號的概率密度曲線也必然會出現(xiàn)如圖4所示的長拖尾現(xiàn)象[8]。

圖4 瞬態(tài)信號所致概率密度曲線的長拖尾現(xiàn)象Fig.4 Heavy-tailed in probability densitycurve induced by transient signal

但是僅依據(jù)長拖尾部分難以判斷瞬態(tài)信號的幅值范圍，因此，可提取信號的極大值點(diǎn)，并計(jì)算極大值點(diǎn)的概率密度曲線p(i,g)。由于極值點(diǎn)幅值間往往具有較強(qiáng)的不連續(xù)性，所以在極值點(diǎn)的幅值概率密度曲線中的長拖尾部分會出現(xiàn)小幅值波動情況，而且這部分的幅值必然與瞬態(tài)沖擊成分有關(guān)，基于此特點(diǎn)，為實(shí)現(xiàn)對瞬態(tài)成分的自動識別，計(jì)算p(i,g)關(guān)于g的導(dǎo)數(shù)

(5)

式中:g表示逆序?yàn)V波之后信號幅值的極大值點(diǎn)。

利用d(i,g)曲線，按以下準(zhǔn)則[10]判斷第i通道里是否含有瞬態(tài)沖擊成分：

① 在d(i,g)曲線經(jīng)歷最小值之后，存在向上過零點(diǎn)，記該點(diǎn)的橫坐標(biāo)為g0。

②ηq為一個(gè)小值，且

(6)

式中:θ是預(yù)先給定的閾值，本文取θ=0.1。

1.4 瞬態(tài)信號篩選提取

依據(jù)上述兩條判斷準(zhǔn)則，可初步提取與瞬態(tài)沖擊成分相關(guān)的極值點(diǎn)，若d(i,g)滿足以上準(zhǔn)則，則令：

(7)

式中:幅值為1的位置很可能為瞬態(tài)沖擊信號的極值點(diǎn)。

為使初步提取的極值點(diǎn)連續(xù)化，并形成波形，可將y2(i,t)和一方波信號做卷積得到y(tǒng)3(i,t)，即

y3(i,t)=y2(i,t)*r(t)

(8)

式中:r(t)表示方波信號。

由于信號在經(jīng)時(shí)頻分解之后，瞬態(tài)沖擊成分會占據(jù)一定的連續(xù)頻帶，且在某些頻率處達(dá)到峰值。為使沖擊成分辨識更加突出明顯，將y3(i,t)沿各通道方向累加求和得到θ1(t)，即

(9)

此外，為明確沖擊成分的頻域分布情況，將y3(i,t)沿著時(shí)間方向累加求和得到θ2(i)，即

(10)

繼而計(jì)算Ω(i,t)，即

Ω(i,t)=θ1(t)θ2(i)

(11)

可知，若Ω(i,t)非零，則說明在第i通道、t時(shí)刻時(shí)存在瞬態(tài)沖擊成分。

為了進(jìn)一步將y3(i,t)中無關(guān)的極值點(diǎn)刪除，使其達(dá)到簡煉的目的，繼而計(jì)算z(i,t)，即

z(i,t)=y3(i,t)Ω(i,t)

(12)

1.5 斷齒瞬態(tài)信號區(qū)分提取

由于輪齒間嚙合拍打、存在側(cè)隙或軸承出現(xiàn)小故障，均會引起瞬態(tài)沖擊成分，但是其往往高頻部分能量較小，且占據(jù)的頻帶較窄，這也意味著會出現(xiàn)與斷齒故障無關(guān)的極值點(diǎn)。

而斷齒導(dǎo)致的沖擊響應(yīng)成分往往能量較大，且會占據(jù)較寬的頻段。因此，為將淹沒在其中的斷齒沖擊響應(yīng)成分準(zhǔn)確提取，對z(i,t)進(jìn)行分頻段統(tǒng)計(jì)，平均分為α個(gè)頻段依次累加分別得到s1(t)、s2(t)、s3(t),…,sα(t)，即

(13)

(14)

(15)

?

(16)

式中:ν=N/α，N為濾波器個(gè)數(shù)。

將sk(t)(k=1,2,3,…,α)兩兩相乘得到alp(t)，即

alp(t)=sl(t)sp(t)(l,p=1,2,3,…,αl≠p)

(17)

按圖5所示的流程原則比較有值區(qū)域和無值區(qū)域，通過選取sk(t)和alp(t)來獲取提取結(jié)果。

圖5 sk(t),alp(t)選擇流程圖Fig.5 Selection flow chart of sk(t) and alp(t)

圖5中

C(t)=count{sk(t)>0,k=1,2,3,…,α}

(18)

表示當(dāng)在t時(shí)刻時(shí)，若第k個(gè)頻段上滿足sk(t)>0，則計(jì)為1，否則為零，各頻段依次將結(jié)果累加為C(t)。

繼而用alp(t)分別和含有瞬態(tài)沖擊成分的濾波信號相乘之后再累加求和得到Tq，即

(19)

最終，Tq(t)即為提取出的只含有由斷齒所引發(fā)的瞬態(tài)沖擊響應(yīng)成分。

2 實(shí)驗(yàn)驗(yàn)證

為驗(yàn)證本文所提方法在齒輪斷齒故障特征提取中的有效性，搭建如圖6所示的二級斜齒圓柱齒輪減速器實(shí)驗(yàn)臺及人為斷齒故障,其主要零部件參數(shù)為：電動機(jī)額定功率750 W,最大轉(zhuǎn)速1 400 r/min；減速器高低速級齒數(shù)分別為13，86，14，85；選用壓電式加速度傳感器拾取振動信號，考慮到?jīng)_擊成分會激發(fā)各階固有頻率，且其固有頻率往往較高，為滿足采樣定理且提高信號的還原度，設(shè)置采樣頻率為10 240 Hz。

圖6 減速器實(shí)驗(yàn)臺和人為齒輪斷齒故障Fig.6 Redactor test rig and artificial gear fracture fault

測得減速器高速軸轉(zhuǎn)速為937 r/min，理論計(jì)算旋轉(zhuǎn)頻率fr=15.62 Hz，嚙合頻率fm=203.3 Hz，圖7為斷齒故障的振動加速度信號時(shí)域波形圖和頻譜圖。在時(shí)域波形圖中雖然可以明顯看出存在瞬態(tài)沖擊成分，但難以分辨出斷齒故障特征；頻譜圖中雖然存在旋轉(zhuǎn)頻率和嚙合頻率及其高次倍頻，但譜線過于密集，不能準(zhǔn)確判斷出是何故障。

圖7 齒輪斷齒故障信號時(shí)域波形及其頻譜Fig.7 Time waveform and frequency spectrum of gear fracture fault

濾波器參數(shù)選擇如下：濾波器個(gè)數(shù)N=200；濾波器階數(shù)n=4；相位φi=0；中心頻率fi按照對數(shù)均勻分布設(shè)置。圖8是經(jīng)過逆GT之后的時(shí)頻圖，從此圖中也難以辨識瞬態(tài)沖擊成分。

圖8 逆GT濾波結(jié)果Fig.8 Filtering results of the inverse GT

圖9 各通道過零點(diǎn)Fig.9 Zero-crossing of each channel

圖10 第196通道濾波結(jié)果Fig.10 The 196th channel’s filtering results

依據(jù)給定的兩條判斷準(zhǔn)則判斷篩選后，各通道與其對應(yīng)的過零點(diǎn)g0對應(yīng)關(guān)系如圖9所示，其中幅值為零表示該通道無向上的過零點(diǎn)，即表明無瞬態(tài)沖擊成分。

以第196通道為例說明。圖10是振動信號經(jīng)逆GT濾波之后的濾波結(jié)果及局部放大圖，在圖10(b)中所指幅值為9.509。

第196通道的概率密度曲線p(g)及其導(dǎo)數(shù)d(g)如圖11和圖12所示，從11(b)中看出p(g)的長拖尾部分出現(xiàn)小幅值波動，且在第三個(gè)波峰處幅值為9.469。在12(b)中，g0=9.469是d(196,g)經(jīng)最小值之后逐漸上升且由負(fù)變正的第一個(gè)幅值點(diǎn)，并且絕大多數(shù)沖擊成分的極值點(diǎn)的幅值都會大于g0點(diǎn)，依此便可確定此通道的瞬態(tài)沖擊信號的幅值范圍。

圖11 第196通道概率密度曲線Fig.11 The 196th channel’s probability density curve

圖12 第196通道概率密度曲線的導(dǎo)數(shù)Fig.12 The 196th channel’s derivative of probability density

在10 (b)中幅值為9.509的極大值點(diǎn)與圖12 (b)的第一個(gè)過零點(diǎn)g0點(diǎn)相對應(yīng)。繼而根據(jù)上述理論依據(jù)，計(jì)算θ1(t)、θ2(t)如圖13所示，從圖13 (a)中已明顯的看出周期性瞬態(tài)沖擊成分，但依然有其他沖擊成分。

圖13 各通道累加結(jié)果Fig.13 Accumulative results of each channel

4個(gè)不同頻段的sk(t)(k=1,2,3,4)分布情況如圖14所示。由圖可知低頻段s1(t)為零，中低頻段s2(t)已出現(xiàn)明顯的周期性，且有值區(qū)域均與高頻段s4(t)相對應(yīng)。為此，經(jīng)sk(t)兩兩相乘得到alp(t)如圖15所示，因s1(t)為零，圖中只給出a23(t)、a24(t)、a34(t)的分布情況。

經(jīng)流程原則選擇篩選后，a24(t)=s2(t)s4(t)即為最終保留的篩選結(jié)果，最終提取波形如圖16(a)所示，從圖中可以看出其周期性非常明顯，每兩個(gè)沖擊之間的平均時(shí)間間隔為Δt=0.422 4 s,經(jīng)計(jì)算得中間軸轉(zhuǎn)速r=2.360 6 r/s，亦即中間軸大齒輪每轉(zhuǎn)一圈耗時(shí)為0.423 6 s，和提取的瞬態(tài)沖擊時(shí)間間隔相一致，從而表明系統(tǒng)發(fā)生斷齒故障，且發(fā)生在中間軸的大齒輪上，與人為制造的斷齒位置相匹配。

圖14 不同頻段sk(t)圖Fig.14 sk(t) diagram in different frequency bands

圖16(b)是提取結(jié)果的局部放大圖。從圖中可以看出，在經(jīng)一次嚙合過程中出現(xiàn)兩次沖擊(嚙入嚙出各一次)，對于其余三個(gè)嚙合碰撞也出現(xiàn)同樣的結(jié)果，并且在圖16(c)所示的實(shí)測信號的對應(yīng)位置波形振動趨勢和提取結(jié)果局部放大圖相一致。

圖15 alp(t)圖Fig.15 alp(t) diagram

圖16 最終信號提取結(jié)果Fig.16 Final signal extraction result

3 結(jié) 論

本文基于聽覺模型和極值點(diǎn)概率密度，提出了一種針對斷齒故障特征的提取方法，試驗(yàn)驗(yàn)證結(jié)果表明：

(1) 極值點(diǎn)概率密度經(jīng)求導(dǎo)后，利用其是否存在向上過零點(diǎn)即可自適應(yīng)的明確各濾波通道中是否存在瞬態(tài)成分及其對應(yīng)的幅值范圍。

(2) 該方法能夠?qū)帻X故障所誘發(fā)的瞬態(tài)沖擊成分與其他無關(guān)沖擊相區(qū)分，只提取與斷齒有關(guān)的沖擊響應(yīng)成分，且提取結(jié)果便于人和計(jì)算機(jī)識別，為設(shè)備斷齒故障檢測帶來便捷。

(3) 輪齒嚙合過程中，在斷齒處會產(chǎn)生兩次沖擊，可作為斷齒故障的又一特征。

[1] 唐貴基，龐爾軍，王曉龍. 基于EMD的齒輪箱齒輪故診斷的研究[J]. 機(jī)床與液壓，2013，41(13)：188-190.

TANG Guiji，PANG Erjun，WANG Xiaolong. Research on gear fault diagnosis based on EMD[J]. Mchine Tool & Hydraulics，2013，41(13)：188-190.

[2] 馮偉. 基于振動分析的齒輪斷齒故障研究[J]. 廣州航海:高等?？茖W(xué)校學(xué)報(bào)，2007，15(1)：13-16.

FENG Wei. Investigation on gear fracture failure based on vibration analysis[J]. Journal of Guang Zhou Maritime College，2007，15(1)：13-16.

[3] 馬銳，陳予恕. 齒輪傳動系統(tǒng)斷齒故障的機(jī)理研究[J]. 振動與沖擊，2013，32(21)：47-51.

MA Rui，CHEN Yushu. Fault mechanism of a gear system with tooth broken[J]. Journal of Vibration and Shock，2013，32(21)：47-51.

[4] JENA D P，SAHOO S，PANIGRAHI S N. Gear fault diagnosis using active noise cancellation and adaptive wavelet transform[J]. Measurement，2014，47:356-372.

[5] 栗茂林，梁霖，王孫安，等. 基于連續(xù)小波系數(shù)非線性流形 學(xué)習(xí)的沖擊特征提取方法[J].振動與沖擊，2012，31(1)：106-111.

LI Maolin，LIANG Lin，WANG Sunan，et al. Mechanical impact feature extraction method based on nonlinear manifold learning of continuous wavelet coefficients[J]. Journal of Vibration and Shock，2012，31(1)：106-111.

[6] WANG Shibin，CAI Gaigai，ZHU Zhongkui，et al. Transient signal analysis based on Levenberg-Marquardt method for fault feature extraction of rotating machines [J]. Mechanical Systems and Signal Processing，2015，54(55)：16-40.

[7] 嚴(yán)?？担茗P星. 一種基于形態(tài)提升的自適應(yīng)軸承微沖擊提 取方法[J]. 振動與沖擊，2013，32(24)：198-203.

YAN Baokang，ZHOU Fengxing. Adaptive weak impulse extraction method of rolling bearings based on morphological lifting wavelet[J]. Journal of Vibration and Shock，2013，32(24)：198-203.

[8] GUO W，TSE P W，Djordjevich A. Faulty bearing singnal recovery from large noise using a hybrid method based on spectral kurtosis and ensemble empirical mode decomposition [J]. Measurement，2012，45:1308-1322.

[9] KAYA E M，ELHILALI M. A temporal saliency map for modeling auditoryattention [C]//46th Annual Conference on Information Sciences and Systems. Princeton：IEEE，2012:1-6.

[10] 李允公，張金萍，戴麗. 基于極值點(diǎn)概率密度和聽覺模型的瞬態(tài)信號提取方法研究[J]. 振動與沖擊，2015，34(21)：37-53.

LI Yungong，ZHANG Jinping，DAI Li. A method extracting transient signals based on probability density of extreme points and on auditory model[J]. Journal of Vibration and Shock，2015，34(21)：37-53.

[11] LI Yungong，ZHANG Jinping，DAI Li，et al. Auditory-model-based feature extraction method for mechanical faults diagnosis[J]. Chinese Journal of Mechanical Engineering，2010，21(3)：391-397.

A method extracting fault features of gear teeth fractures based on an auditory modeland probability density of extreme points

WU Wenshou, LI Yungong, WANG Bo, LI Guomeng, SHI Yuehong

(School of Mechanical Engineering & Automation，Northeastern University，Shenyang 110819，China)

The collision and impact at gear teeth fractures are the most important features in meshing process of gear pairs. Considering a human auditory system has an instinctive response to sudden transient acoustic signals, in order to extract transient impulse response components induced by gear teeth fracture faults, a method extracting fault features of gear teeth fractures based on an auditory model and signal probability density of extreme points was proposed. Firstly, band-pass filtering with Gammatone filters , phase adjustment and extreme points extraction for signals were conducted, and then the amplitude probability densities of extreme points were calculated, their derivatives were used to judge if there are transient impact components in all filtered signals. Those extreme points with transient impact components were extracted. Meanwhile, the whole system vibration might produce extreme points being irrelevant to gear teeth fracture impacts. In order to accurately extract impacts of gear teeth fractures, according to transient signals’ frequency band continuity and multi-band distribution characteristics, an appropriate extraction method was designed. The real measured signals showed that the proposed method can accurately extract fault features of gear teeth fractures, and can extract impact components only induced by gear teeth fracture faults in a variety of transient impulse response components, and the extraction accuracy is higher.

gear teeth fracture fault; auditory model; probability density; transient signals; fault diagnosis; feature extraction

國家自然科學(xué)基金資助項(xiàng)目(51275080)

2015-12-31 修改稿收到日期：2016-03-04

吳文壽 男，碩士生，1991年生

李允公 男，博士，副教授，1976年生

TH17

A

10.13465/j.cnki.jvs.2016.19.017