程衛(wèi)東,趙德尊

(北京交通大學(xué) 機(jī)械與電子控制工程學(xué)院,北京 100044)

?

用于滾動(dòng)軸承轉(zhuǎn)頻估計(jì)的EMD軟閾值降噪算法

程衛(wèi)東,趙德尊

(北京交通大學(xué) 機(jī)械與電子控制工程學(xué)院,北京 100044)

摘要:針對(duì)經(jīng)驗(yàn)?zāi)Ｊ椒纸?EMD)互相關(guān)系數(shù)-峭度準(zhǔn)則降噪方法與小波閾值降噪方法的不足,提出EMD與小波軟閾值降噪相結(jié)合的降噪方法.該方法主要包括以下4部分：1)對(duì)原始信號(hào)進(jìn)行EMD分解得到固有模態(tài)函數(shù)(IMF)的集合；2)計(jì)算各個(gè)IMF與原始信號(hào)的互相關(guān)系數(shù)以及各IMF的峭度值；3)利用互相關(guān)系數(shù)-峭度準(zhǔn)則選擇需要降噪的IMF以及需要剔除的IMF；4)對(duì)選定的IMF進(jìn)行閾值降噪后與剩余的IMF相加重構(gòu)信號(hào).利用仿真和實(shí)測(cè)的故障軸承信號(hào)對(duì)所提出算法以及EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪方法進(jìn)行對(duì)比驗(yàn)證.結(jié)果表明：采用EMD軟閾值降噪方法比采用EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪方法對(duì)信號(hào)進(jìn)行預(yù)處理,更能確保軸承振動(dòng)信號(hào)的完整性,突出信號(hào)的故障特征,降低瞬時(shí)轉(zhuǎn)頻估計(jì)的誤差.

關(guān)鍵詞：經(jīng)驗(yàn)?zāi)Ｊ椒纸?EMD);軟閾值降噪；滾動(dòng)軸承；信號(hào)預(yù)處理；瞬時(shí)轉(zhuǎn)頻(IRF)估計(jì)

為了克服相關(guān)硬件設(shè)備安裝空間及安裝成本對(duì)傳統(tǒng)階比分析技術(shù)的限制[1-3],階比分析技術(shù)框架下不依靠轉(zhuǎn)速計(jì)直接根據(jù)時(shí)域振動(dòng)信號(hào)對(duì)瞬時(shí)轉(zhuǎn)頻進(jìn)行估計(jì)的研究引起了廣泛關(guān)注[4-5].Boashash等[6-7]詳細(xì)討論了瞬時(shí)頻率這一概念并從時(shí)頻分布的角度對(duì)其進(jìn)行估計(jì).郭瑜等[8]利用峰值搜索算法對(duì)電動(dòng)機(jī)變轉(zhuǎn)速振動(dòng)信號(hào)時(shí)頻圖中的瞬時(shí)頻率進(jìn)行提取以估計(jì)其瞬時(shí)轉(zhuǎn)頻.趙曉平等[9]則提出了短時(shí)傅里葉與Viterbi擬合算法相結(jié)合的方法對(duì)與臥旋離心機(jī)瞬時(shí)轉(zhuǎn)頻相關(guān)的瞬時(shí)頻率進(jìn)行提取.然而,對(duì)于變轉(zhuǎn)速滾動(dòng)軸承振動(dòng)信號(hào)，基于時(shí)頻譜的瞬時(shí)轉(zhuǎn)頻估計(jì)方法無(wú)能為力.為此，Wang等[10]提出了基于瞬時(shí)故障特征頻率(instantaneous fault characteristic frequency, IFCF)的滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻(instantaneous rotational frequency, IRF)估計(jì)方法.該方法對(duì)于滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻的提取有較好的效果.但實(shí)際工況中滾動(dòng)軸承低轉(zhuǎn)速運(yùn)行[11-12]和早期故障[13-15]普遍存在.此時(shí),故障沖擊信息易被背景噪音覆蓋,較低信號(hào)的信噪比.基于IFCF的IRF估計(jì)方法由于受到噪聲干擾,估計(jì)的瞬時(shí)轉(zhuǎn)頻將會(huì)產(chǎn)生較大的誤差,其中根據(jù)實(shí)驗(yàn)數(shù)據(jù)可知低信噪比情況下誤差可達(dá)到25%,甚至更多,這將嚴(yán)重影響后續(xù)滾動(dòng)軸承的故障診斷.因此,需要對(duì)原始信號(hào)進(jìn)行降噪預(yù)處理以降低滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻估計(jì)的誤差.

