基于離散小波變換和Kalman濾波的直升機(jī)主減智能狀態(tài)預(yù)測

劉立生 ，楊宇航

(1.南京理工大學(xué) 自動化學(xué)院，南京 210094;2.總參謀部陸航研究所 可靠性室，北京 101121)

直升機(jī)主減具有傳動效率高、結(jié)構(gòu)緊湊等特點(diǎn)，起減速、轉(zhuǎn)向及并車的作用。它將一臺或多臺高轉(zhuǎn)速小扭矩的發(fā)動機(jī)功率變成低轉(zhuǎn)速、大扭矩的功率傳遞給旋翼軸，并按轉(zhuǎn)速和扭矩的需要將功率傳遞給旋翼、尾槳及各個(gè)附件，以保證直升機(jī)正常工作。作為旋轉(zhuǎn)機(jī)械設(shè)備中的一種，主減常見故障類型有機(jī)器運(yùn)行失穩(wěn)、發(fā)生異常振動和噪聲、轉(zhuǎn)速發(fā)生變化等。Land［1］對HUMS(Healthy and Usage Monitoring Systems)進(jìn)行研究;Samuel等［2］提出CBM+項(xiàng)目。這些項(xiàng)目的共同點(diǎn)是都對主減的狀態(tài)進(jìn)行監(jiān)控，從而減少了使用與保障(O＆S)費(fèi)用，提高了直升機(jī)的安全性和維修效率。

振動信號常被用來對機(jī)械的故障進(jìn)行診斷［3-4］，利用小波包變換［5-6］或者離散小波變換［7-10］對振動信號進(jìn)行分解是兩種常用方法。在進(jìn)行小波變換時(shí)，小波的選擇是關(guān)鍵，選擇不同的小波可能會得到不同的結(jié)果，Rafiee等［11］比較了324 種母小波，發(fā)現(xiàn)“db44”與振動信號最相似，分解效果最好。大量的研究是針對故障檢測進(jìn)行的，而對故障進(jìn)行預(yù)測的研究則相對比較少，文獻(xiàn)［12］對預(yù)測的方法進(jìn)行了總結(jié)，包括五種方法:基于物理模型、人工神經(jīng)網(wǎng)絡(luò)、先驗(yàn)知識(專家系統(tǒng)、模糊系統(tǒng))、統(tǒng)計(jì)方法(趨勢推斷、ARMA模型)以及隨機(jī)情況(Kalman濾波、Markov模型、Bayesian網(wǎng)絡(luò)、可靠性函數(shù))?，F(xiàn)實(shí)生活中設(shè)備受外界的干擾是隨機(jī)的，故對隨機(jī)情況進(jìn)行研究更具有意義。Kalman濾波［13-17］通過一種算法排除可能的隨機(jī)干擾，提高檢測精度，能夠很好地對數(shù)據(jù)進(jìn)行預(yù)測，廣泛應(yīng)用于機(jī)器人導(dǎo)航、控制、傳感器數(shù)據(jù)融合以及軍事方面的雷達(dá)系統(tǒng)以及導(dǎo)彈追蹤等。

本文研究了離散小波變換、Kalman濾波以及Elman神經(jīng)網(wǎng)絡(luò)組成的智能狀態(tài)預(yù)測系統(tǒng)。文中首先闡述了離散小波變換、Parseval定理、Kalman濾波以及Elman神經(jīng)網(wǎng)絡(luò)等理論知識，在文獻(xiàn)［13］Kalman預(yù)測算法的基礎(chǔ)上提出了一種新的Kalman濾波預(yù)測算法，接著給出了主減智能狀態(tài)預(yù)測系統(tǒng)的預(yù)測流程圖，然后以某型直升機(jī)主減上8路振動傳感器采集的振動數(shù)據(jù)為例，去噪、分解，利用Parseval定理提取各層能量，分別使用這兩種Kalman濾波算法對各路傳感器各個(gè)時(shí)刻的特征向量進(jìn)行預(yù)測，并用Elman神經(jīng)網(wǎng)絡(luò)對預(yù)測值進(jìn)行診斷，最終得出:本文提出的算法組成的智能狀態(tài)預(yù)測系統(tǒng)能夠更好地對故障進(jìn)行預(yù)測。

1 基于離散小波變換的特征提取

