HHT和馬氏距離融合的船用空壓機(jī)故障診斷

胡歡歡，王永堅(jiān)，邱 晨

(集美大學(xué)輪機(jī)工程學(xué)院， 福建 廈門 361021)

0 引言

舶用二級(jí)往復(fù)式空壓機(jī)是船舶壓縮空氣系統(tǒng)的關(guān)鍵設(shè)備，其工作性能的好壞對(duì)船舶正常營運(yùn)將產(chǎn)生影響，開展船用空壓機(jī)故障類型的有效識(shí)別和故障診斷具有良好的現(xiàn)實(shí)意義。船用二級(jí)往復(fù)式空氣壓縮機(jī)最常見的故障是一、二級(jí)缸套-活塞組發(fā)生斷環(huán)，依靠空壓機(jī)性能參數(shù)很難及時(shí)、準(zhǔn)確地查找故障問題，本文提出一種基于振動(dòng)信號(hào)分析的船用二級(jí)往復(fù)式空氣壓縮機(jī)活塞環(huán)斷環(huán)故障診斷方法。由于空壓機(jī)振源多且機(jī)械結(jié)構(gòu)復(fù)雜，振動(dòng)信號(hào)是非線性非平穩(wěn)的，如何從此類信號(hào)中提取故障特征將是一件十分關(guān)鍵的工作。傳統(tǒng)的振動(dòng)信號(hào)處理方法在處理非線性非平衡信號(hào)時(shí)存在缺陷，如短時(shí)傅里葉變換因高低頻要求不同的窗函數(shù)，合適的窗函數(shù)不易選取，Wigner分布會(huì)產(chǎn)生二次混疊，Cohen類很難選擇合適的核函數(shù)，小波變換不具有自適應(yīng)性，難以選擇小波基[1-2]。希爾伯特-黃變換(hilbert-uuang transform，HHT)使用經(jīng)驗(yàn)?zāi)B(tài)分解算法，將信號(hào)分解成一組本征模態(tài)函數(shù)[3-5]，此類函數(shù)具有單量自適應(yīng)性，高時(shí)頻分辨率和良好的時(shí)頻聚集性的優(yōu)點(diǎn)，故可以處理非線性非平穩(wěn)信號(hào)。對(duì)重構(gòu)后的信號(hào)進(jìn)行希爾伯特變換，可以得到信號(hào)瞬時(shí)時(shí)頻分布，進(jìn)而可以得到Hilbert邊際譜[6-7]，它可以準(zhǔn)確地反映出信號(hào)幅值和頻率之間的關(guān)系。馬氏距離依據(jù)數(shù)據(jù)的協(xié)方差距離，具有良好的小樣本分類特性，可以準(zhǔn)確地反映兩個(gè)樣本的相似程度[8]，因此在機(jī)械故障識(shí)別領(lǐng)域得到廣泛應(yīng)用。

綜上，本文提出一種把希爾伯特-黃邊際譜和馬氏距離融合的船用二級(jí)往復(fù)式空壓機(jī)斷環(huán)故障診斷方法，通過搭建實(shí)驗(yàn)平臺(tái)(以船用應(yīng)急空壓機(jī)為實(shí)驗(yàn)對(duì)象)，人為模擬空壓機(jī)三種狀態(tài)類型(正常狀態(tài)、一級(jí)活塞斷環(huán)、二級(jí)活塞斷環(huán))進(jìn)行實(shí)驗(yàn)，驗(yàn)證所提方法的有效性和故障識(shí)別的準(zhǔn)確性。

1 信號(hào)分析方法

經(jīng)驗(yàn)?zāi)B(tài)分解和希爾伯特變換及譜分析兩個(gè)部分組成希爾伯特-黃變換(HHT)。經(jīng)驗(yàn)?zāi)B(tài)分解(EMD)算法得到本征模態(tài)函數(shù)(IMF)必須滿足一下兩個(gè)條件：1)極值點(diǎn)和零點(diǎn)的個(gè)數(shù)相差不多于一個(gè):2)極大值和極小值點(diǎn)的包絡(luò)線的平均值等于零[5]。

