基于量子疊加態(tài)參數(shù)估計的機械振動信號降噪方法

陳彥龍，張培林，王懷光

(軍械工程學(xué)院 七系，石家莊 050003)

機械狀態(tài)監(jiān)測有效方法為振動信號分析，而現(xiàn)場采集的振動信號會受噪聲干擾，直接影響機械設(shè)備狀態(tài)監(jiān)測。小波分析在振動信號降噪中得以廣泛應(yīng)用，通常假設(shè)小波分解系數(shù)滿足高斯分布，但小波壓縮特性使小波系數(shù)概率密度分布較高斯分布在零值位置更尖，并在分布兩端呈明顯拖尾趨勢，故設(shè)小波分解系數(shù)滿足高斯分布所得結(jié)論并不準(zhǔn)確。陶新民等[1]用廣義高斯分布分析軸承振動信號多尺度小波分解系數(shù)的統(tǒng)計特征，實現(xiàn)軸承故障診斷。張志剛等[2]據(jù)機械振動信號有效小波系數(shù)統(tǒng)計分布進(jìn)行拉普拉斯建模，對主減速器中齒輪故障信號進(jìn)行有效降噪。

小波變換中，雙數(shù)復(fù)小波(DTCWT)具有近似平移不變性等性質(zhì)[3]，信號分析效果更好。本文將用DTCWT實現(xiàn)一維機械振動信號的小波變換。而多數(shù)研究均假設(shè)信號雙樹復(fù)小波系數(shù)的實、虛部獨立分布，未考慮信號虛、實部間關(guān)系，若采用二維廣義高斯分布及二維拉普拉斯分布[4]對實、虛部小波系數(shù)建模，則參數(shù)估計復(fù)雜，需大量數(shù)值計算，較難實現(xiàn)。為此，本文基于二維正態(tài)分布提出含可變參數(shù)的聯(lián)合概率密度分布。并受文獻(xiàn)[5]QSP基本理論啟發(fā)，借鑒量子疊加態(tài)概念，對信號進(jìn)行分析、處理：① 對信號DTCWT小波系數(shù)建模，建立含可變參數(shù)的聯(lián)合概率密度分布；② 據(jù)貝葉斯極大后驗估計(MAP)準(zhǔn)則，結(jié)合父-子代小波系數(shù)尺度間相關(guān)性，導(dǎo)出基于疊加態(tài)參數(shù)估計的自適應(yīng)收縮函數(shù)，實現(xiàn)一維機械振動信號降噪。結(jié)果表明，該方法對機械振動信號降噪性能良好。

1 DTCWT小波系數(shù)聯(lián)合概率密度

1.1 二維概率密度函數(shù)

研究信號DTCWT小波系數(shù)分布規(guī)律，建立數(shù)學(xué)模型，有助于實現(xiàn)信號降噪，由含噪信號中估計出原信號。對含加性高斯白噪聲信號y(t)進(jìn)行小波變換，由小波變換線性性質(zhì)，有：

Y=X+N

(1)

式中：Y=Yr+iYi為含噪信號復(fù)小波系數(shù)；X=Xr+iXi為無噪信號復(fù)小波系數(shù)；N=Nr+iNi為噪聲復(fù)小波系數(shù)。

高斯白噪聲經(jīng)小波變換后仍為同方差的高斯白噪聲。用拉普拉斯函數(shù)與廣義高斯函數(shù)建立小波系數(shù)模型，但直接利用二維拉普拉斯分布及廣義高斯函數(shù)計算形式復(fù)雜，模型參數(shù)估計難以實現(xiàn)，通用性不強。為提高雙樹復(fù)小波系數(shù)分布模型通用性，本文基于二維正態(tài)分布提出帶可調(diào)參數(shù)k的概率密度函數(shù)(PDF)，設(shè)無噪信號小波系數(shù)真實概率密度分布為f0(xr,xi)，數(shù)學(xué)模型為

(2)

式中：-∞

1.2 邊緣概率密度函數(shù)

分析聯(lián)合概率密度函數(shù)可獲得邊緣密度概率密度函數(shù)，即雙樹復(fù)小波系數(shù)實、虛部分布規(guī)律，為參數(shù)估計提供依據(jù)。Xr邊緣概率密度函數(shù)為