函數(shù)ψ(t)∈L2(R)，其 Fourier變換(ω)滿足允許條件，則稱ψ(t)為一個(gè)基本小波或者母小波(mother wavelet)，將母小波經(jīng)過縮放和平移之后，就可以得到小波序列(a，b∈R，a≠0)。給定基本小波函數(shù)ψ(t)，信號f(t)的連續(xù)小波變換為:

其中:a為尺度參數(shù);b為平移參數(shù)。

對尺度參數(shù)a和平移參數(shù)b進(jìn)行離散化處理，則離散小波變換為:

式中:a=2-j，b=n2-j。

離散小波變換可以實(shí)現(xiàn)多分辨率分析，它將信號分解為近似部分cAj小波系數(shù)和細(xì)節(jié)部分cDj小波系數(shù)，接著對近似部分進(jìn)一步的分解，如此進(jìn)行迭代，其分解如圖1所示，其J層分解表示為:

圖1 離散小波變換Fig.1 Discrete wavelet transform

離散小波變換各層能量的提取采用Parseval定理［18］:

式中:x(t)是時(shí)域信號，x(f)是離散信號的Fourier變換，N是采樣周期。把式(4)代入式(3)得:

2 Kalman濾波

Kalman濾波是以最小均方誤差估計(jì)為最佳準(zhǔn)則來尋求一套遞推估計(jì)的算法，其基本思想是:采用信號與噪聲的狀態(tài)空間模型，利用前一時(shí)刻的估計(jì)值和現(xiàn)時(shí)刻的觀測值，用狀態(tài)方程和遞推方法來估計(jì)非平穩(wěn)隨機(jī)信號的波形，求出現(xiàn)時(shí)刻的估計(jì)值。

考慮到直升機(jī)主減振動信號的特點(diǎn)，假設(shè)其狀態(tài)模型:

觀測模型:

式中:X(k)=［x(k)(k)］T，x(k)是k時(shí)刻的值(k)是x(k)在k時(shí)刻的變化率，Γ=［T2/2T］T，T為觀測周期，Z(k)是k時(shí)刻的測量值，H=［10］。w(k)和v(k)分別表示過程和測量的噪聲，假設(shè)為高斯白噪聲，協(xié)方差分別為Q，R。

定義:

(k|k-1):由k-1時(shí)刻的值估計(jì)(k);

P(k| k-1)k|k-1)對應(yīng)的估計(jì)誤差協(xié)方差;

(k|k):k時(shí)刻的最優(yōu)化估算值;

P(k| k)(k|k)對應(yīng)的估計(jì)誤差協(xié)方差;

對目標(biāo)進(jìn)行Kalman預(yù)測具體的流程如圖2，實(shí)際中通常無法得到目標(biāo)的初始值，常利用其前兩個(gè)觀測值建立初始估計(jì)，即:

進(jìn)而得到初始估計(jì)的估計(jì)誤差:

從而得到初始估計(jì)誤差協(xié)方差:

圖2 參數(shù)預(yù)測流程Fig.2 The flowchart of parameter estimation

在獲得變化率后，預(yù)測值使用得到，預(yù)測值有兩種算法:

(1)由初始兩個(gè)觀測值Z(1)，Z(2)，利用Kalman濾波求取不同時(shí)刻的變化率，獲得各時(shí)刻預(yù)測值［13］。

(2)由初始兩個(gè)觀測值Z(1)，Z(2)預(yù)測下一時(shí)刻的值(3)，用Z(2)與(3)構(gòu)成新的觀測值，預(yù)測下一時(shí)刻值(4)，再用(3)與(4)構(gòu)成新的觀測值，對下一時(shí)刻進(jìn)行預(yù)測，如此反復(fù)迭代，從而獲得所有時(shí)刻的預(yù)測值，該算法是本文提出的方法。

為了簡潔起見，用“A”表示文獻(xiàn)［13］的算法，“B”表示本文的算法。

3 Elman神經(jīng)網(wǎng)絡(luò)

Elman神經(jīng)網(wǎng)絡(luò)由Elman提出，該模型在前饋網(wǎng)絡(luò)的隱含層中增加一個(gè)承接層，進(jìn)行一步延時(shí)，達(dá)到記憶的目的，從而使系統(tǒng)具有適應(yīng)時(shí)變特性的能力，能直接反映系統(tǒng)的動態(tài)特性。Elman型回歸神經(jīng)網(wǎng)絡(luò)的結(jié)構(gòu)如圖2，它一般由輸入層、隱含層、承接層和輸出層組成，，其中承接層也稱狀態(tài)層或上下文層，隱含層的輸出通過承接層的延遲與存儲，自聯(lián)到隱含層的輸入。這樣就使其對歷史數(shù)據(jù)具有敏感性，增加了網(wǎng)絡(luò)自身處理動態(tài)信息的能力，從而達(dá)到了動態(tài)建模的目的。

