包絡(luò)譜稀疏指標(biāo)盲濾波器的軸承狀態(tài)監(jiān)測(cè)

吳海威，王業(yè)統(tǒng)，李 美，黃福偉

（1.海南科技職業(yè)大學(xué)信息工程學(xué)院，海南 ?？?571126；2. 海南大學(xué)機(jī)電工程學(xué)院，海南 ?？?570228；3.海南大學(xué)計(jì)算機(jī)與網(wǎng)絡(luò)空間安全學(xué)院，海南 ?？?570228）

1 引言

旋轉(zhuǎn)機(jī)械的狀態(tài)監(jiān)測(cè)應(yīng)用已經(jīng)越來越廣泛，先進(jìn)的工業(yè)開始在其設(shè)計(jì)中引入越來越多的專用狀態(tài)監(jiān)測(cè)系統(tǒng)［1］。由于工業(yè)系統(tǒng)不斷提高的計(jì)算能力以及存儲(chǔ)能力，并且進(jìn)一步降低現(xiàn)代IT系統(tǒng)的成本，因此基本狀態(tài)監(jiān)測(cè)作為系統(tǒng)監(jiān)控的保護(hù)手段地位更為凸顯。尤其針對(duì)復(fù)雜的機(jī)械設(shè)備，例如風(fēng)力渦輪機(jī)和風(fēng)電機(jī)組的齒輪箱［2-3］。

由于多領(lǐng)域統(tǒng)計(jì)指標(biāo)的簡(jiǎn)單性，盲跟蹤信號(hào)特征的概念已經(jīng)在工業(yè)系統(tǒng)監(jiān)測(cè)中得到廣泛應(yīng)用，然而，這一盲信號(hào)分析概念仍有進(jìn)一步改進(jìn)的空間，并使其與其他最先進(jìn)的振動(dòng)分析方法更好地結(jié)合起來，最普遍的檢測(cè)方法之一是觀察振動(dòng)信號(hào)的循環(huán)平穩(wěn)行為［4］。盲濾波方法類似于最小熵卷積（MED）濾波所采用的盲反卷積思想，該濾波使用時(shí)間波形的峰度作為統(tǒng)計(jì)量來最大化MED，最初應(yīng)用于地震信號(hào)處理領(lǐng)域［5］。然而，峰度作為最大化的措施有一些缺點(diǎn)。當(dāng)脈沖頻率增加時(shí)，峰度值趨于減小。當(dāng)脈沖相互重疊時(shí)，峰度值甚至可以減小到高斯信號(hào)的峰度值；另一方面，當(dāng)脈沖頻率非常低，以至于測(cè)量中只有一個(gè)脈沖時(shí)，峰度值趨于最大［6］。

近年來，人們關(guān)注的焦點(diǎn)已經(jīng)從時(shí)間波形統(tǒng)計(jì)轉(zhuǎn)向循環(huán)平穩(wěn)特性。文獻(xiàn)［7］提出的方法在盲反褶積過程中采用了二階循環(huán)平穩(wěn)性指標(biāo)，它允許用戶以高效靈活的方式最大化ICS2。然而，該方法仍然需要故障頻率的先驗(yàn)知識(shí)才能工作。為了彌補(bǔ)這個(gè)潛在的問題，稀疏性的概念應(yīng)用于包絡(luò)譜，文獻(xiàn)［8］提出在濾波信號(hào)的包絡(luò)譜上使用范數(shù)，他們擴(kuò)展了源于峭度圖的濾波器組的思想，用矩代替峭度。然而，這種方法仍然不是完全盲的，即最終用戶仍然需要用范數(shù)值檢查濾波器組并選擇適當(dāng)?shù)念l帶。

為了解決上述方法存在的問題，提出了一種包絡(luò)譜稀疏指標(biāo)盲濾波器的齒輪狀態(tài)監(jiān)測(cè)方法。該方法能夠在不需要知道有關(guān)機(jī)械部件的特征故障頻率等先驗(yàn)信息條件下，利用包絡(luò)譜的稀疏性來跟蹤振動(dòng)信號(hào)中具有二階循環(huán)平穩(wěn)特征的故障。軸承數(shù)據(jù)集實(shí)驗(yàn)結(jié)果證明提出方法的有效性。