隨著經(jīng)驗(yàn)?zāi)B(tài)分解(emprical mode decomposition, EMD)方法的提出,近年來(lái)很多學(xué)者嘗試將該技術(shù)及其擴(kuò)展算法用于振動(dòng)信號(hào)的降噪預(yù)處理.蘇文勝等[16]提出了基于互相關(guān)系數(shù)-峭度準(zhǔn)則的EMD 降噪方法并運(yùn)用于滾動(dòng)軸承的早期故障診斷.彭暢等[17]則提出了EEMD(ensemble EMD)度量因子(相似性度量和基于距離的度量因子)的降噪方法.唐貴基等[18]將互相關(guān)系數(shù)-峭度準(zhǔn)則與EEMD相結(jié)合選擇最優(yōu)固有模態(tài)函數(shù)(intrinsic mode function, IMF)分量重構(gòu)信號(hào)以降低滾動(dòng)軸承振動(dòng)信號(hào)中的噪聲成分.該類算法的主要思想是根據(jù)一定準(zhǔn)則篩選出含有有用成分的IMF分量對(duì)信號(hào)進(jìn)行重構(gòu),其他分量則被認(rèn)為是偽分量或者噪聲予以剔除.然而,在低信噪比條件下,比如：早期故障或者低轉(zhuǎn)速運(yùn)行,滾動(dòng)軸承振動(dòng)信號(hào)中含有強(qiáng)烈的背景噪音,往往噪聲起主導(dǎo)作用的IMF分量中也含有故障沖擊成分.因此，單純根據(jù)相關(guān)準(zhǔn)則選擇某幾個(gè)IMF分量重構(gòu)信號(hào)容易丟失原有信號(hào)的有用信息,對(duì)后續(xù)的瞬時(shí)轉(zhuǎn)頻估計(jì)以及階比分析的結(jié)果造成影響.小波變換[19-20]由于基函數(shù)固定、多分辨率恒定等因素的限制缺乏自適應(yīng)性.對(duì)于滾動(dòng)軸承變轉(zhuǎn)速振動(dòng)信號(hào)而言,小波閾值去噪過(guò)程中最優(yōu)分解層次以及小波基函數(shù)不易確定,降噪步驟較為繁瑣.

為此,本文在上述方法的基礎(chǔ)上提出EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪與小波軟閾值降噪相結(jié)合的方法對(duì)滾動(dòng)軸承變轉(zhuǎn)速信號(hào)進(jìn)行預(yù)處理,然后根據(jù)基于瞬時(shí)故障特征頻率的滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻估計(jì)方法對(duì)轉(zhuǎn)頻進(jìn)行估計(jì).EMD軟閾值降噪方法通過(guò)對(duì)原始信號(hào)進(jìn)行EMD分解得到本征模態(tài)函數(shù)；設(shè)定互相關(guān)系數(shù)和譜峭度值閾值,根據(jù)互相關(guān)系數(shù)-峭度準(zhǔn)則剔除偽分量、選出需要進(jìn)行降噪的IMF分量,對(duì)選定的IMF分量進(jìn)行軟閾值降噪；最后將降噪后的IMF分量與剩余的IMF分量相加重構(gòu)信號(hào).

1EMD降噪算法在滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻估計(jì)中的問(wèn)題

EMD是由Huang等[21]提出的一種自適應(yīng)的信號(hào)分解方法,基于信號(hào)的局部特征,無(wú)須選擇基函數(shù)就可將原始信號(hào)分解為許多窄帶分量,即IMF.每個(gè)IMF分量滿足以下2個(gè)條件：1)在整個(gè)數(shù)據(jù)段內(nèi),極值點(diǎn)的個(gè)數(shù)和過(guò)零點(diǎn)的個(gè)數(shù)必須相等或最多相差一個(gè)；2)在任意時(shí)刻，由局部極大值點(diǎn)形成的上包絡(luò)線和由局部極小值點(diǎn)形成的下包絡(luò)線的平均值為零.IMF分量按頻率從高到低依次分布,每個(gè)IMF分量中的頻率成分隨信號(hào)本身變化而變化.不同于傳統(tǒng)的濾波器,EMD適用于非線性或者非平穩(wěn)信號(hào)[21].另外,EMD方法不受Heisenberg 測(cè)不準(zhǔn)原理[21]的限制,能夠獲得較高的分辨率.因此,EMD方法被廣泛應(yīng)用于機(jī)械振動(dòng)和地震分析等領(lǐng)域.

基于互相關(guān)系數(shù)-峭度準(zhǔn)則的EMD降噪方法中,互相關(guān)系數(shù)代表IMF分量與原始信號(hào)的相關(guān)度,其值越大,說(shuō)明IMF分量與原始信號(hào)的相關(guān)度越高.因此,可以根據(jù)互相關(guān)系數(shù)辨別EMD分解過(guò)程中由于插值誤差.端點(diǎn)效應(yīng)以及過(guò)分解等諸多因素影響而產(chǎn)生的虛假分量.互相關(guān)系數(shù)的計(jì)算公式如下：

(1)

式中：E表示求數(shù)學(xué)期望；μx和μy分別為原始信號(hào)x和y的均值；σx和σy分別為原始信號(hào)x和y的標(biāo)準(zhǔn)差.

峭度是描述波形尖峰度的一個(gè)無(wú)量綱參數(shù),其對(duì)信號(hào)的沖擊成分非常敏感,即沖擊成分比重越大,峭度值就越大.健康滾動(dòng)軸承的振動(dòng)信號(hào)峭度值分布接近正態(tài)分布,其值大約為3[16].當(dāng)滾動(dòng)軸承出現(xiàn)故障時(shí),峭度值會(huì)明顯增加,因此可以用峭度值衡量IMF分量中故障沖擊成分的多少.峭度值可由下式計(jì)算：

(2)