圖3 Elman神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)Fig.3 Structure of Elman neural network

4 主減智能狀態(tài)預(yù)測流程圖

基于上面介紹的理論，用 DWT、Kalman濾波和Elman神經(jīng)網(wǎng)絡(luò)相結(jié)合構(gòu)成直升機(jī)主減的智能狀態(tài)預(yù)測系統(tǒng)，其流程圖如圖4。

其具體過程為:使用振動傳感器采集直升機(jī)主減上的振動信號，去噪后進(jìn)行離散小波變換，根據(jù)Parseval定理進(jìn)行特征提取，提取的特征向量進(jìn)行Kalman濾波由t時(shí)刻值預(yù)測t+τ時(shí)刻的值Xt+τ，最后經(jīng)Elman神經(jīng)網(wǎng)絡(luò)進(jìn)行診斷。Elman神經(jīng)網(wǎng)絡(luò)的診斷流程為:訓(xùn)練數(shù)據(jù)是每一時(shí)刻由實(shí)測值獲得的特征向量;測試數(shù)據(jù)是Kalman濾波對應(yīng)時(shí)刻的預(yù)測值，預(yù)測數(shù)據(jù)輸入訓(xùn)練好的神經(jīng)網(wǎng)絡(luò)進(jìn)行故障診斷，從而實(shí)現(xiàn)對直升機(jī)主減的狀態(tài)進(jìn)行智能預(yù)測。

圖4 主減智能狀態(tài)預(yù)測算法流程圖Fig.4 The workflow of MGB intelligent condition prediction algorithm

5 試驗(yàn)驗(yàn)證

振動傳感器監(jiān)測直升機(jī)整個(gè)飛行過程中主減上各傳動齒輪及軸承處的振動水平，該振動傳感器是適用于直升飛機(jī)飛行試驗(yàn)的ICP型傳感器，其工作溫度范圍寬，抗干擾能力強(qiáng)，數(shù)據(jù)采集設(shè)備采集8路振動信號，具體檢測部位和傳感器安裝位置如表1所示。

采集直升機(jī)地面試車其扭矩為10%時(shí)的振動數(shù)據(jù)，其初始采樣頻率為20 kHz，由于篇幅關(guān)系，以L001為例說明其具體過程。圖5為采集的L001信號，橫坐標(biāo)為具體的飛行時(shí)刻，縱坐標(biāo)為振動加速度，共有36 001 264個(gè)數(shù)據(jù)，它包括了一次飛行所有狀態(tài)下的振動數(shù)據(jù)，其包含的狀態(tài)太多，如:爬升、平飛、左側(cè)飛、俯沖、傾斜拉起等縱多飛行狀態(tài)。

表1 振動監(jiān)測部位Tab.1 Location for vibration monitoring and measure

圖5 L001采集的原始信號Fig.5 Sample signals of L001

由于數(shù)據(jù)采集設(shè)備配套的FAMOS軟件實(shí)現(xiàn)不了像小波分析、神經(jīng)網(wǎng)絡(luò)、Kalman濾波等復(fù)雜算法，但它可以對信號進(jìn)行重新采樣，而Matlab軟件一次性處理不了采樣頻率為20 kHz的大量數(shù)據(jù)，故本文利用FAMOS軟件對各路振動信號進(jìn)行重新采樣，然后利用Matlab軟件對新采樣后的信號進(jìn)行分析。通過對各路信號進(jìn)行重新采樣，使得采樣頻率為原頻率的1/20，即1 kHz，然后用語句:data1=Cut(L001，t1，t2)截取時(shí)間t1～t2之間的飛行數(shù)據(jù)，把這些數(shù)據(jù)拷貝到“.TXT”文檔中，供Matlab軟件調(diào)用。

選擇啟發(fā)式閥值對截取的信號進(jìn)行去噪處理，對去噪后的信號使用Rafiee［11］中的“db44”進(jìn)行離散小波分解，分解層數(shù)為9，每個(gè)傳感器得到10個(gè)特征向量，由于這10個(gè)能量數(shù)量級不同，如果直接輸入到Kalman濾波和Elman神經(jīng)網(wǎng)絡(luò)中，容易產(chǎn)生誤判，故需用如下公式對其進(jìn)行歸一化。