2 提出的方法

2.1 稀疏性度量

自從稀疏性大部分由于壓縮感知而成為一個(gè)熱門的研究方向以來，研究人員引入了許多稀疏性度量。每種方法都有其優(yōu)點(diǎn)和缺點(diǎn)，其中一個(gè)廣泛使用的稀疏性指標(biāo)是范數(shù)，屬于范數(shù)族［9］。

另一個(gè)平方包絡(luò)譜的指標(biāo)是譜域中的譜負(fù)熵。譜負(fù)熵源于將信號(hào)中的瞬時(shí)能量流轉(zhuǎn)換為概率分布，因此，負(fù)熵考慮了能量流x的歸一化。

譜負(fù)熵是譜熵的負(fù)，在時(shí)域中的定義為：

式中：∈x(n)—信號(hào)xk(n)的時(shí)域包絡(luò)。在光譜域中的包絡(luò)為：

時(shí)域中的負(fù)熵標(biāo)明了脈沖信號(hào)的內(nèi)容，很像峭度圖。另一方面，譜域中的負(fù)熵量化了重復(fù)或循環(huán)平穩(wěn)的信號(hào)內(nèi)容。

因此，這里研究了譜域中的負(fù)熵，并將其作為盲濾波的第二個(gè)度量指標(biāo)。

這里選擇將Hoyer指數(shù)作為包絡(luò)譜的第三個(gè)也是最后一個(gè)稀疏度度，Hoyer指數(shù)本質(zhì)上是范數(shù)的規(guī)范化版本［10］。當(dāng)且僅當(dāng)光譜僅包含單個(gè)非零分量時(shí)，將Hoyer指數(shù)的計(jì)算結(jié)果作為單位；當(dāng)且僅當(dāng)所有光譜分量相等時(shí)，Hoyer指數(shù)的值為0。

式中：x—樣本向量；n—樣本數(shù)；N—樣本總數(shù)。

2.2 盲濾波

盲濾波即找到一種濾波器，該濾波器最大化從噪聲測(cè)量信號(hào)x開始的信號(hào)：

式中：s—估計(jì)輸入；h—逆濾波器，*—卷積運(yùn)算。注意，向量和矩

陣以粗體字表示，以區(qū)別標(biāo)量，卷積運(yùn)算表示為：

式中：L、N—s和h的樣本數(shù)量。因此平方包絡(luò)∈x的計(jì)算式如下：

上式還可以表示為：

式中：sH—s的Hermitian轉(zhuǎn)置；diag(sH) —一個(gè)對(duì)角矩陣，其對(duì)角線上的值為向量sH的值。平方包絡(luò)譜Ex—該平方包絡(luò)的傅里葉變換：

其中，

其中，n= 0..L-N- 1，k= 0..K- 1 和K- 1 對(duì)應(yīng)于最大感興趣頻率的指數(shù)，F(xiàn)的尺寸( )L-N，K。

平方包絡(luò)譜Ex的范數(shù)的計(jì)算公式如下：

用向量表示上式，公式如下：

分母可用以下恒等式計(jì)算：

式（13）的分子可進(jìn)一步推導(dǎo)，如下所示：

因此，式（13）可以改寫為：

將式（10）代入式（16）中，可以得到：

可用式（17）得到廣義Rayleigh熵，并可通過迭代特征值將其最大化：

Rayleigh熵有一個(gè)性質(zhì)：它相對(duì)于h的最大值等于它的最大特征值λ和相應(yīng)的特征向量［11］。

因此，最大化Rayleigh 熵可以找到相應(yīng)指標(biāo)和濾波器的最大值。

然而，為了得到實(shí)特征值，相關(guān)矩陣RXW1Xh和RXW2Xh需要符合Hermitian矩陣，并且RXW2Xh需要滿足半正定。

如果滿足這些條件，Rayleigh熵是迭代計(jì)算濾波器系數(shù)的有效方法。