EMD 互相關(guān)系數(shù)-峭度準(zhǔn)則降噪方法根據(jù)互相關(guān)系數(shù)-峭度準(zhǔn)則,篩選出相關(guān)系數(shù)較大并且峭度值大于3的IMF分量對(duì)信號(hào)進(jìn)行重構(gòu),其他分量則被認(rèn)為是偽分量或者噪聲予以剔除.然而,噪聲往往表現(xiàn)為高頻信號(hào),滾動(dòng)軸承的高頻共振信號(hào)會(huì)集中分布于前幾個(gè)IMF分量.當(dāng)信號(hào)的信噪比較低時(shí),噪聲起主導(dǎo)作用的IMF分量中也存在故障沖擊成分.因此,單純的根據(jù)互相關(guān)系數(shù)-峭度準(zhǔn)則剔除峭度值小于3而互相關(guān)系數(shù)較大的IMF分量會(huì)造成重構(gòu)信號(hào)丟失原有的故障沖擊成分.

本文對(duì)故障滾動(dòng)軸承變轉(zhuǎn)速振動(dòng)信號(hào)進(jìn)行仿真(仿真信號(hào)的具體參數(shù)見(jiàn)3.1節(jié)),人為添加高斯白噪聲和信噪比rSNB=-8 dB.采用基于瞬時(shí)故障特征頻率的滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻估計(jì)方法對(duì)轉(zhuǎn)頻進(jìn)行估計(jì)(見(jiàn)圖1),圖1圓圈表示未經(jīng)過(guò)降噪獲取的IFCF趨勢(shì)線,虛線代表通過(guò)瞬時(shí)故障特征頻率趨勢(shì)線和故障特征系數(shù)(仿真信號(hào)故障特征系數(shù)C=2)計(jì)算得到的瞬時(shí)轉(zhuǎn)頻,實(shí)線表示理論瞬時(shí)轉(zhuǎn)頻.其中t代表時(shí)間,f為頻率.雖然估計(jì)的瞬時(shí)轉(zhuǎn)頻在一定程度上反映出理論轉(zhuǎn)頻的變化趨勢(shì),但是由于噪聲干擾,估計(jì)的瞬時(shí)轉(zhuǎn)頻出現(xiàn)誤差,其中最大差值為14.77 Hz,平均差值為9.36 Hz,嚴(yán)重影響了基于瞬時(shí)故障特征頻率的滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻估計(jì)方法的適用性.圖2表示經(jīng)過(guò)EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪后獲取的IFCF趨勢(shì)線、計(jì)算的瞬時(shí)轉(zhuǎn)頻以及理論瞬時(shí)轉(zhuǎn)頻,計(jì)算的瞬時(shí)轉(zhuǎn)頻(虛線)隨時(shí)間的增加反而減小,與理論的瞬時(shí)轉(zhuǎn)頻(實(shí)線)相比,變化趨勢(shì)完全相反.經(jīng)過(guò)上述討論可知：在低信噪比條件下,基于IFCF的滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻估計(jì)方法受噪聲干擾嚴(yán)重,需要對(duì)原始信號(hào)進(jìn)行降噪預(yù)處理,而EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪方法容易造成滾動(dòng)軸承振動(dòng)信號(hào)中的故障沖擊成分缺失,出現(xiàn)瞬時(shí)轉(zhuǎn)頻估計(jì)不準(zhǔn)的現(xiàn)象.

圖1　未經(jīng)過(guò)降噪獲取的瞬時(shí)轉(zhuǎn)頻和理論瞬時(shí)轉(zhuǎn)頻Fig.1　Calculated IRF without denoising and theoretical IRF

圖2　基于EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪方法得到的瞬時(shí)轉(zhuǎn)頻和理論瞬時(shí)轉(zhuǎn)頻Fig.2　Calculated IRF with EMD cross-correlation coefficient and kurtosis criterion denoising and theoretical IRF

2EMD軟閾值降噪

小波閾值降噪是由Donohon等[22-24]提出的一種降噪算法.其原理是原始振動(dòng)信號(hào)經(jīng)過(guò)小波變換之后,有用信號(hào)對(duì)應(yīng)的小波系數(shù)較大,而噪聲信號(hào)一般對(duì)應(yīng)較小的小波系數(shù),根據(jù)這一思想,通過(guò)選取合理的閾值對(duì)小波系數(shù)進(jìn)行量化處理,以達(dá)到降噪的目的.小波閾值降噪主要包括以下3個(gè)基本步驟：1)選擇小波基以及分解層次j,對(duì)信號(hào)進(jìn)行分解；2)對(duì)從第1層到第j層的每一層高頻系數(shù),選擇一個(gè)閾值進(jìn)行閾值量化處理；3)根據(jù)小波分解的第j層低頻系數(shù)和經(jīng)過(guò)量化處理后的第1層到第j層的高頻系數(shù),進(jìn)行小波重構(gòu).

小波閾值函數(shù)包括軟閾值函數(shù)和硬閾值函數(shù),其中軟閾值處理信號(hào)相對(duì)平滑,軟閾值函數(shù)[23]為

(3)

其中,N為某個(gè)尺度j的小波系數(shù)個(gè)數(shù),J是小波分解的最大層次數(shù).

然而,小波閾值去噪過(guò)程中最優(yōu)分解層數(shù)不易確定,分解以及重構(gòu)小波基函數(shù)選擇不當(dāng)同樣會(huì)造成降噪效果不理想.滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻估計(jì)方法中運(yùn)用小波閾值降噪,步驟比較繁瑣,不易得到最優(yōu)降噪效果.

