轉軸轉速波動信號的時頻脊線提取方法

龔思琪，李舜酩，庾天翼

（南京航空航天大學 能源與動力學院，南京 210016）

變轉速是旋轉機械設備的常見復雜工況，其核心旋轉零部件在變轉速工況下的故障振動信號具有明顯的非平穩特點［1］。基于振動信號處理的故障診斷方法對設備狀態監測和故障診斷都至關重要［2］。

傳統的基于振動信號處理的故障診斷方法大多以不同狀態下的時域、頻域及時頻域信息作為參考，通過對它們的統計特征指標對比分析設備當下的健康狀態及潛在的故障特征［3－4］。在變轉速下的故障振動信號的信號表征與故障模式之間的對應關系更復雜，導致傳統的傅里葉變換并不能完全揭示變轉速信號的局部細節特征。在近20年來，隨著非平穩信號處理理論的發展和電子硬件計算技術的快速成長［5］，對變轉速下機械設備動態信號處理及故障診斷方法開展了廣泛的研究。研究者對小波變換［6］和高精度的頻譜與時頻分析方法［7－8］等故障診斷技術進行了廣泛研究，并綜合運用統計物理學方面的信號特征提取方法［9］、盲源分離方法［10］、變分模式分解［11］等方法提高診斷準確率。瞬時頻率變化脊線是振動信號的重要特征之一。通常利用峰值檢測法提取時頻脊線［12］，文獻［13］提出了非線性擠壓時頻變換方法提取旋轉機械故障振動信號中的多分量非穩態信號。但在受到強噪聲干擾的變轉速工況下的振動信號時頻脊線提取和瞬時轉速提取時，準確率較低，辨識效果較差［13］，且對于含有多個頻率分量的信號，脊線難以辨別。文獻［14］提出了一種雙向搜索時頻脊融合方法，提取時頻變工況行星齒輪箱高速軸的瞬時轉速，可靠性高，但程序編寫困難，工作量大，而且需要豐富的理論知識。

針對噪聲干擾下的變轉速旋轉機械振動信號時頻脊線提取和沖擊特征識別問題，本文在利用短時傅里葉變換得到轉速波動信號的時頻分布圖；再將其沿時間軸分割，利用提取中心頻率的方法逐段提取時頻脊線；然后，利用卡爾曼濾波技術提取整段信號的最優脊線。通過分析脊線特點，可以準確有效地判斷轉子故障發生的時間。進一步，利用本文所提方法提取含有多個頻率分量信號的時頻脊線特征；最后，利用變轉速故障軸承的振動信號驗證本文所提方法的有效性和準確性。結果表明：所提方法的時頻脊線提取準確率高，操作簡單，易于實現，對轉速波動信號的特征提取有重要的參考價值。

1 轉軸故障振動信號的頻譜分析

1．1 實驗信號采集

本研究所使用的試驗臺如圖1所示，由電機驅動，在600～1 500 r／min的轉速下工作，采用油脂潤滑。驅動電機的額定功率為0．75 kW，額定電壓為380 V，額定轉矩是5．0 kN·m。

圖1 傳動軸變轉速故障振動試驗臺及傳感器放置位置

為了測量轉軸在變轉速下的振動信號，在轉軸驅動端和輸出端各安裝一個振動三向加速度傳感器采集傳動變轉速振動信號，傳感器安裝位置如圖1（a）、（b）所示。本次變轉速以勻加速方式實現，使電機以600 r／min的轉速開始工作，勻加速至1 500 r／min。測試所用的振動傳感器型號和基本參數如表1所示。

表1 振動傳感器的基本參數

在試驗過程中，采用LMS振動信號采集系統采集被測轉軸的故障振動信號。根據轉速和采樣頻率可知信號一個周期，即試驗轉軸旋轉一周所采得的數據點個數為2 048～5 120。

在測試過程中，選擇5 s內的振動信號作為分析目標，測試得到的變轉速轉軸的不平衡故障振動信號波形如圖2所示。

圖2 轉速波動信號時域波形

1．2 信號的頻譜特點

根據式（1）所示的傅里葉變換（fourier trans form，FT）公式，通過改進可將信號通過快速傅里葉變換，得到圖3所示的轉速波動信號的頻譜圖。

從圖3轉速波動信號的頻譜圖可以看到：信號主要集中在小于1 000 Hz的范圍內，并伴有噪聲干擾，可以看到頻率在0～400 Hz出現了最大幅值，且有多個頻率分量。

圖3 轉速波動信號的頻譜圖

2 信號的時頻分析

根據式（2）所示的短時傅里葉變換（short time fourier transform，STFT）的表達式，可得到時頻分布圖。