每次迭代只需計(jì)算最大特征值和對(duì)應(yīng)的特征向量，利用冪函數(shù)法等算法可以有效地實(shí)現(xiàn)這一點(diǎn)。因此可改寫為：

其中，

（1）初始化h；

（2）用式（21）基于h和X計(jì)算RXW1Xh和RXW2Xh；

（3）通過式（22）找出λmax和與所用標(biāo)準(zhǔn)的較高值相對(duì)應(yīng)的新濾波器h；

（4）使用新的h返回第二步，直到達(dá)到收斂或者迭代完成。

用向量表示上式：

經(jīng)過多次重新改寫后，廣義Rayleigh熵可表示為：

譜負(fù)熵ΔIE為譜熵的負(fù)值，對(duì)于包絡(luò)譜的定義為：

用向量表示上式：

其廣義Rayleigh熵為：

譜負(fù)熵也具有一定的作用，即其始終大于零且小于ln(N)：

這對(duì)時(shí)間域ΔI∈中的負(fù)熵同樣有效。

證明：由于譜負(fù)熵始終大于零且小于ln(N) ：

根據(jù)香農(nóng)公式信息熵H可得到下式：

式中：pi—一組獨(dú)立的概率。因此，對(duì)于，可證明下式：

由Kullback-Leibler散度（或相對(duì)熵）可知：

因此，當(dāng)，有：

以上證明過程同樣適用于時(shí)域ΔI∈中的負(fù)熵。

2.3 實(shí)際考慮因素

包絡(luò)譜稀疏性的最大化是基于增加稀疏性的主要成分是潛在故障特征的假設(shè)。因此，首先需要對(duì)信號(hào)進(jìn)行白化，以便在解決最大化問題之前去除確定性因素，因?yàn)檫@些因素通常會(huì)在包絡(luò)頻譜中引入高幅度離散峰值。

相關(guān)矩陣RXW1X和RXW2X的計(jì)算包括使用傅立葉矩陣來計(jì)算包絡(luò)譜。然而，在實(shí)際計(jì)算中，可以用FFT代替，以減少計(jì)算時(shí)間和內(nèi)存需求。

計(jì)算的總持續(xù)時(shí)間將主要取決于信號(hào)長度和選擇的濾波器長度。由于長信號(hào)和濾波器的原因，指標(biāo)可能變得相當(dāng)大，因此需要的乘法運(yùn)算和FFT的點(diǎn)數(shù)可能會(huì)導(dǎo)致計(jì)算時(shí)間較長。

2.4 非平穩(wěn)角速度情況的推廣

該方法的改進(jìn)可以使濾波器相關(guān)矩陣的傅里葉矩陣的核依賴于角度而不是依賴時(shí)間。這樣，瞬時(shí)角速度變化可以直接考慮在濾波器中，可以通過速度同步離散傅里葉變換（VSDFT）輕松實(shí)現(xiàn)。然后傅里葉矩陣可以改寫為：

再將其傳遞至濾波器優(yōu)化。

2.5 單輸入多輸出（SIMO）系統(tǒng)的延伸

這里方法允許一個(gè)相當(dāng)直接的擴(kuò)展到一個(gè)有Q響應(yīng)xq的情況。由hq濾波的xq的每個(gè)響應(yīng)可以疊加在一起，然后返回s的估計(jì)。

使用Rayleigh商的迭代過程仍然可以用于SIMO情況，唯一需要的調(diào)整就是相關(guān)矩陣和濾波器向量。將RXW1X和RXW2X表示為互相關(guān)矩陣：

對(duì)于RqWq，xq的加權(quán)自相關(guān)矩陣和非對(duì)角化矩陣是xq和xq?的加權(quán)互相關(guān)矩陣。然后，過濾器h轉(zhuǎn)化為：

然后，興趣信號(hào)s的第q個(gè)貢獻(xiàn)的計(jì)算公式如下：

Xq是式（7）中定義的xq的Toeplitz 矩陣。因此，整個(gè)興趣信號(hào)s的計(jì)算公式如下：

3 仿真結(jié)果分析