文獻(xiàn)[25]提出將EMD與小波閾值降噪相結(jié)合用于齒輪故障模式識(shí)別與診斷,原理是采用小波閾值降噪處理高頻IMF分量后與低頻的IMF分量相加重構(gòu)信號(hào),最后通過(guò)時(shí)頻分析完成齒輪的故障診斷,并取得了較好的效果.然而,僅根據(jù)文獻(xiàn)[25]無(wú)法確定如何選擇需要降噪的IMF分量以及選擇幾個(gè)高頻IMF分量做降噪處理,小波閾值降噪函數(shù)也沒(méi)有具體指出.另外,由于EMD算法本身的缺陷,容易產(chǎn)生偽分量,將閾值去噪后的IMF分量與其余IMF直接相加也會(huì)忽略低頻噪聲信號(hào)和偽分量對(duì)重構(gòu)信號(hào)的影響.為此,本文將EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪與小波軟閾值降噪相結(jié)合,設(shè)定互相關(guān)系數(shù)閾值θ=0.1,剔除互相關(guān)系數(shù)ρxy<0.1的IMF分量,對(duì)互相關(guān)系數(shù)ρxy≥0.1而峭度值小于3的IMF分量進(jìn)行軟閾值降噪,將降噪后IMF分量與剩余的IMF分量相加得到重構(gòu)信號(hào).對(duì)于預(yù)設(shè)的互相關(guān)系數(shù)閾值,本文經(jīng)過(guò)多次試驗(yàn)對(duì)比確定其值為0.1較為合理,選取互相關(guān)系數(shù)大于0.1的IMF分量進(jìn)行后續(xù)處理,既能去除EMD分解過(guò)程中由于插值誤差和過(guò)分解等諸多因素而產(chǎn)生的虛假分量,也能避免含有軸承信息的IMF分量的缺失.本文選取峭度值3作為分界點(diǎn)(第一章已經(jīng)指出,健康軸承的峭度值成正態(tài)分布,其值大約為3,當(dāng)滾動(dòng)軸承出現(xiàn)故障時(shí),峭度值會(huì)明顯增加).因此峭度值大于3的IMF分量軸承故障沖擊較為明顯,無(wú)須進(jìn)一步降噪處理；而峭度值小于3的IMF分量由于噪聲信號(hào)的干擾軸承的故障特征不明顯,需要進(jìn)行降噪處理以突出故障信息.

根據(jù)式(3)給出本文所用到的軟閾值函數(shù)公式：

(4)

式中：ηj表示第j個(gè)IMF分量的函數(shù)式,γj為第j個(gè)IMF分量的閾值，計(jì)算公式如下：

EMD軟閾值降噪算法的步驟如下,具體的流程圖如圖3所示.1)對(duì)原始信號(hào)進(jìn)行EMD分解得到IMF分量；2)計(jì)算各IMF分量與原始信號(hào)的互相關(guān)系數(shù)以及各IMF分量的峭度值；3)根據(jù)互相關(guān)系數(shù)-峭度準(zhǔn)則選出需要進(jìn)行閾值降噪的IMF分量以及剔除偽IMF分量；4)對(duì)選定的IMF分量進(jìn)行閾值降噪后與剩余的IMF分量相加得到重構(gòu)信號(hào).

圖3　經(jīng)驗(yàn)?zāi)B(tài)分解(EMD)軟閾值降噪流程圖Fig.3　Flowchart of empirical mode decomposition (EMD) soft-thresholding denoising

EMD軟閾值降噪方法的優(yōu)點(diǎn)主要包括：1)具備EMD的自適應(yīng)性、多分辨率、完備性等優(yōu)勢(shì);2)確保噪聲起主導(dǎo)作用并且含故障沖擊成分的IMF分量不被剔除;3)避免了小波最優(yōu)分解層數(shù)的選擇,同時(shí)也避免了小波分解過(guò)程中小波基函數(shù)選擇的繁瑣性,從而更容易對(duì)滾動(dòng)軸承非平穩(wěn)振動(dòng)信號(hào)進(jìn)行降噪處理.

3仿真與實(shí)測(cè)信號(hào)驗(yàn)證

3.1仿真信號(hào)分析

為了驗(yàn)證本文降噪方法在低信噪比條件下滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻估計(jì)中的有效性,首先對(duì)故障軸承變轉(zhuǎn)速振動(dòng)信號(hào)進(jìn)行仿真.故障軸承非平穩(wěn)信號(hào)構(gòu)造公式[26]如下：

sin[wr(t-tm)]}.

(5)

式中：Am表示由故障引起的第m個(gè)沖擊的幅值;β為結(jié)構(gòu)衰減系數(shù);wr為由軸承故障引起的共振頻率;u(t)為單位階躍函數(shù);tm為第m個(gè)沖擊發(fā)生的時(shí)間,可由遞推下式確定：

t1=(1+μ)·1/f(t0)/n,

tm=(1+μ)·1/f(tm-1)/n,m=2,3,…,M.

(6)

式中：f(t)是軸承轉(zhuǎn)頻隨時(shí)間的變化規(guī)律；t0=0；n為每轉(zhuǎn)出現(xiàn)的故障沖擊數(shù).

本文根據(jù)式(6)對(duì)故障軸承低信噪比非平穩(wěn)信號(hào)進(jìn)行仿真,如下式：

gb=g(t)+n(t).

(7)

