陳旭陽(yáng),韓振南,王志堅(jiān)
(1.太原理工大學(xué) 機(jī)械與運(yùn)載工程學(xué)院,太原 030024;2.中北大學(xué) 機(jī)械工程學(xué)院,太原 030051)
由于機(jī)械設(shè)備運(yùn)轉(zhuǎn)時(shí)采集到的振動(dòng)信號(hào)普遍存在噪聲干擾,因此降噪和故障特征提取成為機(jī)械故障診斷的重點(diǎn)和難點(diǎn)。準(zhǔn)確提取故障特征,對(duì)有效地降低機(jī)械設(shè)備振動(dòng)信號(hào)噪聲具有非常重要的意義。
雙樹(shù)復(fù)小波變換(dual-tree complex wavelet decomposition,DTCWT)是傳統(tǒng)離散小波變換的改進(jìn)方法,是由KINGSBURY正式提出的[1],即用奇、偶濾波器組實(shí)現(xiàn)小波的分解與重構(gòu),不僅具有傳統(tǒng)小波的時(shí)頻局部化分析能力,而且有離散小波所不具備的抗頻率混疊性、近似平移不變性、完全重構(gòu)性、有限的數(shù)據(jù)冗余和高效的計(jì)算效率等優(yōu)良性質(zhì)[2-3]。
基于參數(shù)估計(jì)和小波系數(shù)之間的相關(guān)性提出的閾值選取方法,普遍缺點(diǎn)是難以選取合適的閾值,從而影響了實(shí)際應(yīng)用[4-5]。
高階統(tǒng)計(jì)量方法對(duì)多種噪聲都有很好的抑制作用。它不僅對(duì)未知自相關(guān)加性噪聲不敏感,抑制高斯色噪聲,而且對(duì)另一類均勻?qū)ΨQ分布的非高斯有色噪聲也不敏感,因此成為非平穩(wěn)、非高斯信號(hào)和非最小相位系統(tǒng)的主要數(shù)學(xué)分析工具[6-7]。
由于小波分解層數(shù)的多少直接影響信號(hào)的去噪效果,因此采用粒子群算法自適應(yīng)選擇小波的分解層數(shù)。本文提出了自適應(yīng)雙樹(shù)復(fù)小波和高階累積量的信號(hào)降噪方法,利用小波系數(shù)的周期性和小波系數(shù)模震蕩小的特點(diǎn)以及高斯隨機(jī)過(guò)程的高階累積量為零的特性實(shí)現(xiàn)信號(hào)降噪處理。
雙樹(shù)復(fù)小波變換采用二叉樹(shù)結(jié)構(gòu)的兩路濾波器組進(jìn)行信號(hào)的分解和重構(gòu)(見(jiàn)圖1),一樹(shù)生成實(shí)部,一樹(shù)生成虛部,合理設(shè)計(jì)實(shí)、虛部樹(shù)低通濾波器,滿足半采樣延遲條件,具有近似平移不變性。兩樹(shù)濾波器采樣頻率相同,但是它們之間的延遲恰好是一個(gè)采樣間隔,這樣虛部樹(shù)中第一層的二抽取恰好采到實(shí)部樹(shù)中二抽取所丟掉的采樣值。
3層雙樹(shù)復(fù)小波分解和重構(gòu)過(guò)程,如圖1所示。分解時(shí),h0、h1為實(shí)部樹(shù)低、高通濾波器,g0、g1為虛部樹(shù)低、高通濾波器。重構(gòu)時(shí),h'0、h'1為實(shí)部樹(shù)濾波器組,g'0、g'1為虛部樹(shù)濾波器組。
小波函數(shù)表示為如下形式:
ψ(t)=ψh(t)+iψg(t) .
(1)
式中:ψh(t),ψg(t)為兩個(gè)實(shí)小波;i為復(fù)數(shù)單位。
圖1中,虛線上方實(shí)部樹(shù)變換的小波系數(shù)和尺度系數(shù)可由式(2)、式(3)計(jì)算:
(2)
(3)
同理,下方虛部樹(shù)變換的小波系數(shù)和尺度系數(shù)可由式(4)和式(5)計(jì)算:
(4)
(5)
因此,可得到雙樹(shù)復(fù)小波變換的小波系數(shù)和尺度系數(shù):
(6)
(7)
最后,雙樹(shù)復(fù)小波變換的小波系數(shù)和尺度系數(shù)可由式(8)和(9)重構(gòu):
(8)
(9)
雙樹(shù)復(fù)小波變換后的重構(gòu)信號(hào)可表示為:
(10)
近幾年的研究結(jié)果表明,大多數(shù)機(jī)械故障信號(hào)都是非高斯、非平穩(wěn)、非線性過(guò)程,僅用二階統(tǒng)計(jì)量很難全面描述其特性。高階累積量對(duì)高斯噪聲與對(duì)稱分布噪聲中的信號(hào)檢測(cè)、恢復(fù)和特征增強(qiáng)等方面與二階統(tǒng)計(jì)量相比有獨(dú)特之處,特別是非平穩(wěn)、非線性過(guò)程的處理。二階累積量處理的信號(hào),去噪不徹底。
由于旋轉(zhuǎn)機(jī)械周期性的工作性質(zhì),導(dǎo)致了在工作過(guò)程中測(cè)試的信號(hào)呈現(xiàn)周期性,因而采集到的信號(hào)大多數(shù)為對(duì)稱的。對(duì)于對(duì)稱分布的有用的非噪聲信號(hào),其三階累積量也會(huì)為0.但不管是對(duì)稱分布還是非對(duì)稱分布的非高斯信號(hào),其四階累積量都不為0[8].
一般地,k維隨機(jī)變量(x1,x2,…,xk),其聯(lián)合概率密度函數(shù)表示為f(x1,x2,…,xk),則第一聯(lián)合特征函數(shù)被定義為:
(11)
對(duì)應(yīng)的第二聯(lián)合特征函數(shù)為:
ψ(ω1,ω2,…,ωk)=lnφ(ω1,ω2,…,ωk)=lnE{ej(ω1xx+ω2x2+…ωkxk)} .
(12)
因此,k個(gè)隨機(jī)變量(x1,x2,…,xk)的r=r1+r2+…+rk階聯(lián)合矩為對(duì)第一聯(lián)合特征函數(shù)求k次導(dǎo)數(shù)得:
(13)
同樣地,k個(gè)隨機(jī)變量的r階累積量為對(duì)第二聯(lián)合特征函數(shù)求k次導(dǎo)數(shù)得:
(14)
一組隨機(jī)變量的k階累積量可以表示為:
Ck=cum(x1,x2,…,xk) .
(15)
式中:cum()代表聯(lián)合累積量。設(shè){x(n)}為零均值的實(shí)隨機(jī)過(guò)程,k階累積量Ckx(τ1,τ2,…τk-1)定義為:
Ckx(τ1,τ2,…,τk-1)=cum{(x(n),x(n+τ1),…,x(n+τk-1))} .
(16)
實(shí)際中,可以由M-C(矩-累積量)公式直接得到一個(gè)最簡(jiǎn)單的關(guān)于高階累積量的關(guān)系式。
二階累積量:
C2x(τ)=E{x(n)x(n+τ)}=Rx(τ) .
(17)
三階累積量:
C3x(τ1,τ2)=E{x(n)x(n+τ1)·x(n+τ2)}=m3x(τ1,τ2) .
(18)
四階累積量:
C4x(τ1,τ2,τ3)=E{x(n)x(n+τ1)x(n+τ2)·x(n+τ3)}-Rx(τ1)Rx(τ2-τ3)-Rx(τ2)·Rx(τ3-τ1)-…-Rx(τ3)Rx(τ1-τ2) .
(19)
式中,Rx(τ)為{x(n)}的二階矩,即自相關(guān)函數(shù)。如果x為高斯隨機(jī)過(guò)程,則x的四階累積量為0.
不同于進(jìn)化類算法的復(fù)雜遺傳因子的過(guò)程,粒子群優(yōu)化算法采用速度-位移模式進(jìn)行全局搜索優(yōu)化。同時(shí)繼承了人工智能的記憶特性,可以根據(jù)當(dāng)前最好的那個(gè)粒子位置及時(shí)調(diào)整搜索策略,是一種高效的搜索算法[9]。
粒子群算法的思想來(lái)源于鳥(niǎo)群覓食的群體共享機(jī)制模擬。每個(gè)粒子都是解空間中的一個(gè)解,根據(jù)自己和同伴的經(jīng)驗(yàn)隨時(shí)改變自己的速度和位置(即層數(shù))從而到達(dá)目標(biāo)地點(diǎn)。每個(gè)粒子的最優(yōu)解即個(gè)體極值為pb而整個(gè)群體的最優(yōu)解即全局極值為gb,每個(gè)粒子都是通過(guò)pb和gb隨時(shí)調(diào)整自己的位置和速度。改變位置和速度的公式為:
Vi(t+1)=ωVi(t)+c1rand(pbi-xi(t))+c2rand(gb-xi(t)) .
(20)
Xi(t+1)=xi(t)+Vi(t) .
(21)
具體步驟為:
1) 初始化各粒子的速度和位置。
2) 確定一個(gè)適應(yīng)度函數(shù),將降噪后的信噪比和降噪前信噪比差值的絕對(duì)值作為適應(yīng)度函數(shù)。比較粒子個(gè)體與種群的適應(yīng)度值,個(gè)體中信噪比值最大的粒子為個(gè)體最優(yōu)值,種群中信噪比值最大的粒子為種群最優(yōu)值。
3) 通過(guò)式(20)和式(21)更新粒子的速度和位置。i表示第i個(gè)粒子,t表示迭代次數(shù),rand為[0,1]上的隨機(jī)數(shù),c1,c2為學(xué)習(xí)因子。種群規(guī)模t取10次,視問(wèn)題的復(fù)雜程度而定;ω取0.4~1;和取1.0.
4) 計(jì)算粒子和種群的適應(yīng)度函數(shù)值,通過(guò)與前一步中的個(gè)體和種群最優(yōu)值進(jìn)行比較,更新個(gè)體和種群最優(yōu)值。判斷新種群是否達(dá)到終止條件,達(dá)到則執(zhí)行步驟5),否則繼續(xù)執(zhí)行步驟4).
5) 退出迭代,輸出全局最優(yōu)值(即層數(shù))。
假設(shè)信號(hào)中含有未知功率密度譜的零均值高斯噪聲,信號(hào)模型可以表示為:
x(k)=f(k)+n(k) .
(22)
假設(shè)xj是尺度j上的小波系數(shù),因?yàn)樾〔ㄊ蔷€性變換,而且信號(hào)相互獨(dú)立,所以
xj=fj+nj.
(23)
式中,n是高斯噪聲,f是信號(hào)成分。利用高斯隨機(jī)過(guò)程的四階累積量為0,則:
(24)
可以得到:
(25)
含信號(hào)成分的小波系數(shù)的四階累積量大于0,而噪聲的4類累積量等于0,因此去掉了小波系數(shù)的噪聲,保留小波系數(shù)的有用信號(hào)。估計(jì)噪聲方差值[7]:
(26)
在尺度j上,噪聲的方差是:
(27)
其中,G是高、低通濾波器增益。
對(duì)各層小波系數(shù)進(jìn)行處理時(shí)采用如下原則:
(28)
1) 利用粒子群算法,選擇雙樹(shù)復(fù)小波分解最優(yōu)分解層數(shù)。
2) 對(duì)含噪信號(hào)按粒子群所尋到的最優(yōu)分解層數(shù)進(jìn)行雙樹(shù)復(fù)小波分解。
3) 將各層雙樹(shù)復(fù)小波系數(shù)進(jìn)行反變換,分別重構(gòu)出各層系數(shù)。
4) 對(duì)重構(gòu)后的高頻小波系數(shù)進(jìn)行式(27)的門(mén)限降噪處理。
5) 將降噪后的信號(hào)與傳統(tǒng)的閾值處理雙樹(shù)復(fù)小波系數(shù)降噪得到的信號(hào)進(jìn)行對(duì)比,得出結(jié)論。具體流程如圖2所示。
圖2 實(shí)施流程圖Fig.2 Implementation flow chart
構(gòu)造如下齒輪局部故障仿真信號(hào)驗(yàn)證本文方法的可行性和有效性,調(diào)幅-調(diào)頻信號(hào)的表達(dá)式:
x(t)=0.2[1+cos(2π×30t)]+[1+cos(2π×30t)]×cos(2π×120t)+[1+cos(2π×30t)]×cos[2π×150t+cos(2π×5t)] .
(29)
由圖3可以看出用粒子群尋優(yōu)時(shí),當(dāng)分解層數(shù)為3層時(shí)所得適應(yīng)度值最大。
圖3 適應(yīng)度值隨分解層數(shù)的變化Fig.3 Fitness values varies with the number of decomposed layers
圖4為純凈信號(hào)的波形圖和頻譜圖。在以上純凈信號(hào)中加入隨機(jī)白噪聲,如圖5所示為含噪信號(hào)波形圖和頻譜圖。30 Hz頻率被淹沒(méi)在噪聲中,無(wú)法識(shí)別出。
采用本文方法進(jìn)行降噪處理,降噪效果如下圖6所示,可以清晰地識(shí)別出各特征頻率成分,降噪效果比較理想。
圖7所示為用硬閾值方法降噪后的信號(hào)的波形及頻譜圖。從圖中可以看出30 Hz的特征頻率被淹沒(méi)在噪聲中,降噪效果不理想。
圖4 純凈信號(hào)波形圖及頻譜圖Fig.4 Pure signal waveform diagram and spectrogram
圖5 含噪信號(hào)波形圖及頻譜圖Fig.5 Noise signal waveform diagram and spectrogram
圖6 本文降噪信號(hào)波形圖及頻譜圖Fig.6 This article noise reduction signal waveform diagram and spectrogram
圖8所示為用軟閾值方法降噪后的信號(hào)的波形及頻譜圖。從圖中可以看出30 Hz的特征頻率被淹沒(méi)在噪聲中,降噪效果不理想。
通過(guò)對(duì)比從波形圖中可以明顯、清晰地看出經(jīng)本文方法處理后的信號(hào)降噪效果比較理想且從頻譜圖可以識(shí)別出各特征頻率成分,和軟、硬閾值法相比更能有效地降低信號(hào)噪聲,提高信噪比。
圖7 硬閾值降噪信號(hào)波形圖及頻譜圖Fig.7 Hard threshold noise reduction signal waveform diagram and spectrogram
圖8 軟閾值降噪信號(hào)波形圖及頻譜圖Fig.8 Soft threshold noise reduction signal waveform diagram and spectrogram
實(shí)驗(yàn)信號(hào)是通過(guò)齒輪疲勞強(qiáng)度實(shí)驗(yàn)臺(tái)通過(guò)傳感器采集到的振動(dòng)信號(hào)。如圖9所示,加速度傳感器安裝在齒輪箱端蓋上,電動(dòng)機(jī)轉(zhuǎn)速保持在1 200 r/min左右,扭矩為800 N左右。采樣點(diǎn)數(shù)為1 024,采樣頻率為1 000 Hz.試驗(yàn)齒輪的傳動(dòng)比為1∶1,采取半齒嚙合。
圖9 齒輪疲勞強(qiáng)度實(shí)驗(yàn)臺(tái)Fig.9 Gear fatigue strength test bench
在齒輪嚙合的實(shí)測(cè)振動(dòng)信號(hào)中,由于在齒輪嚙合過(guò)程中故障信號(hào)(齒輪點(diǎn)蝕)被淹沒(méi)在電動(dòng)機(jī)轉(zhuǎn)動(dòng)和散熱大風(fēng)扇轉(zhuǎn)動(dòng)等強(qiáng)背景噪聲下,因而難以從采集到的信號(hào)中直接判斷出齒輪箱是否存在故障。采用本文的降噪方法處理含噪信號(hào),使之能夠有效降噪,從而能更好地提取故障特征。圖10為無(wú)故障的實(shí)驗(yàn)信號(hào)波形圖及其頻譜和包絡(luò)譜圖,從圖中可以看出嚙合頻率的半頻55 Hz和一倍頻110 Hz處雖然比較明顯,但有干擾頻率存在。
圖10 無(wú)故障實(shí)驗(yàn)信號(hào)及其譜圖Fig.10 No trouble experimental signal and its spectrum
圖11為降噪后的無(wú)故障信號(hào)的波形圖、頻譜圖和包絡(luò)譜圖,從包絡(luò)譜圖中可以看出半頻和一倍頻非常清晰、明顯且干擾頻率成分明顯減少。
圖11 降噪后的信號(hào)及其譜圖Fig.11 Signal after noise reduction and its spectrum
圖12是經(jīng)過(guò)軟閾值處理后的信號(hào)波形圖和幅值譜、包絡(luò)譜圖。從包絡(luò)譜和幅值譜圖中可以看到在100~500 Hz之間干擾頻率較多。
圖12 軟閾值降噪后的信號(hào)及其譜圖Fig.12 Signal after soft-threshold noise reduction and its spectrum
通過(guò)圖10—圖12的對(duì)比,可以清晰地看到本文方法的降噪效果要優(yōu)于傳統(tǒng)的軟閾值降噪方法,從兩種方法的包絡(luò)譜圖中可以看出,本文方法比軟閾值法更加有效的消除噪聲信號(hào),減少了頻率干擾。
因?yàn)闊o(wú)法知道降噪前信號(hào)的信噪比,因此將降噪后信號(hào)的信噪比絕對(duì)值的大小作為適應(yīng)度函數(shù)。選取降噪后信號(hào)信噪比絕對(duì)值的最大值所對(duì)應(yīng)的小波分解層數(shù)為最優(yōu)分解層數(shù)。由圖13可以看出在第二層之后,降噪后信號(hào)的信噪比基本趨于平穩(wěn),因此本文選擇的小波最優(yōu)分解層數(shù)為3層。
圖13 信噪比的大小隨分解層數(shù)的變化Fig.13 Size of the signal to noise ratio varies with the number of decomposed layers
圖14所示為加速度傳感器安裝在齒輪箱端蓋上,電動(dòng)機(jī)轉(zhuǎn)速保持在1 200 r/min左右,扭矩為1 000 N左右。采樣點(diǎn)數(shù)為1 024,采樣頻率為1 000 Hz,采集到的是有故障的實(shí)驗(yàn)信號(hào)。因轉(zhuǎn)速控制有偏差因此計(jì)算得到的嚙合頻率存在偏差。
通過(guò)傳感器和信號(hào)采集儀采集到的實(shí)驗(yàn)信號(hào)波形圖以及包絡(luò)譜圖如圖14所示。由于故障信號(hào)比較微弱,而且淹沒(méi)在強(qiáng)背景噪聲下,從幅值譜圖和包絡(luò)譜圖中看出有邊頻成分和很?chē)?yán)重的頻率干擾20 Hz以及其他干擾頻率成分。從信號(hào)的包絡(luò)譜圖中不能準(zhǔn)確地反映出故障特征頻率。
圖14 故障信號(hào)波形及其頻譜圖Fig.14 Fault signal waveform and its spectrum
為了能更好地降低噪聲干擾,診斷出故障信號(hào),故使用本文方法進(jìn)行降噪,降噪后的信號(hào)波形圖和包絡(luò)譜圖如圖15所示。從圖中可以很清晰地看出半頻50.55 Hz處和嚙合頻率107.4 Hz與正常齒輪(圖12)對(duì)比發(fā)現(xiàn),在半頻和嚙合頻率周?chē)咁l成分明顯,由此可看出存在故障特征。與圖14對(duì)比看出,波形清晰,干擾信號(hào)幅值有很大程度降低,干擾頻率成分被消除,去噪效果明顯;且與無(wú)故障圖11相比邊頻成分明顯,幅值不均一。
圖15 故障信號(hào)降噪后的波形及其頻譜圖Fig.15 Waveform of the fault signal after noise reduction and its spectrum
1) 將四階累積量方法引入雙樹(shù)復(fù)小波變換降噪中,根據(jù)信號(hào)和噪聲的統(tǒng)計(jì)特性進(jìn)行信噪分離。取得了很好的降噪效果。
2) 因?yàn)殡p樹(shù)復(fù)小波分解不同分解層數(shù)會(huì)影響降噪效果,分解層數(shù)少則降噪效果不理想,而分解層數(shù)太大則容易丟失有用信號(hào)。即使分解層數(shù)大,降噪效果很好,但丟失有用信號(hào)后,信噪比反而會(huì)降低。利用信噪比作為適應(yīng)度函數(shù),通過(guò)粒子群優(yōu)化選擇信噪比最大的分解層數(shù)。
3) 仿真信號(hào)降噪結(jié)果表明,和雙樹(shù)復(fù)小波變換的傳統(tǒng)軟、硬閾值法相比,該方法在不同信號(hào)和噪聲水平下均能表現(xiàn)出良好的自適應(yīng)性和降噪效果。
4) 實(shí)驗(yàn)信號(hào)處理結(jié)果表明,該方法能夠有效抑制齒輪箱振動(dòng)信號(hào)中的強(qiáng)背景噪聲,信號(hào)處理中能更好地減少干擾頻率,有效提取故障頻率。