為了驗(yàn)證所提出的方法，首先考慮了一個(gè)仿真案例：具有高斯分布振幅的周期脈沖。為了增加一些對(duì)比參考，將該方法與基于最小熵反褶積（MED）的時(shí)域盲反褶積濾波（ICS2）和基于循環(huán)平穩(wěn)源信號(hào)的時(shí)域譜負(fù)熵濾波的性能進(jìn)行了比較［12-13］。

周期脈沖信號(hào)s1與IRFg1，s和高斯噪聲g1，n進(jìn)行卷積。生成的信號(hào)，如圖1所示。但是，不是只分析單個(gè)仿真，而是將加性高斯噪聲的信噪比從-80DB 變化到0DB，以評(píng)估濾波器的趨勢(shì)行為。整個(gè)信號(hào)的總方差保持不變。對(duì)于時(shí)域盲濾波器，使用40個(gè)樣本的濾波器長度，而基于稀疏性的譜濾波器使用20個(gè)樣本的濾波器長度，最大迭代次數(shù)設(shè)置為50。

圖1 模擬情況1輸入信號(hào)Fig.1 Simulate the Input Signal of Case 1

從-80dB到0dB SNR的結(jié)果指示器變化曲線，如圖2所示?？梢钥闯觯谶@種特殊情況下，ICS2濾波器的性能優(yōu)于MED和時(shí)域負(fù)熵濾波，基于稀疏度指標(biāo)的濾波器顯示在大約-38dB處顯著增加，類似于ICS2濾波器。

圖2 仿真案例中不同指標(biāo)的趨勢(shì)分析Fig.2 Trend Analysis of Different Indicators in Simulation Case

為了進(jìn)一步驗(yàn)證，平方包絡(luò)光譜的顏色圖，如圖3所示。時(shí)域?yàn)V波器（ICS2、MED和時(shí)域負(fù)熵）的包絡(luò)特性表明，對(duì)于更高的信噪比，故障頻率在4Hz處具有明顯的諧波。相比之下，稀疏濾波器的包絡(luò)譜在4Hz處顯示出一次故障諧波，但同時(shí)也顯示出一個(gè)較低的頻率分量，其振幅較高。這是由于模擬信號(hào)的包絡(luò)譜包含非常低的頻率調(diào)制，這是由于周期脈沖的高斯分布振幅造成的。此外，這種低頻調(diào)制增加了包絡(luò)頻譜的稀疏性。

圖3 模擬情況下歸一化平方包絡(luò)譜的趨勢(shì)分析Fig.3 Trend Analysis of Normalized Square Envelope Spectrum in Simulation

濾波信號(hào)，如圖4所示。信噪比為-19db。信號(hào)的振幅在-1和+1之間進(jìn)行標(biāo)準(zhǔn)化。結(jié)果表明，時(shí)域?yàn)V波器對(duì)輸入信號(hào)的周期脈沖序列進(jìn)行去卷積，而基于包絡(luò)譜稀疏性的濾波器提取周期脈沖序列用IRF反褶積后的周期脈沖序列。然而，脈沖序列產(chǎn)生的包絡(luò)譜顯示出多個(gè)諧波，與單個(gè)離散峰值相比，減少了稀疏性。因此基于稀疏性的方法在與IRF卷積后過濾掉脈沖序列，這導(dǎo)致稀疏性介于純正弦調(diào)制和脈沖序列之間，因此，在通過跟蹤故障進(jìn)行故障檢測(cè)的意義上，不必去除與IRF 的卷積包絡(luò)譜中的峰值。在有無IRF 卷積的情況下，包絡(luò)譜的特征將有所不同。此外，由于假設(shè)主要的二階循環(huán)平穩(wěn)分量來自故障源，因此無論IRF卷積如何，都可以通過包絡(luò)譜稀疏性來跟蹤故障的演化。

圖4 仿真情況1的-19dB信噪比下的時(shí)間波形Fig.4 Time Waveform Under -19dB SNR in Simulation Case 1

4 實(shí)驗(yàn)分析和結(jié)果

