ICEEMD和HD的單向閥早期故障信號降噪方法

錢恩麗,黃國勇*,,2,何冬,李鍶宇

(1. 昆明理工大學 信息工程與自動化學院,昆明 650500;2. 昆明理工大學 民航與航空學院,昆明 650500)

往復式高壓隔膜泵是礦漿管道輸送的核心動力設備,單向閥是高壓隔膜泵中運動最頻繁的核心部件之一,這使單向閥較泵的其他部件容易出現故障[1]。實際工程中采集的單向閥振動信號伴隨大量背景噪聲,從而導致單向閥早期故障振動特征難以有效提取。由于單向閥故障沖擊成分能量占振動總能量比例小,單向閥早期故障難以檢測和診斷,故對單向閥早期故障分析時需先進行降噪處理。

對機械故障信號降噪時多采用自適應分解降噪方法,其中經驗模態(tài)分解(Empirical mode decomposition, EMD)在機械故障診斷領域得到廣泛應用,但EMD存在模態(tài)混疊等缺點,大量EMD的改進方法被提出,文獻[2]將集合經驗模態(tài)分解(Ensemble empirical mode decomposition, EEMD) 與相關系數分析相結合,對滾動軸承進行故障診斷,EEMD對經EMD獲得的模態(tài)進行平均,進一步改善EMD存在模態(tài)混疊問題,但EEMD的完備性較差。Torres等[3]提出自適應噪聲完全集合經驗模態(tài)分解(Complete ensemble empirical mode decomposition with adaptive noise, CEEMDAN),在EEMD分解的每個階段添加特定噪聲,并計算唯一殘差以獲得每種模式。該方法克服EEMD的不足,但CEEMDAN仍存在分解信號時會出現虛假分量的缺點。Colominas等[4]提出ICEEMD,對CEEMDAN進一步改進,該方法不僅保留EMD方法的優(yōu)點,且其模式中殘留噪聲小,重構信號誤差小且避免在分解過程中產生偽模態(tài)。ICEEMD將一個復雜信號分解為多個從高頻到低頻順序分布的IMF分量,大量噪聲被優(yōu)先分解到高頻IMF分量中,利用相關系數等方法選取有效IMF分量進行信號重構時,由于早期故障信號中含有大量噪聲致使含有大量噪聲的IMF分量(噪聲IMF分量)與原始信號的相關性更大,因此易優(yōu)先篩選到噪聲IMF分量,若直接把分解后的前若干個IMF分量當作噪聲分量剔除,又過于主觀,缺乏嚴謹的理論框架支撐,且很難恰當地分離噪聲分量。故需找到能較好地將ICEEMD分解得到IMF分量中的噪聲IMF分量分離的方法。基于幾何度量的豪斯多夫距離(HD)已被證明在分量選取方面具有良好的應用效果[5]。文獻[6]將變分模態(tài)分解與豪斯多夫距離相結合,證明其方法降噪的有效性且成功對管道泄漏信號降噪。文獻[7]將EMD與HD結合對信號去噪,通過仿真信號和實測陀螺靜態(tài)漂移數據驗證了方法的有效性。

鑒于HD在分量選取方面的優(yōu)越性,將其與ICEEMD相結合,彌補ICEEMD缺乏分離噪聲IMF分量依據的缺點,同時達到對單向閥原始信號降噪的目的。利用HD將噪聲IMF分量從ICEEMD分解得到的IMF分量中分離,以峭度為指標在剩余IMF分量中選取敏感IMF分量重構,將重構信號進行希爾伯特包絡解調,提出一種ICEEMD和HD的單向閥早期故障信號降噪新方法。

1 基本理論

1.1 改進的完備集合經驗模態(tài)分解

ICEEMD基于EMD分解原理,在提取第n個IMF分量時,加入特殊的噪聲En(w(i)),抑制模態(tài)混疊及分解過程中產生偽模態(tài)的缺陷,增強抗干擾性。ICEEMD具體分解步驟[4]如下:

1) 對原始信號x加入零均值的單位協(xié)方差的高斯白噪聲w(i)(i=1,2,…,Ι),得

x(i)=x+β0E1(w(i))

(1)