EMD分解實(shí)測(cè)信號(hào)x(t)的步驟如下。

步驟一應(yīng)用三次樣條算法求振動(dòng)信號(hào)x(t)所有的局部極大值、極小值點(diǎn)的包絡(luò)線，即上、下包絡(luò)線，然后求取這兩個(gè)包絡(luò)線的平均值，記為m1，將x(t)減去m1，得到h1=x(t)-m1。

步驟二如果h1滿足上述的條件，那么h1就是第一個(gè)本征模態(tài)函數(shù)(IMF1)；如果h1不滿足條件，則把h1作為原始數(shù)據(jù)，重復(fù)步驟一，即先得到上、下包絡(luò)線的平均值m11,h1減去m11得到h11，再判斷是否滿足條件，如果不滿足，則重復(fù)循環(huán)k次，直到h1k滿足條件，并記c1=h1k為第一個(gè)本征模態(tài)函數(shù)，h1k=h1(k-1)-m1k。

步驟三把第一個(gè)本征模態(tài)函數(shù)c1從原始信號(hào)x(t)中分離出來，從而得到r1=x(t)-c1。

步驟四將r1當(dāng)作原始數(shù)據(jù)，重復(fù)上面的三個(gè)步驟，得出x(t)的第2個(gè)滿足條件的IMF分量c2(IMF2)。如此重復(fù)n次，得出x(t)的n個(gè)IMF分量，這樣就有如下式子：

當(dāng)余下的信號(hào)函數(shù)不再滿足提取IMF分量的條件時(shí)，循環(huán)結(jié)束，余下的信號(hào)函數(shù)稱為殘余項(xiàng)rn。

于是原始信號(hào)x(t)可以表示為n個(gè)IMF分量和一個(gè)殘余項(xiàng)rn之和，這樣就重構(gòu)了原始信號(hào)，即為：

于是求得每一個(gè)IMF分量的幅值和相位函數(shù)分別為：

φi(t)=arctan [H[ci(t)]/ci(t)]。

進(jìn)而求得瞬時(shí)頻率ωi(t)=(1/2π)dφi(t)/dt。

由于殘余項(xiàng)rn對(duì)重構(gòu)原始信號(hào)影響較小，故省略不計(jì)，從而得到重構(gòu)的原始信號(hào)為：

Re為實(shí)部，n為IMF分量的個(gè)數(shù)。由此式可以看出幅值ai(t)和頻率ωi(t)都是時(shí)間t的函數(shù)，從而得到Hilbert譜，記為：

(1)

進(jìn)而得在某一頻段ω1～ω2的能量函數(shù)S為：

(2)

2 實(shí)驗(yàn)數(shù)據(jù)采集

本實(shí)驗(yàn)臺(tái)選用應(yīng)急船用往二級(jí)復(fù)式空壓機(jī)，額定轉(zhuǎn)速為2800 r/min,額定功率為7.5 kW,額定壓力12.5 MPa;儲(chǔ)氣罐額定壓力為12.5 MPa；采用東華測(cè)試公司的信號(hào)采集儀DH-5922，以及相關(guān)的DHDAS信號(hào)測(cè)試分析系統(tǒng)軟件；傳感器采用加速度傳感器，最大可測(cè)加速度為50g m/s2，靈敏度為100.8，激勵(lì)電壓為24 V，諧振頻率為30 kHz。

將4個(gè)加速度傳感器分別安裝于一、二級(jí)氣缸蓋頂端和一、二級(jí)缸套側(cè)面中間的位置。搭建的實(shí)驗(yàn)平臺(tái)如下圖1所示。本次采用人為模擬故障問題實(shí)驗(yàn)，設(shè)正常狀態(tài)、一級(jí)活塞斷環(huán)、二級(jí)活塞斷環(huán)3種故障，如圖2所示。