(3)

由于

(4)

于是

(5)

令

(6)

則有

(7)

即

(8)

同理，Xi的邊緣概率密度函數(shù)為

(9)

1.3 可調(diào)參數(shù)k確定

由以上推導(dǎo)發(fā)現(xiàn)，k為小波系數(shù)概率密度函數(shù)的關(guān)鍵參數(shù)，其取值關(guān)系概率密度函數(shù)形狀及擬合精度。據(jù)各尺度中雙樹復(fù)小波系數(shù)二維概率密度函數(shù)與邊緣概率密度函數(shù)最小均方差共同確定k值：

k1=E[f0(xr,xi)-f(xr,xi)]2

(10)

k2=E[f0Xr(xr)-fXr(xr)]2

(11)

k3=E[f0Xi(xi)-fXi(xi)]2

(12)

(13)

式中：f0(xr,xi)，f0Xr(xr)，f0Xi(xi)為實際分布；f(xr,xi),fXr(xr)，fXi(xi)為曲線擬合分布。

擬合過程中，據(jù)雙樹復(fù)小波邊緣概率密度函數(shù)fXr(xr)，fXi(xi)的PDF最大值為Pr，Pi求出σr，σi：

(14)

(15)

k的誤差精度為0.1時，雙樹復(fù)小波系數(shù)分布擬合精度不夠；k的誤差精度為0.001時，分布擬合精度與k的誤差精度0.01時接近。綜合考慮分布擬合精度與計算速度，設(shè)k值的誤差精度為0.01。

2 量子疊加態(tài)參數(shù)估計模型

2.1 量子疊加態(tài)

(16)

式中：α，β為一對復(fù)數(shù)，稱量子態(tài)概率幅，且滿足

(17)

2.2 貝葉斯MAP估計器

貝葉斯理論在小波降噪中得到較好應(yīng)用。通過貝葉斯理論，可實現(xiàn)對小波系數(shù)收縮，降低噪聲。在貝葉斯統(tǒng)計理論中，最大后驗估計(MAP)為常用方法。本文據(jù)含噪信號雙樹復(fù)小波系數(shù)Y，獲得使后驗概率密度函數(shù)PX|Y取最大值的無噪小波系數(shù)X。此處以實部邊緣概率密度函數(shù)為例：

(18)

據(jù)貝葉斯原理：

(19)

設(shè)噪聲小波系數(shù)服從均值為0、標(biāo)準(zhǔn)差為σn、相關(guān)系數(shù)ρn=0的高斯分布。據(jù)式(5)及噪聲模型，式(19)可寫為

(20)

式中：f(Xr)=lnpXr(Xr)。

對式(20)求Xr的導(dǎo)數(shù)，并令其為0，得：

(21)

據(jù)Yr與Xr的邊緣概率密度函數(shù)，經(jīng)推導(dǎo)，Xr的估計式為

(22)

(23)

2.3 量子疊加態(tài)參數(shù)估計

基于量子疊加態(tài)原理展開分析，給出基于量子疊加態(tài)的參數(shù)估計算法，并對式(22)中關(guān)鍵參數(shù)進(jìn)行估計。信號、噪聲在小波域奇異特性不同，信號突變點處小波系數(shù)幅值明顯增大，而噪聲小波系數(shù)幅值將隨尺度的增加迅速減小[7]。因此可據(jù)相鄰尺度的小波系數(shù)乘積減少噪聲。

對雙樹復(fù)小波系數(shù)而言，若父系數(shù)模較大，則對應(yīng)位置的子系數(shù)同樣具有較大模值。父-子代小波系數(shù)模的乘積表達(dá)式為

(24)

(25)

(26)

(27)

(28)

當(dāng)前位置j量子疊加態(tài)信號方差σ估計式為

(29)

式中：Wj為以當(dāng)前系數(shù)Y(s,j)為中心的子窗口；Mj為窗口中系數(shù)個數(shù)。取第s尺度上鄰域窗W(j)寬度為Ms=2s+1-1。

當(dāng)前位置j量子疊加態(tài)信號期望u估計式為