式中:En(·)代表EMD分解出的第n個IMF分量;為在所添加噪聲和添加了噪聲的殘差之間求得所需要的信噪比而選擇常數βn-1=εn-1std(rn-1),當n=1時,β0=ε0std(x)/std(E1(w(i))),其中ε0為第一個添加的噪聲和分析信號之間所需信噪比的倒數,一般取ε0=0.2。

2) 通過EMD對x(i)迭代i次計算x(i)的局部均值,從而得到第一個殘差r1和IMF1分量為:

r1=〈M(x)(i)〉

(2)

IMF1=x-r1

(3)

式中:〈·〉是求平均值運算符;M(·)是求分析信號的局部平均值運算符。

3) 第二個殘差是r1+β1E2(w(i))的平均局部均值,因此,r2和IMF2為:

r2=〈M(r1+β1E2(w(i)))〉

(4)

IMF2=r1-r2

(5)

4) 設n=3,4,…,m時,計算第n個殘差為

rn=〈M(rn-1+βn-1En(w(i)))〉

(6)

5) 計算第n個IMFn分量為

IMFn=rn-1-rn

(7)

ICEEMD能自適應地將一個復雜信號分解為按高頻到低頻順序分布且具有不同特征的多個IMF分量,其中包括噪聲分量,有部分噪聲的分量以及有效信號分量[8]。由于大部分隨機噪聲表現為高頻噪聲,故大部分噪聲會被優(yōu)先分解到高頻分量中。單向閥早期故障信號中含有大量噪聲,因此需要先將ICEEMD分解得到的噪聲IMF分量從所有IMF分量中分離。

1.2 豪斯多夫距離

豪斯多夫距離是描述兩個點集之間距離的方法,該方法常用來衡量同一個空間內兩個點集間相似程度[9]。其具體定義如下:

A={a1,…,am}和B={b1,…,bn}分別為有限點的數據集合。令d(a,b)為兩個點之間的歐氏距離,那么點a到數據集B的距離為

(8)

數據集A到數據集B的直接距離為

(9)

同樣,數據集B到數據集A的直接距離為

(10)

因此,兩個數據集之間的豪斯多夫距離定義為

HD(A,B)=max(D(A,B),D(B,A))

(11)

表示數據分布形狀的概率密度函數(Probability density function,PDF),可反映兩種信號的差異[10]。對原始信號求PDF可以保留原始信號完整的信息,可以根據輸入原始信號的PDF與各IMF分量的PDF之間的相似性來選擇相關IMF分量。故利用豪斯多夫距離度量ICEEMD分解得到IMF分量與原始信號的概率相似性,從而將噪聲IMF分量分離,具體式為

hd(i)=HD[PDF(x(t)),PDF(IMFi(t))]

(12)

計算每個IMF分量概率密度函數與原信號概率密度函數之間的豪斯多夫距離,以HD為依據分離噪聲IMF分量。噪聲IMF分量和有效IMF分量的分界點定義為:定位得到第一個豪斯多夫距離的局部最大值,在其后豪斯多夫距離開始減小時所對應的IMF分量即為分界點[7]?；谝陨戏治?利用HD定位出噪聲IMF分量與有效IMF分量的分界點,以此分離ICEEMD分解出IMF分量中的噪聲IMF分量。

1.3 峭度

峭度對沖擊脈沖及脈沖類故障信號十分敏感,適用于表面損傷類故障,尤其是早期故障診斷,其計算公式[11-12]為

(13)

式中:K為信號x的峭度指標;N為信號x的長度;μ為信號x的均值;σ為信號x的標準差。

單向閥發(fā)生早期故障時會產生沖擊成分,因此以峭度為指標,從去除噪聲IMF分量的剩余IMF分量中篩選峭度值較大的敏感IMF分量進行重構。目前IMF分量選擇個數沒有統(tǒng)一標準,經多次試驗分析,本文以峭度為指標選取剩余IMF分量個數的一半進行重構。

1.4 希爾伯特包絡譜

設x(t)為時域信號,其Hilbert變換定義為

(14)

則原始信號x(t)和它的Hilbert變換信號xh(t)可以構成一個新的解析信號xa(t),即

xa(t)=x(t)+jxh(t)

(15)

其幅值為

(16)

解析信號的幅值為實信號的包絡,對包絡進行快速傅里葉變換得到希爾伯特包絡譜[13]。求解重構信號的希爾伯特包絡譜,得到信號的頻率信息。

2 ICEEMD和HD的單向閥早期故障信號降噪方法