對Kalman濾波，分別采用A，B兩種算法對8路傳感器每個(gè)時(shí)刻的各個(gè)特征向量進(jìn)行預(yù)測，預(yù)測的總時(shí)間為30 s，各個(gè)時(shí)刻的誤差采用如下公式:

其中，xi，j為第j個(gè)傳感器的第i個(gè)測量值;為第j個(gè)傳感器的第i個(gè)預(yù)測值;N為傳感器個(gè)數(shù)。

采用A算法時(shí)，各時(shí)刻誤差組成的誤差曲線如圖6，橫坐標(biāo)為預(yù)測時(shí)間，縱坐標(biāo)為誤差值，可以看出采用無迭代方法A時(shí)，預(yù)測時(shí)間越長，其對應(yīng)的誤差就越大，而且誤差的變化率也在逐漸變大;采用B算法時(shí)，各個(gè)時(shí)刻對應(yīng)的誤差曲線如圖7，可以發(fā)現(xiàn)隨著預(yù)測間隔的變大，誤差曲線上下波動，而不是單調(diào)上升。比較圖7和圖6，可以清楚地發(fā)現(xiàn)，相同時(shí)刻算法B的誤差值比算法A小得多，部分時(shí)刻的具體值如表2，通過表2和對應(yīng)的圖可以知道，雖然圖7中波動比較大，但是它的縱坐標(biāo)很小，從而說明算法B更適合于對各個(gè)時(shí)刻進(jìn)行預(yù)測，它與實(shí)際測量值比較接近。

考慮主減常見的五種故障狀態(tài):健康、運(yùn)動失穩(wěn)、異常振動、噪聲、轉(zhuǎn)速變化，分別用(1 0 0 0 0);(0 1 0 0 0);(0 0 1 0 0);(0 0 0 1 0);(0 0 0 0 1)表示。

為了研究Elman神經(jīng)網(wǎng)絡(luò)各個(gè)時(shí)刻的逼近精度，選取健康狀態(tài)下的數(shù)據(jù)進(jìn)行訓(xùn)練和測試，未經(jīng)過Kalman濾波的各個(gè)時(shí)刻的特征數(shù)據(jù)用于訓(xùn)練神經(jīng)網(wǎng)絡(luò)，經(jīng)過Kalman濾波后獲得的預(yù)測值用來對網(wǎng)絡(luò)進(jìn)行測試，使用newelm構(gòu)建兩層Elman神經(jīng)網(wǎng)絡(luò)，隱含層函數(shù)為tansig傳遞函數(shù)，輸出層為purelin傳遞函數(shù)，算法A與算法B部分時(shí)刻的逼近情況如表3所示，這里的逼近精度是相對于各個(gè)具體時(shí)刻的實(shí)際測量值而言的。

表2 不同預(yù)測時(shí)刻的誤差值Tab.2 Error value at different prediction horizons

表3 不同預(yù)測時(shí)刻神經(jīng)網(wǎng)絡(luò)的逼近情況Tab.3 The approach accuracies of neural network at different prediction horizons

從上表可以看出，隨著預(yù)測時(shí)間的變大，算法A逼近精度在不斷地下降，如L008傳感器中，逼近精度在1 s時(shí)為 0.977 9，5 s時(shí)為 0.884 0，15 s時(shí)為0.690 7，28 s時(shí)為0.437 3;而 B 算法逼近精度相對都比較高，并且我們發(fā)現(xiàn)，隨著時(shí)間的推移，其相對逼近精度穩(wěn)定在了一個(gè)值上。這進(jìn)一步說明本文提出的算法更適合用于直升機(jī)主減智能狀態(tài)預(yù)測系統(tǒng)的構(gòu)建。

6 結(jié)論

直升機(jī)主減是傳動裝置中結(jié)構(gòu)最復(fù)雜、尺寸最大、重量最重的一個(gè)部件，它運(yùn)轉(zhuǎn)的狀況直接影響著直升機(jī)的安全性和可靠性。針對其故障特點(diǎn)，文中提出了一種基于離散小波變換、Kalman濾波以及Elman神經(jīng)網(wǎng)絡(luò)組成的直升機(jī)主減智能狀態(tài)預(yù)測系統(tǒng)，通過理論分析和實(shí)驗(yàn)驗(yàn)證得到如下結(jié)論:

(1)離散小波變換中，母小波選擇“db44”能很好地對信號進(jìn)行分解，有利于應(yīng)用Parseval定理對其各層進(jìn)行特征提取。

(2)本文的 Kalman濾波預(yù)測算法更適用于對主減的特征向量進(jìn)行預(yù)測，它能夠?qū)?路信號的80個(gè)特征向量進(jìn)行預(yù)測，并將各路特征向量的逼近精度穩(wěn)定在一個(gè)比較高的固定值上。

(3)應(yīng)用DWT、Kalman濾波以及Elman神經(jīng)網(wǎng)絡(luò)構(gòu)建直升機(jī)主減的智能狀態(tài)預(yù)測系統(tǒng)是可行的，有效的，它能對主減的狀態(tài)進(jìn)行精確的預(yù)測，將為未來HUMS、PHM、CBM和CBM+系統(tǒng)的進(jìn)一步開發(fā)提供新的技術(shù)參考。

［1］ Land J E.HUMS-the benefits-past，present and future［J］.IEEE，2001，6(s):3083-3094.

［2］ Samuel P D，Pines D J.A review of vibration-based techniques for helicopter transmission diagnostics［J］.Journal of Sound and Vibration，2005，282:475-508.

［3］Blunt D M，Keller J A.Detection of a fatigue crack in a UH-60A planetgearcarrierusing vibration analysis［J］.Mechanical Systems and Signal Processing，2006，20:2095-2111.

［4］Rafiee J，Rafiee M A，Tse P W.Application of mother wavelet functionsforautomatic gearand bearing fault diagnosis［J］.Expert Systems with Applications，2010，37:4568-4579.

［5］Wang D Y，Zhang W Z，Zhang J G，et al.Fault bearing identification based on wavelet packet transform technique and artificial neural network［J］.International Conference on System Science，EngineeringDesign and Manufacturing Informatization，2010，2(s):11-14.

［6］ Wu J D，Liu C H.An expert system for fault diagnosis in internal combustion engines using wavelet packet transform and neural network［J］.Expert Systems with Applications，2009，36:4278-4286.

［7］ Wu J D，Liu C H.Investigation of engine fault diagnosis using discrete wavelet transform and neural network ［J］.Expert Systems with Applications，2008，35:1200-1213.

［8］ Saravanan N，Ramachandran K I.Incipient gear box fault diagnosis using discrete wavelet transform(DWT)for feature extraction and classification using artificial neural network(ANN)［J］.Expert Systems with Applications，2010，37:4168-4181.

［9］Wu J D，Huang C K，Chang Y W.Fault diagnosis for internal combustion engines using intake manifold pressure and artificial neural network ［J］.Expert Systems with Applications，2010，37:949-958.

［10］ Wu J D，Kuo J M.An automotive generator fault diagnosis system using discrete wavelet transform and artificial neural network［J］.Expert Systems with Applications，2009，36:9776-9783.

［11］ Rafiee J，Rafiee M A，Prause N，et al.Application of mother wavelet functionsforautomatic gearand bearing fault diagnosis［J］.Expert Systems with Applications，2010，37:4568-4579.

［12］ Sikorska J Z，Hodkiewicz M，Ma L.Prognostic modelling options for remaining useful life estimation by industry［J］.Mechanical Systems and Signal Processing，2011，25:1803-1836.

［13］ Wu S L，Bechhoefer E，He D.A practical regime prediction approach for HUMS applications［J］.American Helicopter Society 63rd Annual Forum，Virginia Beach，VA，2007:1-8.

［14］Peel L，Driven D.Data driven prognostics using a kalman filter ensembleofneuralnetwork models ［J］. International Conference on Prognostics and Helth Management，2008:1-6.

［15］Lall P，Lowe R，Goebel K.Prognostics using kalman-filter models and metrics for risk assessment in BGAs under shock and vibration loads［J］. Electronic Components and Technology Conference，2010:889-901.

［16］ Carr M J， Wang W B. An approximate algorithm for prognostic modelling using condition monitoring information［J］.European Journal of Operational Research，2011，211:90-96.

［17］ Matej G，Dani J，Pavle B，et al.Model-based prognostics of gearhealth using stochastic dynamicalmodels ［J］.Mechanical Systems and Signal Processing，2010，25(2):537-548.

［18］ GaingZ L. Wavelet-basedneuralnetworkforpower disturbance recognition and classification ［J］. IEEE Transactions on Power Delivery，2004，19(4):1560-1568.