u(Yj)=u(Xj+Nj)=u(Xj)+u(Nj)

(30)

由于

u(Nj)=0

(31)

故

u(Xj)=u(Yj)

(32)

本文所提基于量子疊加態(tài)的小波系數(shù)估計方法，能較好利用小波系數(shù)尺度間相關(guān)性。將式(28)、(29)代入式(22)知，對小波系數(shù)Y(s,j)，若其量子疊加態(tài)表明有用信號出現(xiàn)概率大，收縮因子自適應(yīng)變大；反之，收縮因子自適應(yīng)變小，該自適應(yīng)收縮函數(shù)有利于噪聲抑制。

3 降噪算法步驟

本文的信號降噪算法步驟為:

(1) 利用DTCWT實現(xiàn)含噪信號y(t)的小波分解；

(2) 據(jù)各尺度高頻小波系數(shù)二維概率密度函數(shù)及邊緣概率密度函數(shù)，擬合確定k值；

(4) DTCWT逆變換獲得降噪信號。

4 仿真與實際信號驗證

4.1 仿真信號

分別采用HeaviSine、Blocks、Bumps信號[9]驗證本文模型降噪的可行性及有效性。信號長度2 048點。對3種信號分別加入高斯白噪聲，使含噪信號的信噪比分別為13 dB，14 dB，15 dB，用本文方法降噪，并用5層分解及重構(gòu)。對去噪效果，采用信噪比(SNR)進(jìn)行評價。若在原始信號y(t)中加入不同噪聲n(t)，降噪后所得信號為d(t)，單位為dB，則信噪比定義為

SNR=10lg(Pd/Pn)

(33)

式中：Pd為原始信號功率；Pn為噪聲功率。信號噪聲由原始信號y(t)減去降噪信號d(t)獲得：

(34)

用傳統(tǒng)軟、硬閾值降噪法進(jìn)行效果對比，閾值λ為

(35)

表1為不同算法降噪結(jié)果。由表1看出，基于量子疊加態(tài)參數(shù)估計的信號降噪模型降噪效果最好，極大提高了信號的信噪比。表2為不同雙樹復(fù)小波基下本文降噪方法降噪效果。隨小波基變化，仍能取得良好降噪效果，表明本文方法適用性良好。

表1 不同算法降噪結(jié)果

表2 不同雙樹復(fù)小波基降噪結(jié)果

4.2 實際信號

采用某新型機械設(shè)備的綜合傳動裝置進(jìn)行研究，在軸承內(nèi)圈加工1×0.2 mm劃痕。采集加速度振動信號，實測轉(zhuǎn)動速度1 830 r/min，在3檔檔位上測量，傳感器安裝在對應(yīng)軸承位置箱蓋上方。內(nèi)圈故障理論頻率為158 Hz。軸承故障見圖1。采樣頻率12 kHz，采樣時間1 s。為便于比較，取2 048點進(jìn)行分析。故障信號波形見圖2，由頻譜中看出，故障頻率淹沒在噪聲頻率中。

圖1 機械綜合傳動裝置

不同降噪方法結(jié)果比較見圖3～圖5，運算時間見表3。

(1)波形分析。由圖形看出，軟、硬閾值降噪方法去除大量噪聲的同時也去除有用信號，使降噪后脈沖周期不準(zhǔn)確；而基于量子疊加態(tài)的參數(shù)估計降噪方法所得脈沖周期準(zhǔn)確，且故障脈沖較軟、硬閾值降噪保留更多，有益故障判斷。通過測量，降噪后信號明顯存在間隔6.3 ms沖擊，對應(yīng)于內(nèi)圈故障頻率158 Hz，可確定為滾動軸承內(nèi)圈出現(xiàn)故障。因此，量子降噪方法使噪聲大部分被消除的同時，較好保留了軸承故障沖擊特征，脈沖波形完整。

圖2 含強噪聲故障信號