這里方法的數(shù)據(jù)集是Cincinnati 大學(xué)智能維修系統(tǒng)中心（IMS）提供的來自美國宇航局預(yù)測(cè)數(shù)據(jù)倉庫的振動(dòng)信號(hào)，該數(shù)據(jù)集是在一個(gè)軸承發(fā)生外圈故障的運(yùn)行到故障試驗(yàn)期間測(cè)量的。軸承試驗(yàn)臺(tái)，如圖5所示。由一根軸連接到一臺(tái)轉(zhuǎn)速為2000RPM（轉(zhuǎn)/分）的交流電機(jī)上。施加2700kg的徑向載荷。四個(gè)Rexnord ZA-2115軸承與高端加速計(jì)一起安裝在軸上。總共進(jìn)行了三次故障檢修的實(shí)驗(yàn)，采樣速率為20kHz，每10min記錄一次。這里描述的方法在第二個(gè)數(shù)據(jù)集（包含984個(gè)測(cè)量樣本）上進(jìn)行了測(cè)試，其中軸承1由于外圈缺陷而發(fā)生故障。

圖5 軸承試驗(yàn)設(shè)備Fig.5 Bearing Test Equipment

然后比較了模擬使用的六種濾波器，為比較各種方法的性能提供一些參考。外圈的特征故障頻率約為236Hz，因此用作ICS2濾波器的輸入。

過濾器長度設(shè)為20 個(gè)樣本，最大迭代次數(shù)設(shè)為50。理論上，ICS2濾波器直接跟蹤故障頻率調(diào)制，因此它跟蹤故障發(fā)生時(shí)的信號(hào)特征變化的性能最好。

結(jié)果似乎證實(shí)了上述結(jié)論，如圖6所示。在第4天開始時(shí)，ICS2過濾器的ICS2指標(biāo)開始上升。MED濾波器和時(shí)域負(fù)熵濾波器的性能很差，因?yàn)槲墨I(xiàn)［14］已經(jīng)證明峰度或其他跟蹤非高斯性的指標(biāo)在這個(gè)特定的數(shù)據(jù)集上表現(xiàn)不佳。事實(shí)上，從測(cè)試開始信號(hào)就已經(jīng)是非高斯信號(hào)了。因此，在這種情況下，這兩個(gè)濾波器在有效跟蹤故障方面表現(xiàn)不佳。

圖6 在IMS數(shù)據(jù)集上用六種不同的方法濾波后變化Fig.6 Changes After Filtering with Six Different Methods on IMS Data Set

另一方面，基于包絡(luò)譜稀疏性的濾波器表明，由于故障，信號(hào)的主要循環(huán)平穩(wěn)性增加。在第5天開始時(shí)，所有三個(gè)過濾器的指示值都會(huì)突然增加。軸承劣化主要表現(xiàn)為循環(huán)平穩(wěn)性的顯著波動(dòng)，并不是非高斯性的波動(dòng)。濾波從噪聲中剔除包絡(luò)譜中的故障頻率，濾波前后包絡(luò)譜的差異，如圖7所示。從而簡(jiǎn)化并加速了對(duì)故障類型的信號(hào)的檢查。

圖7 測(cè)量用Hoyer指數(shù)濾波器濾波后的包絡(luò)譜Fig.7 Envelope Spectrum After Filtering by Hoyer Exponential Filter

5 結(jié)果分析

實(shí)驗(yàn)結(jié)果表明，在不預(yù)先知道故障頻率的情況下，利用振動(dòng)信號(hào)的盲濾波方法進(jìn)行狀態(tài)監(jiān)測(cè)是可行的。仿真結(jié)果表明，基于包絡(luò)譜稀疏性的濾波方法優(yōu)于傳統(tǒng)的基于時(shí)間指示器的濾波方法，如MED。尤其是當(dāng)故障特征在信號(hào)中引入循環(huán)平穩(wěn)性時(shí)，這一點(diǎn)尤為明顯，而包絡(luò)譜可以有效地檢測(cè)到信號(hào)。針對(duì)故障脈沖重復(fù)率高且在下一個(gè)脈沖前沒有完全衰減的情況，說明基于峭度濾波具有一定的局限性。在這種情況下，由于缺乏強(qiáng)烈的非高斯性質(zhì)，峰度只會(huì)受到輕微影響，當(dāng)軸承頻率在千赫茲范圍內(nèi)時(shí)，可能會(huì)發(fā)生這一現(xiàn)象。