采集正常狀態(tài)下的振動(dòng)信號(hào)，具體過程如下：開啟DH-5922信號(hào)采集儀，再打開DHDAS信號(hào)測(cè)試分析系統(tǒng)軟件，將采集頻率調(diào)節(jié)到12 kHz,采樣點(diǎn)數(shù)為25 600個(gè)，數(shù)據(jù)格式為TXT，后啟動(dòng)空壓機(jī)，當(dāng)空壓機(jī)處于工作狀態(tài)時(shí)，按下信號(hào)采集儀采集信號(hào)鍵，直到信號(hào)采集點(diǎn)數(shù)完成，停止采集信號(hào)，得到正常狀態(tài)振動(dòng)信號(hào)。將一級(jí)斷裂的活塞環(huán)安裝于對(duì)應(yīng)的活塞上，重復(fù)上述采集信號(hào)的步驟，得到一級(jí)活塞斷環(huán)的故障振動(dòng)信號(hào)數(shù)據(jù)；同樣可得二級(jí)活塞斷環(huán)狀態(tài)下的振動(dòng)信號(hào)數(shù)據(jù)。本次實(shí)驗(yàn)每種狀態(tài)各采集10組數(shù)據(jù)，共40組數(shù)據(jù)，三種狀態(tài)共120組數(shù)據(jù)。

3 數(shù)據(jù)分析

將采集到的數(shù)據(jù)使用MATLAB2013b軟件進(jìn)行處理，得到三種狀態(tài)下信號(hào)原始圖，如圖3所示。

由圖3能夠得出，由于背景噪聲和環(huán)境等因素的影響，所獲取的三種狀態(tài)下的信號(hào)存在著較多的干擾，這將嚴(yán)重影響故障診斷識(shí)別的準(zhǔn)確性。因此，本文采用一維離散小波降噪方法對(duì)三種信號(hào)進(jìn)行處理，小波降噪?yún)?shù)為：sqtwolog小波閥值，db6小波，6層分解[9-10]。降噪后的信號(hào)如圖4所示。對(duì)比這兩幅圖片可以發(fā)現(xiàn)降噪后的信號(hào)混疊較小，故降噪后的信號(hào)噪聲被明顯抑制。

分別將正常狀態(tài)，一級(jí)活塞斷環(huán)，二級(jí)活塞斷環(huán)三種狀態(tài)下降噪后的信號(hào)數(shù)據(jù)使用EMD分解。列出三種狀態(tài)下的前6個(gè)IMF分量函數(shù)圖，正常狀態(tài)如圖5所示。

由于EMD算法的特性，不可避免地會(huì)出現(xiàn)虛假的IMF或者與原始信號(hào)關(guān)系不大的IMF分量，如何剔除這些分量關(guān)系著故障識(shí)別和診斷的準(zhǔn)確性。由于每一個(gè)IMF分量都正交于原始信號(hào)，故可得每一個(gè)IMF分量與原始信號(hào)的相關(guān)性，這里用相關(guān)系數(shù)ρ來表示相關(guān)性，ρ越大表示ρ所對(duì)應(yīng)的IMF分量與原始信號(hào)相關(guān)性越大。由概率論相關(guān)系數(shù)計(jì)算公式可得ρ：

其中i表示第i個(gè)IMF分量。

表1列出了三種狀態(tài)下的前6個(gè)IMF分量與原函數(shù)的相關(guān)系數(shù)。由表1可知三種狀態(tài)下的前5個(gè)IMF分量與各自原信號(hào)相關(guān)度較高，第6個(gè)相關(guān)度較低。由此可知第6及其以后的IMF函數(shù)由于相關(guān)度很低，故可以舍棄不計(jì),于是得到5個(gè)敏感的IMF分量。

表1 前6個(gè)IMF分量的相關(guān)系數(shù)