根據上述理論分析,本文提出了ICEEMD和HD的單向閥早期故障信號降噪方法。其降噪方法流程如圖1所示。

圖1 降噪方法流程圖

具體步驟如下:

1) 對所采集單向閥早期磨損故障數據進行ICEEMD分解,并獲得一系列IMF分量。

2) 計算每個IMF分量及原始信號的概率密度函數,計算每個IMF分量概率密度函數與原信號概率密度函數之間的豪斯多夫距離,并定位豪斯多夫距離的局部最大值,在其后豪斯多夫距離開始減小時所對應的IMF分量為分界點,分界點之前的IMF分量視為噪聲IMF分量,將其分離。

3) 計算所有IMF分量的峭度值,在剩余IMF分量中選取一半峭度值較大的IMF分量進行重構。

4) 利用希爾伯特包絡譜對重構信號進行包絡解調,提取單向閥運行基頻及倍頻,觀察降噪效果。

3 仿真信號分析

為驗證所提方法的有效性,構造仿真信號,其表達式為

s(t)=Asin(2πf0t)+ε(t)

(17)

式中:信號幅值A=0.02;信號頻率f0=1 Hz;ε(t)為SNR=4 dB的高斯白噪聲,采樣頻率為2 560 Hz,采樣點數為10 240。

仿真信號的時域波形和頻譜如圖2所示,特征頻率完全被噪聲淹沒,采用本文所提方法對仿真信號進行降噪。

圖2 仿真信號的時域波形和頻譜圖

仿真信號經ICEEMD自適應分解得到的IMF分量如圖3所示,仿真信號被自適應分解為13個從高頻到低頻排列的IMF分量和一個殘差分量。

圖3 仿真信號IMF分量時域波形圖

計算每個IMF分量概率密度函數與仿真信號概率密度函數之間的豪斯多夫距離,如圖4所示。定位到HD第一個局部最大值IMF9,其后HD開始減小時所對應的IMF分量為IMF10,以IMF10分量為界,IMF1～IMF9分量為要分離的噪聲分量。

圖4 仿真信號PDF與各IMF分量PDF間HD值

ICEEMD分解仿真信號得到IMF分量的峭度值如表1所示,從IMF10分量開始到IMF13分量選取峭度值較大的2個即IMF10、IMF12分量重構。

表1 仿真信號各IMF分量峭度值

利用希爾伯特包絡譜對重構信號進行包絡解調,降噪效果如圖5所示。仿真信號含有大量噪聲,致使降噪后提取仿真信號的特征頻率存在誤差,在誤差允許范圍內,經本文方法降噪后可有效提取出仿真信號被噪聲淹沒的特征頻率1 Hz的近似頻率1.09 Hz,且消除了仿真信號中大量噪聲,驗證了本文降噪方法的有效性。

圖5 仿真信號希爾伯特包絡解調圖

4 試驗分析

4.1 試驗說明

利用所提降噪方法對采集的云南大紅山管道三號高壓隔膜泵站單向閥早期故障信號進行降噪處理。高壓隔膜泵為TZPM系列的三缸曲軸驅動活塞式隔膜泵,最高工作壓力為24.44 MPa。將PCB-ICP型振動加速度傳感器安裝在同一泵組的三對進、出料單向閥上,同時用型號為PCB-394C06的加速度校準器校準檢測加速度信號,傳感器安放位置及故障單向閥如圖6所示,用PXIe-3342高精度8通道動態(tài)數據采集卡采集振動信號數據,采樣頻率為2 560 Hz,通過電腦實時監(jiān)測采集的信號。本文對采集的單向閥早期磨損故障的振動信號進行降噪處理,其采樣數據長度為10 240。由現場實測數據可知高壓隔膜泵的沖程數為每分鐘30～31次,其正常運行時頻率為0.5～0.517 Hz。文獻[14]進一步計算得出單向閥運行的基頻為1 Hz～1.034 Hz。

圖6 傳感器安放位置及故障單向閥

4.2 試驗驗證

單向閥早期磨損故障信號的時域波形和頻譜如圖7所示。

圖7 早期磨損故障信號的時域波形和頻譜圖