式中：STFT（t，w）表示STFT分析結果；f（t）為待分析信號；g（t）表示STFT窗函數。

窗口函數g（t）的選擇決定了分辨率的高低。根據不同窗函數的特點可知：具有帶通特性的漢寧窗或漢明窗雖然對信號濾波有更好的效果［15］，但是會使能量幾乎集中在主瓣。而高斯窗函數可以取得時間和頻率分辨率的最佳折衷［16］。因此，選取高斯窗函數對信號進行STFT。

對變轉速故障振動信號進行STFT分析，得到如圖4所示的二維時頻圖。可以看出：頻率在0～1 000 Hz和4 000 Hz附近的能量高，頻率隨時間基本保持不變，它們是由轉子的固有振動模態引起的。但在6～15 s的時間范圍內，如圖4中的頻率線2和3以及它們附近的多個頻率出現了隨時間變化而不斷上升的趨勢。這表示隨著采樣時間的增加，轉軸的轉動頻率也逐漸增加，這是由于轉子在6～8 s產生波動導致的頻率沿采樣時間的增加而升高。

圖4 轉速波動信號的二維時頻圖

由于轉子轉速在600～1 500 r／min波動，轉子轉頻最高在25 Hz附近，通過圖3的轉速波動信號的頻譜圖可以判斷轉速波動信號的頻率主要集中在0～400 Hz。從圖5可以看到：頻率不隨時間變化的線4，以及在時間大于8 s時，頻率隨時間變化的線5、線6。其中，線4的頻率是56 Hz，大約是轉子頻率的3～4倍。

圖5 時頻圖4的局部放大圖

3 信號的脊線提取

對信號進行STFT后，一些與機械系統零部件的動力學特性相關的曲線會清晰或隱約地出現，這些曲線在三維空間內的分布形態與“山脊”很像，是脊線標志或沿脊頂延伸的線。因此，在時頻圖中把三維空間內沿功率最高點延伸的線稱之為時頻脊線［17］，表現在時頻面上，是明亮程度和趨勢各異的線條，它能夠表征信號頻率變化的各種模式［18］。

信號在時域和頻域的能量可表示為：

能量分布還滿足邊緣性質：

從時頻表示上直接提取的時頻脊線受噪聲影響大，脊線特征不明顯，從峰峰值上得到的脊點，噪聲難以消除，導致所需提取的脊線淹沒在噪聲的峰值點中，如圖6所示。

圖6 轉速波動信號提取時頻脊線

從時頻能量最高點上提取信號的中心主頻率，遍歷時間軸，會得到圖7中的線1、2、3。在0～6 s的時間段內，轉速波動信號的頻率主要在56 Hz附近有小范圍的波動；在6～9 s時間段內，中心頻率與時間的關系變得很復雜；在9～12 s時間段內，中心頻率集中在二階頻率分量處有上升趨勢；在12～20 s的時間段內，中心主頻率在25 Hz附近隨著時間的增加而升高；在12 s處頻率有突變，如圖7中區域3所示。導致這種突變產生的主要原因是在6～8 s時間段內產生故障后，信號的頻率分量與能量的對應關系模糊，使得部分頻率分量的能量突然間提高，而其他頻率分量的能量降低。

圖7 轉速波動信號脊線圖

時頻面內脊線附近局部區域內的時頻點能量突出，具有相對較高的信噪比，機械系統運行過程中狀態一般不會出現突變，出現的時頻脊線具有連續平滑的特點。因此，直接通過能量最高點來提取時頻脊線的方法失效。但該方法表示出了該轉速波動信號的能量變化，為后面沿時間域分割信號的分割點提供了依據。

根據脊線過渡平滑的特點，以及圖7直接提取的中心頻率的劇烈變化，將提取區域按時間劃分為3個部分：0～6 s，6～9 s，9～20 s。目標頻率范圍設置為（f－Δf，f＋Δf），Δf＝5 Hz由于轉子的轉速為600～1 500 r／min，將圖7中線2的12～20 s所示的脊線即表示一階頻率作為即將要提取的目標脊線，通過對時間軸分段提取，將提取的頻率閾值設置在接近一階中心頻率附近，再根據能量最高點提取頻率范圍內的曲線。

圖8 轉速波動信號目標脊線

由圖7的信號脊線圖可以看出：56 Hz的頻率大概是轉速信號的3～4倍頻。由此可推算出：一階頻率在14 Hz附近，在0～6 s將目標頻率f設置為14 Hz，6～9 s頻率變化范圍大，且初步確定是轉速增加的，故將目標頻率設置為f＝18 Hz。9～20 s也是一個轉速上升的過程，但相比第2段，其過渡較為平穩，最高轉速1 500 r／min，頻率為25 Hz，將其目標頻率設置為f＝25 Hz，得到如圖8線1所示的脊線。