式中：Am與時(shí)間成線性變化,軸承瞬時(shí)轉(zhuǎn)頻隨時(shí)間的變化規(guī)律設(shè)定為f(t)=10(t)+30,系統(tǒng)的共振頻率wr=3 500Hz,衰減系數(shù)β=8,采樣頻率Fs=5 000 Hz,采樣點(diǎn)數(shù)為20 000,故障特征系數(shù)C=2；n(t)為高斯白噪聲,信噪比rSNB=-8 dB.仿真信號(hào)的時(shí)域波形圖如圖4所示,U為電壓幅值.

圖4　仿真信號(hào)的時(shí)域波形圖Fig.4　Simulated vibration signal in time domain

圖5　前6個(gè)IMF分量(IMF1-6)Fig.5　Intrinsic mode functions (IMF1-6)

根據(jù)重構(gòu)信號(hào)對(duì)滾動(dòng)軸承的瞬時(shí)轉(zhuǎn)頻進(jìn)行估計(jì),具體步驟如下：1)通過(guò)重構(gòu)信號(hào)的快速譜峭度圖,確定帶通濾波器參數(shù)進(jìn)行帶通濾波；2)對(duì)濾波信號(hào)進(jìn)行Hilbert變換以及STFT變換獲取包絡(luò)時(shí)頻圖,利用基于幅值累加的峰值搜索算法從包絡(luò)時(shí)頻圖中提取瞬時(shí)故障特征頻率并對(duì)其擬合得到瞬時(shí)故障特征頻率趨勢(shì)線；3)瞬時(shí)故障特征頻率趨勢(shì)線與故障特征系數(shù)相除即可得到瞬時(shí)轉(zhuǎn)頻.圖6(a)給出了獲取的瞬時(shí)故障特征頻率趨勢(shì)線、計(jì)算的滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻與理論瞬時(shí)轉(zhuǎn)頻的對(duì)比,理論瞬時(shí)轉(zhuǎn)頻f(t)=10(t)+30,估計(jì)的瞬時(shí)轉(zhuǎn)頻與其相吻合,其中兩者的最大差值為2.00 Hz,平均差值為1.03 Hz.分別采用本文降噪方法、EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪方法以及未采用降噪方法估計(jì)的滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻如圖6(b)所示.其中方塊表示經(jīng)EMD軟閾值降噪后估計(jì)的瞬時(shí)轉(zhuǎn)頻,虛線表示經(jīng)過(guò)EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪后獲取的瞬時(shí)轉(zhuǎn)頻,圓圈表示未降噪直接從原始信號(hào)中估計(jì)的瞬時(shí)轉(zhuǎn)頻,實(shí)線為理論瞬時(shí)轉(zhuǎn)頻.通過(guò)綜合對(duì)比可知,單純的采用EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪方法對(duì)原始信號(hào)進(jìn)行預(yù)處理容易造成信號(hào)中的故障沖擊成分的缺失,出現(xiàn)瞬時(shí)轉(zhuǎn)頻估計(jì)不準(zhǔn)現(xiàn)象；EMD軟閾值降噪優(yōu)于EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪方法,保證信號(hào)完整性的同時(shí),降低了噪聲干擾,突出了信號(hào)的故障信息,進(jìn)而瞬時(shí)轉(zhuǎn)頻的估計(jì)結(jié)果更接近理論值.

表1IMF1-IMF7與原始信號(hào)的互相關(guān)系數(shù)及其各自峭度值和IMF′的峭度值

Tab.1Cross-correlation coefficients and kurtosis values of IMF1-IMF7 and IMF′ kurtosis values

IMF序號(hào)互相關(guān)系數(shù)KIMFIMF'IMF10.72442.54313.8331IMF20.47663.3816—IMF30.30892.96164.6204IMF40.21633.0943—IMF50.15142.87054.4870IMF60.11192.80393.0255IMF70.07283.1745—

圖6　本文算法獲取的轉(zhuǎn)頻及其與傳統(tǒng)算法的綜合比對(duì)Fig.6　Calculated IRF by proposed algorithm and comparison with traditional algorithms

3.2實(shí)測(cè)信號(hào)分析

為進(jìn)一步驗(yàn)證EMD軟閾值降噪在實(shí)際工況低信噪比非平穩(wěn)條件下估計(jì)滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻方法中的實(shí)用性,利用北京交通大學(xué)機(jī)電學(xué)院現(xiàn)代制造技術(shù)綜合實(shí)驗(yàn)中心ST-5000A型多功能柔性轉(zhuǎn)子試驗(yàn)臺(tái)(見(jiàn)圖7)進(jìn)行實(shí)測(cè)驗(yàn)證.該試驗(yàn)臺(tái)由速度調(diào)節(jié)器控制電機(jī)轉(zhuǎn)速,加速度傳感器安裝在離軸承較近的位置以準(zhǔn)確的測(cè)量其振動(dòng)信號(hào),轉(zhuǎn)速計(jì)安裝在最右側(cè)用于測(cè)量轉(zhuǎn)速,采集裝置為YE6231采集卡及其配套的采集軟件.采集卡以相同采樣率同時(shí)采集兩路信號(hào),一路為滾動(dòng)軸承故障信號(hào),另一路為故障信號(hào)相對(duì)應(yīng)的轉(zhuǎn)速脈沖信號(hào).利用電火花切割凹槽模擬軸承外圈早期裂紋故障,實(shí)驗(yàn)軸承型號(hào)為6 000.該軸承的滾珠個(gè)數(shù)為7，滾動(dòng)體直徑為4.8 mm,節(jié)圓直徑為17.65 mm,接觸角為0.根據(jù)其幾何參數(shù)計(jì)算出外圈故障特征系數(shù)Co=2.548.實(shí)測(cè)信號(hào)的時(shí)域波形圖如圖8所示,采樣率為24 000 Hz,采樣時(shí)長(zhǎng)為2 s,瞬時(shí)轉(zhuǎn)頻隨時(shí)間增加而逐步增加.根據(jù)圖8可知,由于滾動(dòng)軸承的故障輕微,故障沖擊淹沒(méi)在背景噪音中,表現(xiàn)不明顯.