由時域波形圖可知單向閥早期磨損時運行特征頻率被噪聲完全淹沒;由頻譜圖可知單向閥發(fā)生早期磨損時能量主要集中在0～250 Hz,且信號含有大量噪聲,不能直接提取到單向閥運行的基頻及其倍頻,需對信號進行降噪處理。

單向閥早期磨損信號經ICEEMD自適應分解得到的IMF分量如圖8所示,實驗中添加200組零均值的單位協(xié)方差的高斯白噪聲,其標準差為原始信號標準差的0.2倍。原始信號被自適應分解為11個從高頻到低頻排列的IMF分量和一個殘差分量,高頻IMF分量中含有的噪聲成分高于低頻IMF分量,且含有用信息少,故需先將噪聲IMF分量分離。

圖8 早期磨損故障信號IMF分量時域波形圖

據前文所述降噪過程,計算每個IMF分量概率密度函數與原始信號概率密度函數之間的豪斯多夫距離,如圖9所示。

圖9 ICEEMD分解各IMF分量PDF與原信號PDF間HD值

定位到HD第一個局部最大值IMF3,其后HD開始減小時所對應的IMF分量為IMF4,以IMF4分量為界,IMF1～IMF3分量為要分離的噪聲分量。ICEEMD分解得到IMF分量的峭度值如表2所示,從IMF4分量開始到IMF11分量選取峭度值較大的4個即IMF6、IMF7、IMF9、IMF10分量重構。

表2 早期磨損故障信號IMF分量峭度值

利用希爾伯特包絡譜對重構信號進行包絡解調,降噪效果如圖10所示。從圖10可以明顯看到20 Hz以后的大量噪聲被消除,同時保留了20 Hz以前的大量故障沖擊成分。進一步觀察0～20 Hz的頻率成分,可以清楚得到與高壓隔膜泵運行基頻相近的頻率成分0.625 Hz,與單向閥運行基頻相近的頻率成分0.937 5 Hz及其2～9倍頻。因圖片大小限制,在所有希爾伯特包絡解調圖中僅標出前3個頻率具體數值,其他頻率均只標記出所在位置。

圖10 所提方法希爾伯特包絡解調圖

4.3 對比分析

為驗證本文降噪方法的有效性及優(yōu)越性,在原始數據相同的前提下,將其與以下幾種降噪方法進行比較:

1) 為驗證ICEEMD分解方法的有效性,將本文方法與VMD+HD+峭度方法進行對比。將本方法中的ICEEMD用VMD替換,其它降噪步驟與本文方法相同。VMD分解時需要確定分解模態(tài)數和懲罰因子,借鑒文獻[15]的選取方法,本文選取分解模態(tài)數為6,懲罰因子為2000。

由于VMD將原始信號分解為從低頻到高頻的一系列IMF分量,噪聲優(yōu)先分解到高頻IMF分量中,因此通過HD定位噪聲IMF分量分界點時,從高頻IMF分量開始定位HD第一個局部最大值,其后HD開始減小時所對應的IMF分量即為分界點,將噪聲IMF分量分離。VMD分解得到的每個IMF分量的概率密度函數與原始信號的概率密度函數之間的豪斯多夫距離如圖11所示,HD第一個局部最大值為IMF5分量,其后HD開始減小時所對應的IMF分量為IMF4,以IMF4分量為分界點,IMF5和IMF6分量為噪聲IMF分量。

圖11 VMD分解各IMF分量PDF與原信號PDF間HD值

在剩余IMF分量中選取峭度值較大的兩個IMF分量即IMF3和IMF4分量進行重構,利用Hilbert包絡譜對重構信號進行解調,其降噪效果如圖12所示。圖中20 Hz以后還存在大量噪聲的干擾。進一步觀察0～20 Hz的頻率成分,該方法降噪后只提取到了高壓隔膜泵運行的近似基頻0.625 Hz和單向閥運行近似基頻0.937 5 Hz及其2、3、4、5、7、9倍頻,且在0～20 Hz之間仍存在噪聲成分的干擾。故ICEEMD對原始信號具有較好分解降噪效果,且優(yōu)于VMD分解方法。

圖12 VMD+HD+峭度希爾伯特包絡解調圖