從圖8中的線1可以看到：在0～5 s目標脊線在這一時間段的頻率隨時間波動比較劇烈，結合圖7中該段時間內的能量集中在56 Hz附近，因此在圖8中該段時間內在該處能量比較分散會引起較大的波動。目標脊線在5～10 s之間的過渡十分復雜，時間在10～20 s之間的目標脊線波動較小，但頻率隨時間增加的趨勢也更為明顯。

4 最優脊線提取

卡爾曼濾波器是一種最優化自回歸數據處理算法，對于解決很大部分問題都會是最優、效率最高甚至是最有用的［20］方法。因此，本文采用卡爾曼濾波在時頻面內提取目標脊線的最優脊線［21］。

Kalman濾波器的系統測量值、狀態運動方程和觀測方程分別為［22］：

式中：X（k）是k時刻的系統狀態；U（k）是k時刻對系統的控制量；A是狀態轉移矩陣；Z（k）是k時刻的觀測值；H是觀測矩陣；W（k）和V（k）表示上述兩個方程中的噪聲。

在實際問題中，通常假設A（k）＝A，H（k）＝H。設W和V為噪聲（white gaussian noise），協方差分別為Q，R（這里假設它們不隨系統狀態變化而變化）。首先，利用狀態方程預測下一狀態的系統。假設當前時刻為k，根據系統的模型，可以基于系統的上一狀態而預測出現在狀態：

式中：X（k｜k－1）為利用上一狀態的預測結果；X（k－1｜k－1）為上一狀態的最優結果；U（k）為現在狀態的控制量，若沒有控制量，它可以為0。

對應于X（k｜k－1）的協方差，即上面預測值相對于真實狀態的誤差方差陣更新公式為：

式中：P（k｜K－1）是X（k｜k－1）對應的協方差；P（k－1｜k－1）是X（k－1｜k－1）對應的協方差；A′表示A的轉置矩陣；Q是系統過程的協方差。

當第k時刻的觀測值Z（k）到達以后，利用預測值中沒有的信息去修正第k個時刻的狀態預測值。結合預測值和測量值，可以得到狀態（k）時最優化估算值X（k｜k）：

式中，Kg為卡爾曼增益（Kalman Gain），可通過式（13）估計：

為了使Kalman濾波器運行直到系統過程結束，需要更新k狀態下X（k｜k）的協方差：

式中，I為單位矩陣，對于單模型單測量，I＝1。當系統進入k＋1狀態時，P（k｜k）就是式（11）的P（k－1｜k－1）。

通過以上卡爾曼濾波器的5個核心公式不斷進行迭代優化運算，可以完整實現Kalman濾波。將原始脊線的初始坐標作為卡爾曼濾波的初始狀態輸入，將P的初值設為1，保證卡爾曼濾波能夠穩定運行。通過多次試驗可得協方差Q與R分別為10－8，10－5時，Kalman濾波的濾波效果最佳。根據脊線平滑且表示信號的頻率隨時間變化趨勢的特點，得到圖9所示的轉速波動信號的最優目標脊線。圖9中藍色曲線表示經過沿時間軸分割后提取的原始脊線，紅色曲線表示經過卡爾曼濾波技術后得到的最優脊線，橙色曲線表示兩曲線之間的誤差。為了使脊線信息更加明顯，將縱坐標范圍縮小至0～45 Hz。相比原始脊線，最優脊線更加平滑，且其表達的頻率隨時間變化的趨勢沒有改變。在0～4 s時間段內，可以看到經卡爾曼濾波器提取的脊線頻率更加集中，這也更加真實地反映了在這段時間內轉子正常運轉。而誤差曲線主要集中在0 Hz附近，說明最優脊線與實際信號誤差小，最優脊線能反映實際信號波動情況的真實信息。

圖9 轉速波動信號最優脊線

5 結論

1）在STFT時頻分析的基礎上，結合能量最高點提取中心主頻率的方法，提取轉速波動信號的時頻脊線，可以觀察到信號頻率在不同時刻發生的變化，并可以判斷轉子系統是否發生故障以及故障發生的時間。相比直接提取時頻脊，其所受噪聲影響更小。

2）通過對不平穩信號的時頻脊線進行卡爾曼濾波，可以有效地找到目標脊線的最優解，得到表達信號頻率隨時間變化的趨勢的平滑的時頻脊線。

3）通過對轉軸試驗臺實測信號變換結果的觀察分析可知：從故障振動信號的時頻譜圖提取時頻脊線，并利用卡爾曼濾波估計最優脊線，可以有效地濾除背景噪聲，剔除干擾分量，得到平滑且清晰的時頻脊線。