機械振動信號降噪技術(shù)分析研究設(shè)計
機械振動信號降噪技術(shù)分析研究設(shè)計,機械振動,信號,技術(shù),分析研究,設(shè)計
工程壓縮機信號S3在小波和濾波器降噪中為了不與其他數(shù)據(jù)沖突改為SS3;
dt=0.0004。
機械振動信號降噪技術(shù)分析研究
摘 要
旋轉(zhuǎn)的機械在其運行的過程中振動中包含了很豐富的運行狀態(tài)信息,所以對振動信號進行診斷監(jiān)測就顯得十分必要。然而在實際的獲取信號過程中,由于測試環(huán)境的不同、測試儀器的不同以及不可預測的人為因素等諸多干擾的情況下,被測信號往往不僅僅包括有用振動信號還參雜著大量的隨機的干擾信號,也就是噪聲。若不對被測信號進行一定的降噪處理,之后無論是提取的特征還是對數(shù)據(jù)的處理都會有不可忽略甚至影響研究方向的誤差,這在工程應(yīng)用中顯然是不允許的,所以對振動信號進行有效的去降噪處理是一個迫切需要完成的任務(wù)。機械振動信號的降噪技術(shù)分析顯然能為我們解決此類問題所帶來的困擾,因為它能凈化采集的的振動信號,即增大信噪比,進而使我們獲得更準確的振動特征參數(shù)等,以便能了解利用或消除振動。為消除噪聲從而獲取較為純凈的振動信號所提出了很多方法,有時頻分析法、傅里葉變換、最優(yōu)估計法、最優(yōu)濾波法、自適應(yīng)濾波、小波分析和經(jīng)驗?zāi)B(tài)分解等方法。本文首先通過對于振動信號降噪的重要性進行說明后,提出了為什么會選用濾波器分析、小波分析和經(jīng)驗?zāi)B(tài)分解三種方法的原因并簡要說明了其各自的發(fā)展,進而通過對濾波器濾波信號降噪技術(shù)、小濾變換降噪信號技術(shù)和經(jīng)驗?zāi)B(tài)分解信號降噪技術(shù)原理的分析研究,根據(jù)其各自降噪機理和特點利用Matlab 軟件對簡單數(shù)字信號的進行模擬降噪過程來闡述以上三種方法的實用性,最后對一個通過壓縮機信號的降噪處理來驗證振動信號降噪在實際工程信號處理中的應(yīng)用。
關(guān)鍵詞:振動信號降噪;信噪比;小波降噪;濾波器降噪;經(jīng)驗?zāi)B(tài)分解降噪
- I -
ABSTRACT
In the course of its operation, the vibration of the rotating machine contains a wealth of operational information, so it is necessary to diagnose and monitor the vibration signal . However, in the actual process of signal acquisition, due to different testing environment and testing instruments of different and unpredictable human factors such as interference, the measured signal is often not only includes the useful vibration signal also mixed with random interference signal is a lot of noise. If the measured signal noise reduction processing, after both feature extraction and processing of data will have negligible error and even affect the research direction, which is obviously not allowed in engineering application, so the vibration signal effectively to noise reduction processing is an urgent need to complete the task. Mechanical vibration signal analysis and noise reduction technology can obviously bring us to solve such problems, because it can purify the vibration signal acquisition, which increase the SNR is large enough, then we get the vibration parameters more accurately, so as to understand the use of or eliminate vibration. In order to eliminate the noise and vibration signal acquisition is pure proposed a lot of methods, sometimes frequency analysis, Fu Liye transform, optimal estimation, optimal filtering, adaptive filtering, wavelet analysis and empirical mode decomposition method. This will be through the filter, filter noise reduction technology transform noise reduction technology and EMD denoising technology principle and development and utilization of the Matlab software to simulate the simple digital signal denoising process to illustrate the above three methods practicality, at the end of a compressor by the noise reduction signal to verify the vibration and noise reduction the signal in the signal processing in practical engineering application.
Key words:Vibration signal denoising; signal to noise ratio; wavelet denoising; filter noise reduction; empirical mode decomposition
- II -
目 錄
摘 要 I
ABSTRACT II
1 緒論 1
1.1 論文研究的意義 1
1.2 論文國內(nèi)外目前的研究成果 2
1.3 論文研究的發(fā)展方向 3
2 濾波器的信號降噪實現(xiàn) 5
2.1 噪聲分析 5
2.2 濾波器降噪 6
2.3 FIR 濾波器的信號降噪模擬 8
2.4 IIR 濾波器的信號降噪模擬 10
2.5 濾波器的工程信號降噪實例 13
3 小波降噪 16
3.1 振動信號降噪的小波分析原理 16
3.2 振動信號降噪的小波包分析原理 19
3.3 小波分析的降噪模擬 20
3.4 閾值選擇對小波分析影響的降噪模擬 22
3.5 小波包分析的降噪模擬 24
3.6 小波及小波包分析的工程信號降噪實例 26
4 經(jīng)驗?zāi)B(tài)分解降噪 28
4.1 經(jīng)驗?zāi)B(tài)分解原理 28
4.2 EMD 方法的信號噪降模擬 31
4.3 經(jīng)驗?zāi)B(tài)分析的工程信號降噪實例 33
結(jié) 論 35
參 考 文 獻 36
附錄 A 外文資料譯文 38
附錄 B 外文資料原文 46
致 謝 49
- III -
1 緒論
1.1 論文研究的意義
機械在運行過程中產(chǎn)生會產(chǎn)生振動,提取這些包含豐富運行狀態(tài)的振動信號有助于我們對振動機械進行監(jiān)測和診斷[1]。然而在實際的獲取信號過程中,由于被測的振動信號中不可避免的會受到大量干擾信號的的影響。若不對被測的受干擾信號進行除噪處理,分離出有用的振動信號,其特征提取和數(shù)據(jù)分析的可信度就不高,故而影響后續(xù)的研究,這在工程應(yīng)用中顯然是不允許的,所以對振動信號進行有效的去降噪處理是一個迫切需要完成的任務(wù)。目前已經(jīng)被提出的降噪方法有時頻分析法、傅里葉變換、最優(yōu)估計法、最優(yōu)濾波法、自適應(yīng)濾波、小波分析和經(jīng)驗?zāi)B(tài)分解等方法[2]。
濾波器的信號降噪分析,能有效的控制噪聲,對所選用的頻率范圍的幅值及伴隨的特性進行平衡,因此濾波器在信號處理中發(fā)揮著至關(guān)重要的作用。滿足研究條件的振動信號在濾波器中均可以視作無數(shù)的正弦波疊加而成并且正弦波越多,就能越精確的表示出被測信號。這些頻率不同的振動信號組成成分波一般統(tǒng)稱之為諧波成分。由于濾波器是根據(jù)信號的頻率的不同來對信號進行過濾,所以可以通過濾波器的這一特性對已知頻率不同的振動信號和噪聲信號進行分離。
小波分析,利用其時頻分析能力強,可同時采用變化的時頻窗來對振動信號進行多分辨分析,小波分析還能觀察到信號的概貌,即能看到信號細節(jié)[3],于是越來越廣泛的應(yīng)用在了非平穩(wěn)信號的非平衡信號的特征提取及信號處理中。小波分析在非平穩(wěn)信號處理領(lǐng)域內(nèi)的作用巨大,不同尺度下它的處理過程不同,使信號的整體特性和局部特性都能得到較好的反應(yīng)。小波分析的精度可以調(diào)節(jié),這一性質(zhì)既能使其定位到被測信號中的高頻成分也能使其定位到被測信號中的低頻成分,它能做到在時域上的任何位置均有分辨率,這也使小波變換提取的信號信息更加詳細,小波分析所用依據(jù)的小波基緊支性會集中信號的能量于少數(shù)小波系數(shù),達到集中有用信號的目的,噪聲在經(jīng)過小波分析后一般都極小,所以小波分析僅需要大致了解信號的類型去套用標準的降噪方法即可,這一簡單實用的優(yōu)點促使小波分析在信號降噪中得到廣泛應(yīng)用[4]。
- 37 -
經(jīng)驗?zāi)B(tài)分解的全稱是 Empirical Mode Decomposition,簡稱 EMD。是近二十年來新興的信號處理方法。EMD 方法的自適應(yīng)能力強,無需對信號進行預先分析情況下就可以實現(xiàn)信號成分的分離,所以其對機械振動信號的降噪處理的靈活性相較于其他降噪方法有明顯的優(yōu)勢。
1.2 論文國內(nèi)外目前的研究成果
濾波器的研究的重點是:在最優(yōu)情況下如何制造出最合適的濾波器。1940s,現(xiàn)代濾波理論的創(chuàng)使人維納創(chuàng)造了維納濾波器。一,輸入的振動信號是廣義上是平穩(wěn)的;二,輸入的信號必須已知其統(tǒng)計特性是維納濾波器必須滿足的兩個條件。這使維納濾波器不適用于不平穩(wěn)和統(tǒng)計特性未知的信號中[6]。60 年代初卡爾曼創(chuàng)立了卡爾曼濾波理論讓濾波器可以用于非平穩(wěn)信號處理,卡爾曼通過采用狀態(tài)量模型估計非平穩(wěn)和多輸入輸出的隨機信號。但設(shè)計卡爾曼濾波器需預先知道信號和噪聲的統(tǒng)計特性。由于維納濾波和卡爾曼濾波均要已知輸入信號與噪聲統(tǒng)計特性的弊端限制了濾波器在信號處理中應(yīng)用。自適濾波器的出現(xiàn)解決了這一問題,它自 1967 被提出后就隨大規(guī)模集成計算機和電路技術(shù)的快速發(fā)展而發(fā)展,自適應(yīng)降噪就是其應(yīng)用分支。同時自適用濾波理論也與優(yōu)化理論、估計理論與檢測和信息論等密切相關(guān)。
傳統(tǒng)的信號分析方法是傅里葉分析,由于它只能將振動信號在頻域內(nèi)展開就使其不能了解到振動信號隨時間的狀態(tài)信息。這雖然在一些情況下是很有必要的,但由于時域信息對一些信號的重要性,人們將傅里葉分析方法的應(yīng)用范圍擴大,并提出大量時頻分析方法。其中有:短時傅里葉變換、Rand-Wigner、Gabor 等,也包括本次論文重點要重點使用的方法之一:小波變換。
LIN 人等根據(jù)對含有大量噪聲的齒輪箱振動信號的研究處理過程中,提出了一種新的以 morlet 小波基為基礎(chǔ)的連續(xù)小波變換的降噪方法,在時間域內(nèi)連續(xù)小波變換實現(xiàn)機械振動信號降噪方法,實現(xiàn)了周期性沖擊信號的提取。陳志新等用復小波塊閥值降噪法提取弱故障特征信息的方法,取得了不錯的效果。N.G.Nikolaou 用小波變換提取故障特征來對軸承故障進行診斷,效果顯著。何正嘉等在回顧對振動機械的監(jiān)測和診斷領(lǐng)域內(nèi)關(guān)于小波技術(shù)的大量研究成果,總結(jié)出了小波降噪技術(shù)技術(shù)在往復機械、齒輪和軸承等
機械工程方面應(yīng)用中所取得的進展總結(jié)出許多有價值的通用性的特點和規(guī)律為一些未知領(lǐng)域提供了先驗知識。
自 1998 年被發(fā)明到現(xiàn)在短短不到 20 年的時間里,EMD 迅速的發(fā)展成為了一個被廣泛應(yīng)用的信號處理方法。國際上,Huang 等在建立了 HHT 的基本框架上,根據(jù) HHT 的基本原理,發(fā)明了固有模式函數(shù)概念,他在提出 EMD 方法的同時定義了希爾伯特譜及邊際譜的概念,分析了關(guān)于 EMD 的完備性條件與正交性的相關(guān)問題。大量比較和分析 EMD 與小波等的信號分析的差別,Yang 和 Salvino 將此方法用于設(shè)備的診斷領(lǐng)域[7-9];Rai 等用 EMD 提取軸承振動信號的 IMF,并且結(jié)合 FFT 的技術(shù)提取了 IMF 的頻譜特征并且取得了不錯的效果[10];Fandrin[11-12]與Wu zhaohua[13-14]分別用分數(shù)高斯噪聲和白噪聲驗證了 EMD 實質(zhì)上是與小波分解類似并具有品質(zhì)因數(shù)的恒定的帶通濾波器, 這為實現(xiàn)濾波器的 EMD 化設(shè)計打下了基礎(chǔ)。國內(nèi),鐘佑明從理論角度對 HHT 局部乘積定理的認證,并在 HHT 提供了的統(tǒng)一理論依據(jù)的基礎(chǔ)上,提出了邊際譜和傅立葉譜之間存在本質(zhì)區(qū)別這一新的觀點[15]。沈國際提出了 EMD 分解多頻信號的一個必要條件,將其用于齒輪箱故障振動信號分離,取得了良好效果[16]。不過作為一種新興信號分析方法EMD 方法的強大優(yōu)勢雖然已經(jīng)在很多領(lǐng)域體取得了成果,但該方法屬于“經(jīng)驗性”算
法,所以無論在理論基礎(chǔ)或者是實際應(yīng)用中的諸多環(huán)節(jié)依然存在著不足,例如虛假模式、端點效應(yīng)等問題,這極大影響了信號特征的精確性。
1.3 論文研究的發(fā)展方向
近年來提出了不同尺度的元素級聯(lián)合而成的的開閉與組合的廣義形態(tài)濾波器,這種濾波方式能通過對傳統(tǒng)濾波器由于采用相同尺寸結(jié)構(gòu)元素易造成輸出統(tǒng)計變形的缺點而改進來實現(xiàn)振動信號降噪的目的[17]。僅由加減、極大值和極小植運算構(gòu)成的廣義形態(tài)濾波器使降噪運算即簡單以高效,應(yīng)用價值較高。
將小波分析方法和其它如:小波與相關(guān)分析、小波與包絡(luò)譜、小波與獨立分量分析、小波與奇異值分解等相結(jié)合的特征提取方法,小波結(jié)合支持向量機和神經(jīng)網(wǎng)絡(luò)的診斷方法等都得到了研究和應(yīng)用[18]。
經(jīng)驗?zāi)B(tài)分解的研究主要體現(xiàn)在關(guān)于對端點效應(yīng)的處理、模態(tài)混疊、停止準則以及均值曲線擬合這四個方面 。端點處理是對信號篩選過程中的端點的處理,使端點處的發(fā)
散現(xiàn)象減輕;模態(tài)混疊是時間尺度的分解與特征模態(tài)函數(shù)的不對應(yīng),這會破壞各個被分解信號的所包含的物理意義;停止準則就是確定經(jīng)驗?zāi)B(tài)分解篩選過程的是否停止或繼續(xù)的判斷依據(jù),以使信號在滿足指定要求的情況下減少算法的運算量和復雜性;均值曲線擬合是由于上下包絡(luò)線在擬合的過程中的不準確所帶來的“過沖”和“欠沖”進行造成特征模態(tài)的多余或者丟失。EMD 方法也有了和其它信號處理方法結(jié)合的發(fā)展趨勢。比如 EMD 與小波、EMD 與 SVM、EMD 與 ICA 以及 EMD 與神經(jīng)網(wǎng)絡(luò)等方法進行綜合特征提取與診斷方法。
2 濾波器的信號降噪實現(xiàn)
2.1 噪聲分析
2.1.1 機械振動信號噪聲分析
振動的采集是指測量在選定點上進行位移、速度和加速度的加速度的大小。但是由于機械振動信號采集過程中由于測試儀器的精度,測試環(huán)境和其它一些人為的因素都會使被測振動信號與實際振動信號存在一定的偏差,這些干擾被測振動信號的偏差就稱為“噪聲”。
機械振動信號通常有許多噪聲源,且一般隨機性比較強,出現(xiàn)的頻域范圍也比較廣泛而不穩(wěn)定,通常噪聲信號的頻率比機械振動信號的頻率高。噪聲信號的波形一般波動大而不線性,因而難以用經(jīng)驗來預測。
機械振動信號噪聲的危害為機械振動信號的作用所決定,由于機械振動信號的噪聲再利用在現(xiàn)階段較難解決。振動信號的測量可以實現(xiàn)對機械故障進行分析和運行設(shè)備的, 但是夾雜著噪聲的振動信號會使振動信號的可信賴程度降低,甚至可能完全使原始的振動信號被淹沒,這顯然會影響振動信號的后續(xù)分析,所以信號的降噪是伴隨振動信號測量的一個不可或缺的部分。
3.1.2 機械振動噪聲分析
機械振動噪聲一般分為以下三種:
(1) 機械性噪聲:由于固體振動而產(chǎn)生的如齒輪傳動部件、曲柄連桿部件、鏈傳動部件、軸承部件、液壓系統(tǒng)部件等多種運動部件產(chǎn)生的噪聲;
(2) 氣體動力性噪聲:由于氣體振動而產(chǎn)生,如風機、活塞式發(fā)動機、空壓機(115dB)等所產(chǎn)生的噪聲,以及因氣流激發(fā)而產(chǎn)生的噪聲,如噴注噪聲、排氣噪聲、卡門旋渦噪聲;
(3) 電磁性噪聲:因高頻諧磁場的相互作用,引起電磁性振動而產(chǎn)的噪聲機械振動降噪的危害主要體現(xiàn)在以下三個方面:
(1) 對人休心理的危害:噪聲會擾亂人們的心情,影響人們的正常通話交流,分散人的注意力而有發(fā)生其它危險的可能等;
(2) 對人體生理的影響:強度過大和時間過長的噪聲會干擾人的正常睡眠,使人疲勞,影響人的聽力,嚴重的還會導致失聰?shù)龋?
(3) 對生產(chǎn)活動是影響:噪聲可能會影響建筑物,比如強度很大的噪聲會使玻璃振碎,振動噪聲會干擾信號的傳遞,不利于對器械的運行狀態(tài)時行判斷,噪聲疲勞還能導致飛機失事這樣的嚴重事故。
2.2 濾波器降噪
2.2.1 濾波器的原理
利用信號特性的不同來對被研究信號進行篩選,從中獲取信噪比高的有用信號的步驟叫做濾波,以濾波為目的而被使用的的工具統(tǒng)稱為濾波器。經(jīng)典的濾波器可以被分為低通濾波器、高通濾波器、帶通濾波器和帶阻濾波器來實現(xiàn)不同的使用要求,每一種濾波器又根據(jù)所處理信號在時間域內(nèi)是否連續(xù)而被分為數(shù)字濾波器和模擬濾波器。模擬濾波器的優(yōu)點在于處理的速度快, 帶寬大同時不需要 ADC(analog to Digital Converter 模擬數(shù)字轉(zhuǎn)換器) 和 DAC (Digital to Analog Converter 數(shù)字模擬轉(zhuǎn)換器),但是由于它可重復性不強及抗干擾能力弱的原因,模擬濾波器在信號處理中的應(yīng)用受到了限制。而數(shù)字濾波器是利用有限精度算法來實現(xiàn)的、將輸入的數(shù)字信號轉(zhuǎn)換為分析所需要的輸出數(shù)字信號的 LTI 系統(tǒng),其可以通過硬件和軟件兩個部分來實現(xiàn)。從硬件的角度上看,數(shù)字濾波器可以被表示為數(shù)字加法器、數(shù)字乘法器和數(shù)字延時單元經(jīng)過一定運算后構(gòu)成的信號過濾設(shè)備;而從軟件的角度看,數(shù)字濾波器則是一個算法或者計算程序。所以雖然數(shù)字濾波器在硬平臺(如 DSP 芯片)上運行,但濾波器的性能不由硬件決定,它只提供程序運行的基礎(chǔ),濾波器的性能是濾波器的數(shù)字系數(shù)和其算法決定的。
數(shù)字濾波器的頻率響應(yīng)在理想條件下–也就是其幅值只有在通帶頻域中能完全保
留,而在阻帶頻率中被完全過濾:
2.1
2.2
2.2.2 FIR 濾波器
FIR 濾波器是限定時間內(nèi)在原始信號上加入固定值的脈沖信號而反饋輸出信號的一種濾波方式。所以其始終是穩(wěn)定的并能夠用硬件實現(xiàn),由于 FIR 濾波器非線性相位的特性使其可以用于沒有明顯相位信息的信號上,在處理普通的數(shù)字信號時,使用有限長的脈沖響應(yīng)單元能使計算的時間延遲控制在較小的范圍內(nèi),這是瞬態(tài)信號處理的關(guān)鍵。窗函數(shù)法是濾波器進行信號濾波最常用的方法。
公式 2.1 與公式 2.2 對應(yīng)的理想單位脈沖響應(yīng)hd (n) 為:
1
hd (n) =
wc e- jn0we jnw dw = sin(wc (n - n0 ))
ò
-w
2.3 2.3
2p c
p (n - n0 )
從公式中可以看出低通濾波器的單位脈沖響應(yīng) 由無窮的非因果關(guān)系序列構(gòu)成。當把一個具有物理可實現(xiàn)意義的 N 層因果關(guān)系的相位線性相關(guān)的濾波器當作設(shè)計目的時, 要使被測脈沖響應(yīng)hd (n) 被截取成長度是 N 的一段,將起點的位置取為n = 0 處,其中心位置即為n0 = (N -1) / 2 ,進而確保所截取序列對(N-1)/2 是對稱的。設(shè)截取的偶對稱脈沖響應(yīng)表示為 h(n),將它作為實數(shù)濾波器的系數(shù)向量。截取即相當于乘以一個矩形窗:
h(n) = hd (n)RN (n)
2.4
上式中的 RN (n)是寬度不限而長度為 N 的矩形矩陣,其中心位置與hd (n) 對準。實際
N -1
濾波器系統(tǒng) H(z)是 H (z) = ?h(n)z -n
n=0
用有限長的序列 h(n)去代替hd (n) ,這顯然會帶
,
來誤差,在頻率上的表現(xiàn)就是吉布斯效應(yīng)。這個效應(yīng)會使通帶和阻帶內(nèi)有波動性,它的產(chǎn)生原因是由于hd (n) 將頻率帶強制分離,因此被叫做截斷效應(yīng)。脈沖響應(yīng)的截斷使實際頻率特性和理想的特性有誤差,這同時也是濾波器設(shè)計中主要所關(guān)心和要解決的問題, 通過定量分析后選擇適當?shù)拇昂瘮?shù)可以盡可能的減少這種截斷效應(yīng),從而得到特定技術(shù)要求的 FIR 線性相位濾波器。式 2.5 可以表示對振動信號進行傅里葉變換后再根據(jù)卷積定理后過程:
H (e jw )= 1 p H (e jq ) A
(e j (w -q ) )dq
2p ò-p d RN
2.5
在進行信號降噪處理過程中所用函數(shù)fir(1 )的核心思想是采用窗函數(shù)來實現(xiàn)所需功能的 FIR 濾波器的構(gòu)建。它通過選擇合適的窗函數(shù)w(n) ,單位沖激響應(yīng)表示為h(n) ,進而求解濾波器的系數(shù)b(n) = w(n)h(n) ,其中1 £ n £ N 。常用的窗函數(shù)包括:
(1) 矩形窗(Rectangular window);
(2) 漢寧窗(Hanning window);
(3) 海明窗(Hamming window);
(4) 布萊克曼窗(Blackman window)。
2.2.3 IIR 濾波器
IIR 濾波器是無限長時間內(nèi)在原始信號上持續(xù)的加入固定值的脈沖信號而反饋輸出信號的一種濾波方式,其功能實現(xiàn)的遞推關(guān)系式如式 2.6:
N M
y[n] = -?ak y[n - k] + ?bx x[n - k]
2.6
k =1 k =0
M
?bk z
N
-k
IIR 函數(shù)濾波器的系統(tǒng)函數(shù)為:
H (z) = k =0
2.7 2.7
1 - ? ak k =1
z -k
可以采用相應(yīng)函數(shù)分別來設(shè)計巴特沃茲濾波器、切比雪夫 I 型濾波器、切比雪夫Ⅱ 型濾濾器、橢圓數(shù)字濾波器 和遞歸型的數(shù)字濾波器。
功率譜估計只能通過一部分已測得的數(shù)據(jù)來對信號進行處理。它能表達隨機信號在頻域上的統(tǒng)計特性,并且具有明顯的實際物理含義,在振動信號的降噪處理過程中有重要的意義,由于實際得到的信號長度一般都有限,用這種有限長度信號所得功率譜僅代表了隨機信號真實功率譜的估計,所以叫做功率譜估計。經(jīng)典的功率譜估計方法包括周期圖法、Nuttall 法、Welch 法、Bartlett 法和最大熵譜估計法等,這些方法均是以傅立葉變換為基礎(chǔ)的。
2.3 FIR 濾波器的信號降噪模擬
以公式 3.8 為原始信號來對信號進行降噪模擬處理,原始信號主要有頻率分別 5Hz、15Hz 和 40Hz 的正弦信號以及標準方差為 1 的白噪聲信號。
x = 4 ′sin(2 ′p ′ 5′ t) + sin(2′p ′15′ t) + sin(2′p ′ 40′ t) + randn(size(t))
2.3.1 降噪過程
(1) 選擇濾波器的階數(shù)為 20 的帶通濾波器,通頻帶是 0.01-20Hz;
(2) 選擇 fir1 函數(shù)設(shè)計濾波器;
2.8
(3) filterz() 函數(shù)對信號信號降噪濾波,并將其和原始信號圖示之;
(4) 對原始信號和降噪信號分別進行傅里葉變換并圖示頻率譜。
2.3.2 噪降結(jié)果
如圖 2.1 通過對比可以明顯的看出,降噪后的信號明顯則顯的更加平緩,沒有了幅值最高的尖點和一些突兀且沒有規(guī)律的波形;
圖 2.1 FIR 濾波器降噪波形圖
從圖 2.2 中原始信號和降噪后信號兩個頻率譜圖的對比中可以看出,降噪后的信號有效的去除了信號中高頻的部分和隨機的部分。
圖 2.2 FIR 濾波器降噪頻譜圖
2.3.3 噪降代碼
clc; clear all; close all; fs=100;
t=0:1/fs:1-1/fs;
x=4*sin(2*pi*5*t)+5*sin(2*pi*15*t)+5*sin(2*pi*40*t)+0.5*randn(size(t));
n=20;
Wn=[0.01 20]*2/fs;
b=fir1(n,Wn);
y=filter(b,1,x); figure; subplot(2,1,1);
plot(t,x);
title('原始信號');
xlabel('時間t/s');
subplot(2,1,2);
plot(t,y);
title('濾波后的信號');
xlabel(' 時 間 t/s'); set(gcf,'position',[200,200,800,300]); figure(2);
subplot(2,1,1); fftx=fft(x); halffftx=abs(fftx(1:50)); df=0:1:50-1;
plot(df,halffftx); subplot(2,1,2); ffty=fft(y); halfffty=abs(ffty(1:50)); df=0:1:50-1;
plot(df,halfffty); set(gcf,'position',[200,200,800,300])
2.4 IIR 濾波器的信號降噪模擬
為了便于對比不同濾波器的降噪效果,選取與 FIR 濾波器實現(xiàn)所用相同的原則信號進行模擬降噪信號,信號由式 2.8 給出。
2.4.1 模擬降噪過程
(1) 確定濾波器的通帶頻率和阻帶頻率;
(2) 通過 cheb1od 函數(shù)來確定 1 型切比雪夫濾波器的階數(shù);
(3) 函數(shù) cheby1 函數(shù)設(shè)計濾波器;
(4) 使用函數(shù) filter 對 cheby1 設(shè)計好的濾波器進行降噪并將其結(jié)果與原始信號圖示之;
(5) 分別對原始信號和降噪后的信號進行傅里葉變換,用頻率譜圖示之;
2.4.2 降噪結(jié)果
圖 2.3 IIF 濾波器(Ⅰ型切比雪夫濾波器)降噪波形圖
圖 2.3 所示與 FIR 濾波器的效果類似,I 型切比雪夫濾波器也過濾了大量無規(guī)律的白噪聲及高頻信號,這從圖中每個周期內(nèi)僅有兩個波峰和兩個波谷,而且沒有了最鋒利的波形可以看出。
圖 2.4 IIF 濾波器(Ⅰ型切比雪夫濾波器)降噪頻譜圖
從圖 2.4 中降噪前后頻率譜圖能看出明顯的降噪效果,降噪前的信號在 5Hz、15Hz、40Hz 都比較集中,且從 0 到 50Hz 的采樣頻率范圍內(nèi)也一直存在,這顯然會影響對真實的振動信號的處理和分析,經(jīng)過降噪后的信號明顯便加集中在了 5Hz 及 15Hz 的部分。
2.4.3 降噪代碼
clc; clear all; fs=100;
t=0:1/fs:1-1/fs; x=4*sin(2*pi*5*t)+5*sin(2*pi*15*t)+5*sin(2*pi*40*t)+0.5*randn(size(t)); Wp=20*2/fs;
Ws=25*2/fs; Rp=1; Rs=50;
[n,Wn]=cheb1ord(Wp,Ws,Rp,Rs);
[b,a]=cheby1(n,Rp,Wn); y=filter(b,a,x);
figure; subplot(211); plot(t,x);
title('原始信號'); subplot(212); plot(t,y);
title(' 濾 波 后 的 信 號 '); set(gcf,'position',[200,200,800,300]); figure;
subplot(211); fftx=fft(x); halffftx=abs(fftx(1:50)); df=0:1:50-1;
plot(df,halffftx); subplot(212);
ffty=fft(y); halfffty=abs(ffty(1:50)); df=0:1:50-1;
plot(df,halfffty); set(gcf,'position',[200,200,800,300])
對同一信號采用不同的濾波器的降噪實例中可以看出不同濾波器在MATLAB中的的具體實現(xiàn)過程和降噪效果。實現(xiàn)相同性能指標的信號降噪濾波器時,由于IIR濾波器的線性相位不相關(guān),但FIR濾波器的階數(shù)比IIR濾波器的階數(shù)高5到10倍,這會導致它的成本高、信號時延大若要求線性相位,IIR濾波器則要加全通網(wǎng)絡(luò)來校正相位,這也會很大的提高其階數(shù)。而此時的IIR濾波器則可以得到嚴格的線性相位。
另外,IIR模擬濾波器這一本質(zhì)是其無法避開的不足;而FIR濾波器的使用則靈活很多而且還可以用于一些IIR濾波器中切比雪夫等逼近不能設(shè)計出指定標準的情況。因此在實際應(yīng)用過程中要考慮的因素很多。若在相位要求不高的情況下IIR濾波運算量小、速度快且成本低的優(yōu)點就得以發(fā)揮;而對線性相位要求較高的場合,F(xiàn)IR濾波器更加合適。機械振動信號的處理中,因其對信號的的相位要求不高,所以IIR是較為合理的選擇。
2.5 濾波器的工程信號降噪實例
三種降噪方法均采用同一個壓縮機泄漏故障的信號進行降噪處理,以下的工程信號降噪均默認此條件。
2.5.1 FIR濾波器的降噪應(yīng)用實例
圖 2.5 和圖 2.6 是原始的壓縮機信號和經(jīng)過 FIR 濾波器降噪后的壓縮機信號的時域?qū)Ρ葓D,本此降噪選擇對信號大于 2000Hz 頻率的部分進行降噪,圖 2.5 可以明顯看出其波形的幅值的變化,一些鋒利的尖點明顯的減少或消失了。圖 2.6 以 0~5000Hz 為例來分析降噪效果能看出能看出信號的降噪的主要在 2000 到 3000Hz,尤其是 2500 到3000Hz 信號幾乎全部被消除,它基本實現(xiàn)了 FIR 濾波器的降噪目的。
圖 2.5 FIR 濾波器的工程信號降噪時域圖
圖 2.6 FIR 濾波器的工程信號降噪頻域圖
2.5.2 IIR 濾波器的降噪應(yīng)用實例
同樣對高于 2000Hz 的噪聲信號進行去噪處理時,圖 2.7 與圖 2.8 的波形與頻譜圖和FIR 濾波器的降噪時頻圖類似。頻譜圖以 0~5000Hz 為例進行降噪效果分析,IIR 濾波器的降噪也主要體現(xiàn)在 2000 到 3000Hz 上,這符合預期降噪結(jié)果。IIR 濾波器降噪和 FIR 濾波器降噪的差異在本例中差別不是很大,也驗證了 FIR 濾波器和 IIR 濾波器在對相位信息不敏感的信號的處理可以達到類似效果的結(jié)論,但是 FIR 濾波器的層數(shù)更高,程序的運行次數(shù)也更大。
圖 2.7 IIR 濾波器的工程信號降噪時域圖
圖 2.8 IIR 濾波器的工程信號降噪時域圖
3 小波降噪
3.1 振動信號降噪的小波分析原理
頻域譜分析是一種以確定各頻率段信號的強弱為目的的定量分析方法;時域分析也廣泛應(yīng)用于工程技術(shù)。但大多數(shù)時候,僅僅使用時域或頻域分析在振動信號處理過程中都是不嚴謹?shù)?于是一種即能進行時域分析又能進行頻域分析的方法就迫切的需要被提出,而小波分析剛好可以實現(xiàn)這樣的要求。
FFT 是傅里葉變換中最先被使用來加入時域信息的時頻分析方法,它的前提是需要假設(shè)在很短的時間窗內(nèi)的信號是平穩(wěn)且連續(xù)的,通過分隔時間窗的方法進行頻域分析的過程進而實現(xiàn)對局部時間的頻域信息的獲取。但由于時間窗的大小不可改變,這對瞬態(tài)信號的分析來說還顯的力不從心。相較于短時傅立葉變換不能改變時間窗大小的缺點, 小波變換的時間窗可以根據(jù)處理信號的不同甚至信號所包含成分的不同進行改變,當信號不穩(wěn)定時,就選用較小的時間窗來達到精確分解信號的目的,當信號較穩(wěn)定時,選用較大的時間窗來減小運算量從來減小程序的復雜性和冗余性,由于具備這樣的優(yōu)勢使得小波換在信號的降噪處理中取得了不小的成果。
由于小波分析在時域窗和頻域窗所具有靈活性使小波分析在信號降噪的應(yīng)用中也有了很多分支,不過都有三個必要的步驟[19]:
(1) 把信號分解到小波域;
振動信號首先應(yīng)該進行信號分解,可用于分解信號的小波變換方法有連續(xù)小波變換、離散小波變換、多分辨率分析及 Mallat 算法等。
連續(xù)小波變換(CWT)定義是信號乘以小波函數(shù)的尺度與位移在所有時間之和,表示為式 3.1 [20]:
C(scale, position) =
+¥
ò f (t)Y(scale, position)dt
-¥
3.1
經(jīng)過 CWT 的輸出量的是關(guān)于時移和尺度函數(shù)的小波系數(shù),這些小波系數(shù)可以作用在合適的位移和尺度上,從而得到振動信號的小波。
離散小波變換:為了避免在所有尺度上都要求進行大量的小濾系數(shù)計算,只考慮部分尺度和時移來計算是很有意義的。在此背景下,二進尺度和位移就可以使信號分析既十分有效又相當精確,即 DWT。對于多數(shù)的振動信號,低頻成分常常蘊含信號的特征, 而高頻成分則表征了信號的細節(jié)差別。小波分析常用到近似和細節(jié)。近似能表示信號的大尺度,即低頻成分;而細節(jié)則表示信號的小尺度,即高頻成分,所以原始振動信號可以用兩個互補濾波器產(chǎn)生兩個信號。其一階濾波可簡要表示如下圖 3.1:
圖 3.1 信號一階濾波示意圖
對一個振動信號采用上圖 3.1 所示的濾波方法后,理論上會產(chǎn)生兩倍于原始信號的數(shù)據(jù),可以在整個采樣范圍內(nèi)僅在固定間隔內(nèi)采集數(shù)據(jù)減少采集的數(shù)據(jù)量即通過計算DWT 來獲得原始信號的細節(jié)和近似。
根據(jù)要求對所采集的近似振動信號進行連續(xù)分解能使信號在多個頻率段內(nèi)被分析。圖 3.2 的小波樹分解比較清晰的展示了小波分解的過程。圖中 S 表示原始信號 signal,A 表示信號近似 approximate ,D 表示信號細節(jié) detail, 而下標則代表信號的小波分解層數(shù)。
圖 3.2 小波分解樹示意圖
由于小波分解過程是不受其頻域和時域的限制的迭代,所以原則上說小波分解次數(shù)也不受限制。但當信號細節(jié)僅包括單個成分時繼續(xù)劃分會使分解小波失去其物理意義, 所以小波分解一般只根據(jù)對信號的研究需求選擇適當?shù)姆纸鈱訑?shù)。MATLAB 提供的離散小
波分解函數(shù)主要有:dwt(單層一維小波變換函數(shù))、dwt2(單層二維小波變換)、wavedec (多層一維小波分解函數(shù))、wavedec2(多層二維小波分解函數(shù))、dwtmode(設(shè)置離散小波和小波包變換信號或圖像擴展模式)、dwtper(周期化單層一維離散小波變換)、dwtper2(周期化單層二維離散小波變換)。
(2) 對小波系數(shù)作用閾值或掩碼;
在實現(xiàn)機械振動信號的處理過程中首要任務(wù)就是對分解后的各個小波系數(shù)進行適當?shù)拈撝底饔?。所以說閾值的合理選擇與振動信號的降噪效果息息相關(guān),這也吸引了大量學者提出了很多理論和模型,每種模型都有其適應(yīng)范圍。小波變換中,由各層的信噪比來確定所需的閾值。在確定噪聲的強弱程度后就能根據(jù)各層噪聲強度的不同來分別得到各層的作用閾值,一般來說根據(jù)噪聲強弱程度的不同而確定閾值的方式有如下三種:
1) 當采用缺省方式來計算時采用如下公式:
t h r=
3.2 3.2
這種情況下的閾值能通過 MATLAB 中的 ddencmp 指令來計算,其中的n 就是將被分解的振動信號的長度。
2) 由 Birge-Massart 算法來求解閾值時,可按照該步驟實現(xiàn):
首先確定一任意分解層數(shù) j,再依據(jù)確定好的分解層數(shù) j 進行下一步操作來保留更高層系數(shù)
再另選一 i 層(1 £ i £ j ),保存最大值 ni 個系數(shù)的情況下可用式 3.3 來計算:
i
n = M ( j + 2 -1)a
3.3
上式中的參數(shù) M、α 均是均是由以往信息處理所得經(jīng)驗而確定的系數(shù),當選取M = L(1) 為首層信號分解后的小波系數(shù)時,M 通常滿足 L(1) £ M £ 2L(1) ,α 的值依具體信號的特點和要求而定,在對信號進行降噪處理時一般取a = 3 。
3) 由小波包變換來確定的閾值由公式 3.4 計算而得,令 t*為使得函數(shù):
k
crit(t) = -?c2 + 2s 2t(a + log( n / t))
3.4
k £t
當式 3.4 中選擇最小的 t 值,而ck 表示振動信號經(jīng)過小泡包分解后的第 k 大的小波包系數(shù),這種情況下的閾值由公式 3.5 給出:
t h r=| c * |
t 3.5
式 3-11 中的 σ 表示噪聲信號的強弱程度,α 是依據(jù)信號降噪經(jīng)驗而確定的大于1 系數(shù)且必須為實數(shù),α 的取值大小對直接影響降噪信號經(jīng)小波分解后的系數(shù)緊湊程度,α 越大緊湊程度越低,重構(gòu)之后的降噪信號曲線也更趨于光滑,α 一般為 2。
通過適當?shù)姆椒ǐ@得所求閾值后,一般有兩種方式可用來處理分解后的小波, 硬閾值和軟閾值。前者的作用方式是把所有絕對值都不大于的數(shù)值定義為零點閾值;而后者是在前者的條件下把邊界位置的斷點均縮到零,這種閾值的選擇方式能減少信號在分解過程中的中斷。
(3) 重建小波系數(shù),得到處理之后的信號。
將分解并閾值處理后的信號在沒有損失的情況下還原原始信號的過程是小波重構(gòu), 這個步驟也被叫做逆離散小波變換(IDWT)。重構(gòu)過程通過小波系數(shù)來完成。信號降噪過程中的小波分解由濾波與下采樣,小波重構(gòu)濾波和上采樣組成。重構(gòu)過程就是小波分解示意圖的反過程。
3.2 振動信號降噪的小波包分析原理
小波包分析處理信號處理含噪信號時的基本思路及運算過程與小波分析大體相同, 區(qū)別主體體現(xiàn)在小波包分析能提供更靈活的信號分解方法。因為小波包分析中的信號細節(jié)與信號近似同樣能再分解,對于 N 層分解會產(chǎn)生 2N 個路徑,這決定了它在每個頻域段內(nèi)都能獲得比較精確的分析能力,圖 3.3 為小波包分解樹。
圖 3.3 小波包分解樹
與小波分解類似,可供信號的處理而使用的小波包基有很多種。小波包基的選擇方式一般依據(jù)所需處理信號的具體要求而定,其中處理效果最好的小波包基叫做最優(yōu)基。最優(yōu)基是以熵為選擇標準來確定。使用 MATLAB 進行小波包信號降噪處理時選取besttree() 函數(shù)來確定最優(yōu)基進,也稱作計算最優(yōu)樹。小波包在處理含噪聲的振動信號時通常分為以下四個步驟:
(1) 信號的小波包分解;
(2) 確定最優(yōu)小波包基;
(3) 小波包分解系數(shù)的閾值量化;
(4) 信號的小波包重構(gòu)。
與小波分析的核心相似,小波包分析要解決的核心問題也是合理的選取閾值,它一定程度上決定了小波包的振動信號降噪結(jié)果。MATLAB 工具箱中的 wpdencmp()專門用來小波包降噪和壓縮處理的函數(shù)。
3.3 小波分析的降噪模擬
對一個正弦信號 N = 1000;t = 1: N; x = sin(0.03t) 在包含方差為0.5的白噪聲進行降噪處理。
3.3.1 降噪過程
(1) 按要求生成正弦信號并用 randa 函數(shù)加入噪聲信號,把相加結(jié)果放入 ns 中;
(2) 將原始正弦信號,加噪信號以樣本點數(shù)為橫坐標,幅值為縱坐標分別圖示之;
(3) 將正弦信號與加噪信號分析進行快速傅里葉變換并取絕對值從而得出的正弦信號與加噪信號的功率譜圖示之;
(4) 用 wden 通過方法對 ns 一維信號進行小波降噪;使用 db3 型小波,閾值選擇極大極小值閾值選取方法,采用軟閾值的閾值化方法,one 型閾值化分辨率分 5 層對小波進行降噪處理,并將降噪結(jié)果存入 xd。
(5) 分析對降噪后的信號和降噪后經(jīng)過步驟(3)作用的結(jié)果圖示之。
3.3.2 降噪效果
功率譜選擇從 0 到 100 為例如圖 3.4:‘從含噪信號功率譜’和‘消噪信號功率譜’ 圖可以直觀的看出降低噪聲之后的信號更加集中,信號的能量全都集中在頻率為 5Hz
左右,而含噪信號則顯的不那么直觀,雖然它在 5Hz 的部分能量也比較集中,但是在更高的頻率下依然有信息,會造成對振動信號包含信息的干擾。
從‘含噪信號’和‘消噪信號’的對比中可以后看去噪后的信號基本上恢復了原始信號形狀,除去了噪聲干擾,但與‘原始信號’存在區(qū)別,這主要是由選取小波的不同和細節(jié)系數(shù)閾值的選擇相關(guān)。
圖 3.4 含噪正弦信號降噪效果圖
3.3.3 降噪代碼
N=1000; t=1:N;
x=sin(0.03*t); noise=0.5*randn(1,1000); ns=x+noise; subplot(3,2,1);
plot(t,x);
title('原始信號');
px=abs(fft(x));
subplot(3,2,2); plot(px);
title('原始信號功率譜'); subplot(3,2,3);
plot(ns);
title(' 含 噪 信 號 '); pns=abs(fft(ns)); subplot(3,2,4);plot(pns);
title(' 含 噪 信 號 功 率 譜 '); xd=wden(ns,'minimaxi','s','one',5,'db3'); subplot(3,2,5);
plot(xd);
title('降噪信號'); pxd=abs(fft(xd)); subplot(3,2,6); plot(pxd);
title('降噪信號功率譜')
3.4 閾值選擇對小波分析影響的降噪模擬
下面通過一個實例說明閾值選擇對信號的保持及信號降噪效果的影響[21]:
3.3.1 降噪過程
(1) 載入matlab系統(tǒng)的的含噪leleccum信號,取3920個采樣點,將原始噪聲信號圖示之;
(2) 對db1型小波對噪聲信號進行小波分解,分解到分三層,其細節(jié)系數(shù)和近似系數(shù)分別存入相應(yīng)的ca,cd;
(3) 直接對分解的信號進行強制重構(gòu),并將結(jié)果圖示之;
(4) 采用ddncmp函數(shù)獲得信號的默認閾值:‘den’表示降噪,‘wv’表示小波分解,‘s’表示軟閾值的作用方式;并且用wdencmp命令函數(shù)實現(xiàn)信號的降噪過程:‘gbl’表示每層均用同一閾值進行處理,‘c,l’表示將要被處理降噪的輸入信號,‘db1’是選用小波的名稱;將降噪結(jié)果存入s2并圖示之;
(5) 采用wthresh命令對分解細節(jié)系數(shù)向量作用閾值:‘s’表示軟閾值,并將結(jié)果返回給各細節(jié)系數(shù),對作用閾值后的細節(jié)系數(shù)和近似系數(shù)進行重構(gòu)并將其降噪結(jié)果圖示之。
3.3.2 降噪結(jié)果
從圖3.5可以看出給定軟閾值消噪對信號的保留做的較好,它的閾值一般是由經(jīng)驗公式獲得,信號的有用成分不易丟失,可信度比較高;強制消噪的曲線比較圓滑,但容易丟失有用信號,這也與作用閾值的原理相符。
圖 3.5 含噪正弦波去降噪結(jié)果
3.3.3 降噪代碼: load leleccum; s=leleccum(1:3920); ls=length(s); subplot(2,2,1); plot(s);
title('原始信號');grid; [c,l]=wavedec(s,3,'db1');
ca3=appcoef(c,l,'db1',3); cd3=detcoef(c,l,3); cd2=detcoef(c,l,2); cd1=detcoef(c,l,1); cdd3=zeros(1,length(cd3)); cdd2=zeros(1,length(cd2)); cdd1=zeros(1,length(cd1)); c1=[ca3 cdd3 cdd2 cdd1]; s1=waverec(c1,l,'db1');
subplot(2,2,2); plot(s1);
title(' 強 制 消 噪 后 的 信 號 ');grid; [thr,sorh,keepapp]=ddencmp('den','wv',s); s2=wdencmp('gbl',c,l,'db1',3,thr,sorh,keepapp); subplot(2,2,3);
plot(s2);
title('默認閾值消噪后的信號');grid; cd1soft=wthresh(cd1,'s',1.465); cd2soft=wthresh(cd2,'s',1.823); cd3soft=wthresh(cd3,'s',2.768); c2=[ca3 cd3soft cd2soft cd1soft]; s3=waverec(c2,l,'db1');
subplot(2,2,4); plot(s3);
title('給定軟閾值消噪后的信號');grid
3.5 小波包分析的降噪模擬
用小波包來分析對模擬對一個系統(tǒng)內(nèi)置的含噪信號noismima 進行去噪處理的過程。
3.4.1 降噪過程
(1) 裝載 noismima 原始信號并圖示之,采用快速傅里葉變換并取絕對值獲得功率譜并圖示之;
(2) 采用默認閾值,wdencmp()函數(shù)進行去噪同時使用全局閾值的選項對信號進行去噪處理并圖示之,同(1)中的方法作功率譜圖;
(3) 根據(jù)之前的降噪效果,調(diào)節(jié)閾值的大小先進一步進行降噪,圖示并作功率譜圖;
3.4.2 降噪效果
從圖 3.6 可以明顯的看出,調(diào)節(jié)閾值后的小波包去噪更好的去除了原始信號中的噪聲,功率譜圖以 0-500Hz(即一半)的頻率為例說明,調(diào)節(jié)閾值后的消噪信號基本沒有了高頻部分的干擾,所以其信號的能量更加集中,從調(diào)節(jié)閾值后的水口信號圖中也可以體現(xiàn)這一點。
圖 3.6 小波包降噪效果圖
3.4.3 降噪代碼
load noismima; s=noismima(1:1000); figure(1); subplot(3,2,1); plot(s);
xlabel('樣本序號, n'); ylabel('幅值 A');
title('原始信號');
px=abs(fft(s));
subplot(3,2,2); plot(px);
[thr,sorh,keepapp,crit]=ddencmp('den','wp',s); [c,treed,perf0,perfl2]= wpdencmp(s,sorh,3,'db2',crit,thr,keepapp); subplot(3,2,3);
plot(c);
xlabel('樣本序號 n');
ylabel('幅值 A');
title('默認閾值消噪信號'); pns=abs(fft(c));
subplot(3,2,4); plot(pns);
thr=thr+15;
[c1,treed,perf0,perfl2]= wpdencmp(s,sorh,3,'db2',crit,thr,keepapp); subplot(3,2,5);
plot(c1);
xlabel('樣本序號 n');
ylabel('幅值 A');
title('調(diào)節(jié)閾值后的消噪信號'); pxd=abs(fft(c1)); subplot(3,2,6);
plot(pxd);
3.6 小波及小波包分析的工程信號降噪實例
3.6.1 小波降噪的應(yīng)用實例
圖3.7反映了小降噪對信號選擇性強的優(yōu)勢,左面兩圖原始信號和降噪信號的的幅值
收藏