(2)頻譜分析。比較降噪后信號頻譜，小波閾值降噪后，158 Hz頻率突出，但調(diào)制頻率不準(zhǔn)確，不符合軸承內(nèi)圈故障頻率特點，小波閾值降噪方法未能達(dá)到較好降噪效果；采用基于量子疊加態(tài)的參數(shù)估計降噪方法處理后，頻譜中158 Hz頻率突出，且存在明顯邊帶，轉(zhuǎn)頻(30.5 Hz)及二倍頻(61 Hz)、三倍頻(91.5 Hz)明顯，符合軸承內(nèi)圈故障頻譜特點。

圖5 軟閾值降噪消噪結(jié)果

表3 不同算法運算時間

(3) 運算效率。據(jù)表3，量子降噪時間多于軟硬閾值降噪法，但運算時間仍較短，為快速降噪方法。

5 結(jié) 論

(1) 基于二維正態(tài)分布帶可調(diào)參數(shù)k的概率密度函數(shù)考慮雙樹復(fù)小波系數(shù)實、部關(guān)系，k值計算過程兼顧二維概率密度函數(shù)與邊緣概率密度函數(shù)，使該分布模型自適應(yīng)更強。

(2) 基于量子疊加態(tài)的機械振動信號降噪方法，充分利用DTCWT的平移不變性及量子降噪模型可捕獲小波系數(shù)尺度間相關(guān)性特點，可獲得較常規(guī)小波降噪方法更高的信噪比，能有效抑制高斯白噪聲。

(3) 本文所提降噪方法可有效保留軸承沖擊信號特征，據(jù)其周期性沖擊對應(yīng)的特征頻率，可對軸承故障進(jìn)行故障識別，為強背景噪聲下機械故障信息提取提供新途徑。

[1]陶新民,徐晶,杜寶祥,等. 基于小波域廣義高斯分布的軸承故障診斷方法[J]. 機械工程學(xué)報, 2009, 45 (10): 61-67.

TAO Xin-min, XU Jing, DU Bao-xiang, et al. Bearing fault diagnosis based on wavelet-domain generalized gaussian distribution[J]. Journal of Mechanical Engineering, 2009, 45(10): 61-67.

[2]張志剛,周曉軍,宮燃,等. 小波域局部Laplace 模型降噪算法及其在機械故障診斷中應(yīng)用[J]. 機械工程學(xué)報, 2009, 45 (9): 54-57.

ZHANG Zhi-gang,ZHOU Xiao-jun,GONG Ran, et al. Denosing algorithm based on local laplace model in wavelet domain and its application in mechanical fault diagnosis[J]. Journal of Mechanical Engineering, 2009, 45 (9): 54-57.

[3]Nelson J D B, Kingsbury N G. Enhanced shift and scale tolerance for rotation invariant polar matching with dual-tree wavelets[J]. IEEE Trans on Image Processing,2011,20 (3): 814-821.

[4]Eltoft T, Kim T, Lee T W. On the multivariate laplace distribution[J]. IEEE Signal Processing Letters,2006,13(5): 300-303.

[5]Eldar Y C, Oppenheim A V. Quantum signal processing[J]. IEEE Signal Process Mag, 2002, 19 (6): 12-32.

[6]李盼池, 宋考平, 楊二龍. 基于量子門線路的量子神經(jīng)網(wǎng)絡(luò)模型及算法[J]. 控制與決策, 2012, 27 (1): 143-146.

LI Pan-chi, SONG Kao-ping, YANG Er-long. Quantum neural networks model and algorithm based on quantum gates circuit[J]. Control and Decision, 2012, 27 (1): 143-146.

[7]Zhang L, Bao P. Edge detection by scale multiplication in wavelet domain[J]. Pattern Recognition Letters,2002,23(14): 1771-1784.

[8]付曉薇, 丁明躍. 基于量子概率統(tǒng)計的醫(yī)學(xué)圖像增強算法研究[J]. 電子學(xué)報,2010,38(7):1590- 1596.

FU Xiao-wei, DING Ming-yue. Research on image enhancement algorithms of medical images based on quantum probability statistics[J]. Acta Electronica Sinica, 2010, 38 (7): 1590-1596.

[9]Donoho D L, Johnstone I M. Adapting to unknown smoothness via wavelet shrinkage[J]. Journal of American Statistical Association, 1995, 90 (432): 1200-1224.