圖7　滾動(dòng)軸承試驗(yàn)臺(tái)布局Fig.7　Experimental setup of rolling element bearing

圖8　實(shí)測(cè)信號(hào)的時(shí)域波形圖Fig.8　Raw vibration signal in time domain

圖9　本文算法獲取的轉(zhuǎn)頻及其與傳統(tǒng)算法的綜合比對(duì)Fig.9　Calculated IRF by proposed algorithm and comparison with traditional algorithms

根據(jù)基于瞬時(shí)故障特征頻率的滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻估計(jì)方法,從未降噪的滾動(dòng)軸承實(shí)測(cè)信號(hào)中獲取的瞬時(shí)故障特征頻率趨勢(shì)線,使其與故障特征系數(shù)Co=2.548相除得到的滾動(dòng)軸承的瞬時(shí)轉(zhuǎn)頻,以及通過(guò)轉(zhuǎn)速計(jì)實(shí)測(cè)的瞬時(shí)轉(zhuǎn)頻分別如圖9(a)所示.其中,實(shí)測(cè)的瞬時(shí)轉(zhuǎn)頻首先根據(jù)轉(zhuǎn)速計(jì)測(cè)取的脈沖信號(hào)計(jì)算出離散的轉(zhuǎn)速點(diǎn),進(jìn)而通過(guò)曲線擬合所得到.實(shí)測(cè)的瞬時(shí)轉(zhuǎn)頻直接由硬件測(cè)量計(jì)算而來(lái),與振動(dòng)信號(hào)具有同步性,因此理論上可以作為驗(yàn)證本文算法的標(biāo)準(zhǔn).根據(jù)圖9(a)易知,估計(jì)的瞬時(shí)轉(zhuǎn)頻有明顯的誤差,其中最大差值約為46.10 Hz,平均差值為28.05 Hz.圖9(b)表示根據(jù)互相關(guān)系數(shù)-峭度準(zhǔn)則進(jìn)行EMD降噪后得到的瞬時(shí)故障特征頻率趨勢(shì)線、計(jì)算得到的滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻以及實(shí)測(cè)的瞬時(shí)轉(zhuǎn)頻. 可知，估計(jì)的瞬時(shí)轉(zhuǎn)頻與實(shí)測(cè)值的變化趨勢(shì)不符,出現(xiàn)了明顯失真現(xiàn)象,沒(méi)有參考價(jià)值.利用本文方法對(duì)滾動(dòng)軸承實(shí)測(cè)信號(hào)進(jìn)行降噪后獲取的瞬時(shí)故障特征頻率趨勢(shì)線、計(jì)算得到的瞬時(shí)轉(zhuǎn)頻以及實(shí)測(cè)的瞬時(shí)轉(zhuǎn)頻分別如圖9(c)所示,估計(jì)的瞬時(shí)轉(zhuǎn)頻準(zhǔn)確的反應(yīng)了實(shí)測(cè)瞬時(shí)轉(zhuǎn)頻的變化趨勢(shì),數(shù)值與實(shí)測(cè)值十分接近,最大差值約為6.30 Hz,平均差值約為2.60 Hz.圖9(d)表示實(shí)測(cè)瞬時(shí)轉(zhuǎn)頻以及分別利用3種方法估計(jì)得到的瞬時(shí)轉(zhuǎn)頻,其中,方塊表示經(jīng)EMD軟閾值降噪后估計(jì)的瞬時(shí)轉(zhuǎn)頻,虛線表示經(jīng)過(guò)EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪后獲取的瞬時(shí)轉(zhuǎn)頻,圓圈表示未降噪直接從原始信號(hào)中估計(jì)的瞬時(shí)轉(zhuǎn)頻,實(shí)線為實(shí)測(cè)瞬時(shí)轉(zhuǎn)頻.通過(guò)對(duì)滾動(dòng)軸承實(shí)測(cè)振動(dòng)信號(hào)分析進(jìn)一步說(shuō)明：EMD軟閾值降噪算法優(yōu)于EMD互相關(guān)系數(shù)-峭度準(zhǔn)則降噪算法,避免了重構(gòu)信號(hào)的失真,具有較好的降噪效果,適用于變轉(zhuǎn)速條件下故障軸承的瞬時(shí)轉(zhuǎn)頻估計(jì).

4結(jié)語(yǔ)