對(duì)每種狀態(tài)各個(gè)IMF分量使用希爾伯特變換重構(gòu)原函數(shù)后得到其時(shí)間-頻率-幅值三者間關(guān)系的希爾伯特譜，其時(shí)頻三維圖如圖6所示。從圖6可得原函數(shù)頻率幅值較大的部分分布在0～5 000 Hz之間，故信號(hào)的能量也是大部分分布在這個(gè)頻率區(qū)間之中。

時(shí)頻圖所反映的頻率能量分布，當(dāng)空壓機(jī)活塞環(huán)斷裂時(shí)會(huì)激起船舶空壓機(jī)系統(tǒng)的固有頻率，此時(shí)形成共振，故幅值較大，能量較多，因此幅值較大的頻段為固有頻段(f)，且固有頻段是0～5 000 Hz。不同的故障激起因有頻段幅值不同，即能量不同，因此可以將空壓機(jī)固有頻段的能量值作為故障診斷的特征。由式(2)對(duì)希爾伯特譜進(jìn)行時(shí)間(點(diǎn)數(shù))積分便得到的希爾伯特邊際譜，它能準(zhǔn)確地反映信號(hào)幅值隨頻率變化的規(guī)律。

根據(jù)式(2)可以求出各個(gè)狀態(tài)每個(gè)采集樣本的0～5 000 Hz的希爾伯特邊際譜能量值，即得到各種狀態(tài)下固有頻段的能量值。

4 馬氏距離

馬氏距離是多元數(shù)理統(tǒng)計(jì)理論中常用的判別方法之一，用于表示數(shù)據(jù)的協(xié)方差距離。該方法的計(jì)算是建立在總體樣本的基礎(chǔ)上，排除了變量之間相關(guān)性的干擾，是一種衡量兩個(gè)未知樣本集相似程度的有效方法[6]。由于算法簡單，不受量綱的影響，計(jì)算速度快，并適用于小樣本問題的處理，因此，被廣泛應(yīng)用于故障診斷領(lǐng)域。

使用HHT算法將采集到的振動(dòng)信號(hào)進(jìn)行分解，并得到希爾伯特邊際譜，然后得到固有頻段的能量值,以固有頻段的能量值為特征值。選取每種狀態(tài)下每個(gè)加速度傳感器的5組共20組樣本進(jìn)行訓(xùn)練，其余20組進(jìn)行測(cè)試。計(jì)算出每種狀態(tài)下的20組樣本的能量值的均值和方差，記為標(biāo)準(zhǔn)特征值Si和標(biāo)準(zhǔn)方差var(si)。三種狀態(tài)下的標(biāo)準(zhǔn)特征值和標(biāo)準(zhǔn)方差如表2所示。

待檢信號(hào)的固有頻段能量特征值與三種狀態(tài)下的標(biāo)準(zhǔn)特征值之間的馬氏距離由式(3)計(jì)算得到。式(3)可以判斷出兩個(gè)樣本之間的相似性，馬氏距離越小，相似性越大。

(3)

其中：sx是待測(cè)樣本數(shù)據(jù)的能量特征值；si和var(si)為標(biāo)準(zhǔn)特征值和標(biāo)準(zhǔn)方差；di是求得的馬氏距離，d1表示正常狀態(tài)，d2表示一級(jí)活塞斷環(huán)，d3表示二級(jí)活塞斷環(huán)。

表2 標(biāo)準(zhǔn)特征值及其標(biāo)準(zhǔn)方差