2) 為驗證HD對分離ICEEMD分解得到的IMF分量中噪聲IMF分量的有效性,將本文所提方法與ICEEMD+峭度方法進行對比。ICEEMD+峭度方法不再進行噪聲IMF分量分離而直接選取ICEEMD分解出的IMF分量中峭度值較大的IMF分量進行重構,為了保證對比變量的唯一性,選取峭度值較大的4個IMF分量進行重構,然后將重構信號進行希爾伯特包絡解調。由表2中的峭度值可知選取IMF1、IMF2、IMF9、IMF10分量進行重構,其降噪效果如圖13所示。

圖13 ICEEMD+峭度希爾伯特包絡解調圖

由前面分析可知IMF1、IMF2分量為噪聲IMF分量,因此該方法降噪后,在20 Hz以后依然存在大量噪聲的干擾。進一步觀察0～20 Hz的頻率成分,該方法只提取到高壓隔膜泵運行的近似基頻0.625 Hz和單向閥運行頻率的2、3、4、7、9倍頻,單向閥運行基頻被噪聲淹沒。故利用HD可有效分離ICEEMD分解所得IMF分量中的噪聲IMF分量。

3) 為驗證利用峭度指標對剩余IMF分量進行篩選的必要性,將本文所提方法與ICEEMD+HD方法進行對比。ICEEMD+HD方法將分離出噪聲IMF分量后的所有IMF分量進行重構,然后對重構信號進行希爾伯特包絡解調,其降噪效果如圖14所示。該方法降噪后在20～100 Hz之間存在大量噪聲干擾。進一步觀察0～20 Hz的頻率成分,該方法只提取到高壓隔膜泵運行的近似基頻0.625 Hz和單向閥運行頻率的2、3、5、6、7、8、9倍頻。雖然利用HD將含有大量噪聲的IMF分量分離,但剩余的IMF分量中仍然含有部分噪聲,將剩余IMF分量全部重構就會將剩余噪聲全部聚集,增加重構信號噪聲含量,進而出現圖14中存在部分噪聲干擾甚至單向閥運行基頻被噪聲淹沒的現象。故利用峭度指標對剩余IMF分量進行篩選十分必要。

圖14 ICEEMD+HD希爾伯特包絡解調圖

4) 為驗證本文所提降噪方法的優(yōu)越性,將本文方法與文獻[1]方法進行對比。文獻[1]所提方法可提取到單向閥運行的基頻及2倍頻和4倍頻,本文所提方法可有效提取高壓隔膜泵的運行基頻,單向閥的運行基頻及其2至9倍頻,故本文所提方法具有更優(yōu)的降噪效果。

經上述定性分析,本文所提方法較以上3種降噪方法和文獻[1]所提方法具有明顯降噪優(yōu)勢。表3計算了本文所提方法和其他3種方法重構信號的樣本熵[16],進一步定量分析本文降噪方法的優(yōu)越性。

由表3可知,本文所提降噪方法的樣本熵明顯低于其他3種方法,證明所提方法降噪后信號復雜度最低,所含噪聲最少。對比方法3和4可知,方法4的樣本熵低于方法3,說明方法4重構信號所含噪聲低于方法3。由于方法3重構了IMF1、IMF2、IMF9、IMF10分量,由前面分析可知IMF1和IMF2分量為噪聲IMF分量,因此方法3降噪后信號最復雜,所含噪聲最多,樣本熵最大,證明利用HD分離噪聲IMF分量的必要性。通過定量分析,進一步證明本文所提方法對單向閥早期磨損故障信號降噪的有效性和優(yōu)越性。

表3 不同降噪方法樣本熵

5 結論

針對往復式高壓隔膜泵單向閥早期故障振動信號特征信息被噪聲淹沒的問題,提出ICEEMD和HD相結合的單向閥早期故障信號降噪方法。將該方法與其他降噪方法進行對比,從試驗結果和對比分析得到以下結論。

1) 利用本文所提方法對單向閥早期磨損故障信號進行降噪,該方法可以消除信號中大量噪聲成分,保留故障特征信息,有效提取單向閥早期故障信號中單向閥運行的基頻及其倍頻,為單向閥早期故障信號降噪提供一種新方法。

2) 本文方法采取分離噪聲IMF分量以及選取剩余IMF分量中敏感IMF分量的方法進行降噪,這樣會剔除掉部分有用信息,后續(xù)將考慮對所有不同頻段IMF分量依據其不同特性進行具體降噪,不再進行IMF分量篩選,充分保留所有IMF分量中有用信息的同時實現單向閥早期故障信號降噪。