對(duì)仿真信號(hào)分析的結(jié)果表明,該方法能夠保證信號(hào)完整性的同時(shí),突出了故障信息,降低了噪聲干擾,進(jìn)而提高了瞬時(shí)轉(zhuǎn)頻估計(jì)的準(zhǔn)確性.本文對(duì)早期故障滾動(dòng)軸承變轉(zhuǎn)速條件下的實(shí)測(cè)信號(hào)進(jìn)行了分析,采用EMD軟閾值降噪后估計(jì)的瞬時(shí)轉(zhuǎn)頻的誤差明顯減小,進(jìn)一步驗(yàn)證了本文方法在滾動(dòng)軸承瞬時(shí)轉(zhuǎn)頻估計(jì)方法中信號(hào)預(yù)處理的有效性.因此，該方法在低信噪比條件下的滾動(dòng)軸承轉(zhuǎn)頻估計(jì)方面具有較好的應(yīng)用前景.

參考文獻(xiàn)(References):

[1] MING Z, JING L, WANG X F, et al. A tacho-less order tracking technique for large speed variations[J]. Mechanical Systems and Signal Processing, 2013, 40(1): 76-90．

[2] BONNARDOT F,BADAOUI E M, RANDALL R B, et al. Use of the acceleration signal of a gearbox in order to perform angular resampling (with limited speed fluctuation)[J]. Mechanical Systems and Signal Processing, 2005, 19 (4): 766-785．

[3] URBANEK J, BARSZCZ T, ANTONI J. A two-step procedure for estimation of instantaneous rotational speed with large fluctuations [J]. Mechanical System and Signal Processing, 2013, 38:96-102．

[4] COMBET F, ZIMROZ R. A new method for the estimation of the instantaneous speed relative fluctuation in a vibration signal based on the short time scale transform[J].Mechanical System and Signal Processing, 2009,23: 1382-1397．

[5] ZHAO D Z, LI J Y, CHENG W D. Feature extraction of faulty rolling element bearing under variable rotational speed and gear interferences conditions[J]. Shock and Vibration, 2015,2015:1-9．

[6] BOASHASH B. Interpreting and estimating the instantaneous frequency of a signal-Part1: fundamentals [J]. Proceedings of the IEEE, 1992, 80(4):520-538．

[7] BOASHASH B. Interpreting and estimating the instantaneous frequency of a signal-Part 2: algorithms and applications [J]. Proceedings of the IEEE, 1992, 80(4):540-568. [8] 郭瑜, 秦樹人, 湯寶平, 等.基于瞬時(shí)頻率估計(jì)的旋轉(zhuǎn)機(jī)械階比跟蹤[J]. 機(jī)械工程學(xué)報(bào), 2003,39(3):32-36．

GUO Yu, QIN Shu-ren, TANG Bao-ping, et al. Order tracking of rotating machinery based on instantaneous frequency estimation [J].Chinese Journal of Mechanical Engineering, 2003, 39(3):32-36．

[9] 趙曉平,趙秀莉,侯榮濤,等.一種新的旋轉(zhuǎn)機(jī)械升降速階段振動(dòng)信號(hào)的瞬時(shí)頻率估計(jì)算法[J].機(jī)械工程學(xué)報(bào),2011,47(7):103-108．

ZHAO Xiao-ping, ZHAO Xiu-li, HOU Rong-tao, et al.A new method for instantaneous frequency estimation of run-up or run-down vibration signal for rotating machinery[J]. Chinese Journal of Mechanical Engineering, 2011, 47(7):103-108．

[10] WANG T Y, LIANG M, LI J Y, et al. Rolling element bearing fault diagnosis via fault characteristic order (FCO)analysis[J].Mechanical System and Signal Processing, 2014, 45(1): 139-153.[11] DONG G M, CHEN J, ZHAO F G, A frequency-shifted bispectrum for rolling element bearing diagnosis [J]. Journal of Sound and Vibration, 2015, 339:396-418．

[12] WANG J, HE Q B, KONG F R. An improved multiscale noise tuning of stochastic resonance for identifying multiple transient faults in rolling element bearings[J]. Journal of Sound and Vibration, 2014, 333(26):7401-7421．

[13] 馮輔周,司愛(ài)威,饒國(guó)強(qiáng),等. 基于小波相關(guān)排列熵的軸承早期故障診斷技術(shù)[J]. 機(jī)械工程學(xué)報(bào),2012,48(13): 73-79．

FENG Fu-zhou, SI Ai-wei, RAO Guo-qiang, et al. Early fault diagnosis technology for bearing based on wavelet correlation permutation entropy [J]. Journal of Mechanical Engineering, 2012, 48(13): 73-79．

[14] 王仲生, 姜洪開(kāi), 徐一艷. 發(fā)動(dòng)機(jī)轉(zhuǎn)子系統(tǒng)早期故障智能診斷[J]. 航空學(xué)報(bào),2009,30(2):242-246．

WANG Zhong-sheng, JIANG Hong-kai, XU Yi-yan. Early fault intelligent diagnosis of aero-engine rotor systems [J]. Acta Aeronautica et Astronautica Sinica, 2009, 30(2):242-246．

[15] 胥永剛, 馬海龍, 付勝, 等. 機(jī)電設(shè)備早期故障微弱信號(hào)的非線性檢測(cè)方法及工程應(yīng)用[J]. 振動(dòng)工程學(xué)報(bào),2011,24(5):529-538.

XU Yong-gang, MA Hai-long, FU Sheng, et al. Theory and applications of weak signal non-linear detection method for incipient fault diagnosis of mechanical equipments [J]. Journal of Vibration Engineering, 2011, 24(5):529-538.