判斷過程如下：取某一個(gè)測(cè)試樣本的固有頻段的能量值sx，然后分別代入三種狀態(tài)下所對(duì)應(yīng)的標(biāo)準(zhǔn)特征值和標(biāo)準(zhǔn)方差，求得三種狀態(tài)下所對(duì)應(yīng)的三個(gè)馬氏距離(d1，d2，d3)，比較這三個(gè)馬氏距離的大小，其中最小的數(shù)值所對(duì)應(yīng)的狀態(tài)為此樣本所對(duì)應(yīng)的狀態(tài)(即：當(dāng)d1最小時(shí)對(duì)應(yīng)正常狀態(tài)，當(dāng)d2最小時(shí)對(duì)應(yīng)的故障是一級(jí)活塞環(huán)斷裂；當(dāng)d3最小時(shí)對(duì)應(yīng)的故障是二級(jí)活塞環(huán)斷裂)，故而識(shí)別出該樣本所對(duì)應(yīng)的船用空壓機(jī)的狀態(tài)，達(dá)到故障診斷的目的。依照HHT邊際譜和馬氏距離判別方法的各個(gè)步驟，對(duì)三種不同狀態(tài)下的空壓機(jī)60組振動(dòng)信號(hào)進(jìn)行故障診斷研究，診斷結(jié)果表明：在正常狀態(tài)下全部診斷正確，一級(jí)活塞斷環(huán)有一處出現(xiàn)誤判，二級(jí)活塞斷環(huán)出現(xiàn)兩處錯(cuò)誤，故整體診斷正確率為95%。結(jié)果表明該方法可以進(jìn)行船舶二級(jí)往復(fù)式空壓機(jī)故障診斷，且診斷準(zhǔn)確率較高。出現(xiàn)誤判的原因可能是測(cè)量引起誤差。部分診斷結(jié)果如表3所示。

表3 部分馬氏距離診斷結(jié)果

續(xù)表

狀態(tài)Statesd1d2d3診斷結(jié)果Results一級(jí)活塞環(huán)斷裂First piston fracture7.358 11.535 57.259 8一級(jí)活塞環(huán)斷裂First piston fracture一級(jí)活塞環(huán)斷裂First piston fracture7.368 01.852 38.562 9一級(jí)活塞環(huán)斷裂First piston fracture一級(jí)活塞環(huán)斷裂First piston fracture8.996 41.526 17.333 0一級(jí)活塞環(huán)斷裂First piston fracture二級(jí)活塞環(huán)斷裂Second piston fracture32.120 3 10.231 71.365 8二級(jí)活塞環(huán)斷裂Second piston fracture二級(jí)活塞環(huán)斷裂Second piston fracture29.356 2 11.893 51.986 5二級(jí)活塞環(huán)斷裂Second piston fracture二級(jí)活塞環(huán)斷裂Second piston fracture33.235 09.200 12.425 0二級(jí)活塞環(huán)斷裂Second piston fracture二級(jí)活塞環(huán)斷裂Second piston fracture29.556 110.222 43.658 9二級(jí)活塞環(huán)斷裂Second piston fracture二級(jí)活塞環(huán)斷裂Second piston fracture34.207 08.362 41.520 3二級(jí)活塞環(huán)斷裂Second piston fracture

5 結(jié)論

1)對(duì)于船舶二級(jí)往復(fù)式空壓機(jī)振動(dòng)信號(hào)的非線性非平穩(wěn)性，采用EMD將其分解成若干個(gè)IMF分量，求出敏感IMF分量，進(jìn)而對(duì)每個(gè)敏感IMF分量進(jìn)行希爾伯特變換，再重構(gòu)信號(hào)，得出希爾伯特邊際譜，利用空壓機(jī)固有頻段的能量值作為信號(hào)的特征值，以此為判別馬氏距離的依據(jù)。實(shí)驗(yàn)表明，該方法識(shí)別正確率在95%，可以有效進(jìn)行故障識(shí)別。

2)一維離散小波降噪可以明顯降低環(huán)境背景噪聲干擾對(duì)振動(dòng)信號(hào)的影響，提高了故障診斷的準(zhǔn)確性。

3)在船舶往復(fù)式二級(jí)空壓機(jī)斷環(huán)故障診斷過程中，存在誤判是因?yàn)榇嬖跍y(cè)量誤差和系統(tǒng)誤差。本文的研究結(jié)果是基于小樣本實(shí)驗(yàn)，要建立完善的空壓機(jī)活塞環(huán)斷裂的HHT邊際譜的故障診斷還需要大數(shù)據(jù)進(jìn)一步分析。