與傳統(tǒng)的盲反褶積濾波器（如MED）相比，基于稀疏度的濾波器在工作方式上的另一個(gè)顯著區(qū)別是，這里設(shè)計(jì)的濾波器純粹是為了最大化信號(hào)的某一特性。這意味著它們不會(huì)恢復(fù)為源信號(hào)（例如，重復(fù)沖擊）、從混合物中分離信號(hào)或?qū)π盘?hào)進(jìn)行去噪。這與過去的方法不同，過去通常使用盲過濾。因此，用戶在分析結(jié)果時(shí)需要記住這個(gè)屬性，因?yàn)檫^濾后的響應(yīng)可能不能準(zhǔn)確地表示感興趣的信號(hào)，而是最大化所用標(biāo)準(zhǔn)的失真版本。然而，由于包絡(luò)譜是基于振動(dòng)的狀態(tài)監(jiān)測(cè)中最常用的工具之一，因此該方法的泛化性能較高。

實(shí)驗(yàn)結(jié)果表明了這里方法對(duì)循環(huán)平穩(wěn)源的靈敏度較高。雖然稀疏濾波器在檢測(cè)故障時(shí)比ICS2濾波器稍慢，但它們?nèi)匀荒軌蛟诓恢缹?shí)際軸承故障頻率的情況下準(zhǔn)確地跟蹤故障。結(jié)果還表明，簡(jiǎn)單地觀察濾波信號(hào)的包絡(luò)譜稀疏性，可以作為信號(hào)循環(huán)平穩(wěn)性的可靠盲跟蹤度量指標(biāo)。最后，這里的三種稀疏性措施在跟蹤能力方面的表現(xiàn)并不是很突出。

這里采用了一種新的方法來觀察盲濾波和包絡(luò)譜的結(jié)合。然而，幾乎每種方法都存在一些缺點(diǎn)，可能會(huì)對(duì)結(jié)果產(chǎn)生負(fù)面影響。使用這種過濾器的主要難題在于過濾器長度的選擇。這個(gè)長度會(huì)對(duì)結(jié)果和計(jì)算時(shí)間產(chǎn)生較大影響。顯然，較長的濾波器長度會(huì)增加計(jì)算時(shí)間，且也會(huì)導(dǎo)致濾波器系數(shù)的收斂速度變慢或不收斂。如果由于這種不收斂而導(dǎo)致指示值劇烈振蕩，則后者可能特別麻煩。因此，建議縮短濾波器的長度，同時(shí)抑制對(duì)濾波器無用的頻率。

6 結(jié)論

為了解決包絡(luò)譜稀疏問題，提出了一種基于包絡(luò)譜稀疏指標(biāo)盲濾波器的軸承狀態(tài)監(jiān)測(cè)方法。提出了一種盲濾波技術(shù)，利用包絡(luò)譜的稀疏性來跟蹤振動(dòng)信號(hào)中具有二階循環(huán)平穩(wěn)特征的故障的存在。三個(gè)稀疏性度量用于推導(dǎo)適用于基于Rayleigh商迭代的方法的濾波器公式。仿真和實(shí)驗(yàn)結(jié)果證明了如下結(jié)論：

（1）這里方法對(duì)循環(huán)平穩(wěn)源的靈敏度較高，雖然稀疏濾波器在檢測(cè)故障時(shí)比ICS2濾波器稍慢，但它們?nèi)匀荒軌蛟诓恢缹?shí)際軸承故障頻率的情況下準(zhǔn)確地跟蹤故障。

（2）當(dāng)故障信號(hào)不具有強(qiáng)脈沖成分時(shí)，該方法的性能優(yōu)于傳統(tǒng)的盲反卷積方法，提出方法能夠有效解決包絡(luò)譜稀疏問題，實(shí)現(xiàn)有效的濾波器狀態(tài)監(jiān)測(cè)。

（3）較長的濾波器長度會(huì)增加計(jì)算時(shí)間，且也會(huì)導(dǎo)致濾波器系數(shù)的收斂速度變慢或不收斂。