[16] 蘇文勝,王奉濤,張志新,等.EMD降噪和譜峭度法在滾動(dòng)軸承早期故障診斷中的應(yīng)用[J].振動(dòng)與沖擊,2010,29(3): 18-21．SU Wen-sheng, WANG Feng-tao, ZHANG Zhi-xin, et al. Application of EMD denoising and spectral kurtosis in early fault diagnosis of rolling element bearing [J]. Journal of Vibration and Shock, 2010, 29(3):18-21．

[17] 彭暢,柏林,謝小亮. 基于EEMD、度量因子和快速譜峭度圖的滾動(dòng)軸承故障診斷方法[J].振動(dòng)與沖擊2012, 31(20):144-146．

PENG Chang, BO Lin, XIE Xiao-liang. Fault diagnosis method of rolling element bearing based on EEMD, measure-factor and fast kurtogram[J].Journal of Vibration and Shock, 2012, 31(20): 144-146. [18] 唐貴基,王曉龍.基于EEMD降噪和1.5維能量譜的滾動(dòng)軸承故障診斷研究[J].振動(dòng)與沖擊, 2014,33(1):6-10．

TANG Gui-ji, WANG Xiao-long. Fault diagnosis for roller bearings based on EEMD de-noising and 1.5-dimensional energy spectrum [J]. Journal of Vibration and Shock, 2014, 33(1):6-10．

[19] 程正興,楊守志,馮曉霞.小波分析的理論、算法、進(jìn)展與應(yīng)用[M].北京:國(guó)防工業(yè)出版社,2007:332-335．

[20] BADOUR A F, SUNAR M, CHEDED L. Vibration analysis of rotating machinery using time-frequency analysis and wavelet techniques [J].Mechanical Systems and Signal Processing,2011,25(6): 2083-2101．

[21] HUANG N E, SHEN Z, LONG S R, et al. The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis[J]. Proceedings of the Royal Society of London,1998, 454(1971):903-995．

[22] DONOHO D L, JOHNSTONE I M. Threshold selection for wavelet shrinkage of noisy data [J].Engineering in Medicine and Biology Society, 1994, 1:A24-A25．

[23] DONOHO D L. Denoising by soft-thresholding[J]. IEEE Transaction on Information Theory, 1995, 41(3):613- 627．

[24] DONOHO D L, JOHNSTONE I. Wavelet shrinkage asymptotic [J].Journal of Royal Statistical Society, 1995, 57(2):301-369．

[25] 邵忍平,曹精明,李永龍.基于EMD小波閾值去噪和時(shí)頻分析的齒輪故障模式識(shí)別與診斷[J].振動(dòng)與沖擊,2012,31(8):96-101．

SHAO Ren-ping, CAO Jing-ming, LI Yong-long. Gear fault pattern identification and diagnosis using time- frequency analysis and wavelet threshold de-noising based on EMD [J].Journal of Vibration and Shock, 2012, 31(8):96-101．

[26] 王天楊, 李建勇, 程衛(wèi)東.基于改進(jìn)的自適應(yīng)噪聲消除和故障特征階比譜的齒輪噪源干擾下變轉(zhuǎn)速滾動(dòng)軸承故障診斷[J].振動(dòng)與沖擊, 2014,33(18):7-13．

WANG Tian-yang, LI Jian-yong, CHENG Wei-dong. Fault diagnosis of rolling bearing under a variable rotational speed and gear vibration noise based on revised ANC algorithm and FCO spectrum [J]. Journal of Vibration and Shock, 2014, 33 (18):7-13．

DOI:10.3785/j.issn.1008-973X.2016.03.005

收稿日期：2015-04-10.

基金項(xiàng)目：國(guó)家自然科學(xué)基金資助項(xiàng)目(51275030).

作者簡(jiǎn)介：程衛(wèi)東(1967-), 男，副教授，從事制造裝備智能測(cè)控與故障診斷教學(xué)科研工作.ORCID:0000-0001-5085-4758. E-mail: wdcheng@bjtu.edu.cn

中圖分類號(hào)：TH 113.1

文獻(xiàn)標(biāo)志碼：A

文章編號(hào)：1008-973X(2016)03-08-0428

EMD soft-thresholding denoising algorithm for rolling element bearing rotational frequency estimation

CHENG Wei-dong, ZHAO De-zun

(SchoolofMechanicalElectronicandControlEngineering,BeijngJiaotongUniversity,Beijing100044,China)

Abstract:A method based on the empirical mode decomposition (EMD) and wavelet shrinkage was proposed due to the shortcomings of the EMD cross-correlation coefficient and kurtosis criterion denoising and wavelet shrinkage. The method consists of four main steps: (i) IMFs were obtained by decomposing raw signal, (ii) The cross-correlation coefficient between IMFs and the raw signal, and the kurtosis values of the IMFs were calculated, (iii) IMFs with the noise were selected and the false IMFs were removed, (vi) the noise of the selected IMFs was removed by the soft-thresholding denoising method, and then the signal with the rest of the IMFs was reconstructed. The proposed method was tested based on both the simulated and experimental bearing vibration signals. Results show that, compared with EMD cross-correlation coefficient and kurtosis criterion denoising, the method of EMD soft-thresholding denoising can ensure the integrity of the signal, highlight fault features and reduce the error of the instantaneous rotational frequency (IRF).

Key words:empirical mode decomposition (EMD); soft-thresholding denoising; rolling element bearing; signal preprocessing; instantaneous rotational frequency (IRF) estimation