李紀永，李舜酩，陳曉紅，江星星

（1．南京航空航天大學能源與動力學院，南京 210016；2．南京航空航天大學理學院，南京 210016）

基于雙樹復小波變換的轉靜碰摩響應分析

李紀永1，李舜酩1，陳曉紅2，江星星1

（1．南京航空航天大學能源與動力學院，南京 210016；2．南京航空航天大學理學院，南京 210016）

針對碰摩故障診斷中轉靜碰摩響應消噪及特征提取構造分段圓弧平滑閾值函數，利用小波、小波包及雙樹復小波（Dual-tree Complex Wavelet）分別對非平穩(wěn)信號降噪處理，以信噪比為指標判斷降噪效果。結果表明，構造的閾值函數消噪具有高效性。通過計算單轉子局部碰摩動態(tài)響應，用雙樹復小波對位移響應進行特征提取。與小波及小波包相比，雙樹實小波近似解析性及采樣插值具有互補特點，能完整提取動態(tài)響應特征，減少信息丟失。并用仿真及轉子碰摩實驗驗證該方法的有效性。

轉靜碰摩；降噪；動態(tài)響應；閾值函數；雙樹復小波

航空發(fā)動機中轉靜碰摩有兩類，即轉動部分碰摩部位不在轉軸本身如葉片等與靜止部分發(fā)生的碰摩及轉動部分發(fā)生在轉軸本身的碰摩。碰摩故障通常表現為其它故障的間接結果［1－2］，如轉子質量不平衡、熱彎曲、軸系不對中及轉子渦動失穩(wěn)等均可能表現為轉靜碰摩。碰摩發(fā)生時一般為先有局部碰摩，即在一個周期內轉子與定子發(fā)生一次或數次碰撞。沖擊產生接觸力及摩擦力，摩擦力將轉子轉動能量轉化為橫向振動同時伴有較復雜的振動現象，在某些條件下會出現混沌運動。

通過小波及小波包分析碰摩非平穩(wěn)信號［3－4］可提取細節(jié)成分，但仍有3種缺陷限制其應用［5］：①平移敏感性，即輸人信號平移會在小波分解系數中產生不可預計變化；②方向性差，只能提供3個空間方位；③相位信息缺乏。小波變換域的能量越集中說明小波分析效果越好，而該能量集中程度與小波幅頻、相頻特性有關。小波濾波器非理想截止特性及分解過程隔點采樣會致信號分解、重構時產生虛假頻率成分。

雙樹復小波變換（Dual-tree Complex Wavelet Transform，DTCWT）由Kingsbury［6］提出。Salesnick等［7］對雙樹復小波濾波器構造進行深入研究。雙樹復小波變換具有的良好抗頻帶混疊能力及平移不變特性可保證雙樹復小波變換在濾波器特性及周期復合沖擊特征檢測方面優(yōu)于經典離散小波變換、第二代小波變換及經驗模式分解等時頻域分析方法。

降噪可視為經典的非參數回歸，現有函數估計方法包括樣條曲線擬合估計、核函數估計及傅里葉展開估計等［8］。Donolo等［9－10］在小波變換基礎上提出的閾值降噪方法優(yōu)勢在于小波基函數空間多樣性。小波變換稀疏性及去相關性保證真實信號小波變換能量集中于有限的小波系數，而如何選取閾值及閾值函數為實現小波去噪算法中關鍵問題［11］。閾值計算方法［12］有：stein無偏估計、固定閾值及啟發(fā)式閾值；閾值函數有硬閾值、軟閾值。硬閾值含間斷點造成重構信號振蕩，重構信號不光滑；軟閾值所得小波系數與原小波系數存在固定偏差，影響重構信號逼近程度。為此，研究者構建不同閾值函數，包括折衷閾值函數、指數閾值函數等，但大多與真實小波系數有一定偏差。

本文采用stein無偏估計計算閾值，利用構造的閾值函數對雙樹復小波變換后非平穩(wěn)信號進行降噪處理，計算單轉子碰摩響應，并對碰摩位移響應用雙樹復小波進行特征提?。焕脝无D子及多轉子碰摩實驗驗證本文所提方法。

1 雙樹復小波分解

DTCWT由兩實值小波φh（t）與φg（t）構成復小波，即

式中：φg（t）＝H［φh（t）］，H為Hilbert算子。

兩正交小波函數組成Hilbert變換對的充要條件［5］為兩低通濾波器滿足半采樣延遲條件。其幅值與相位存在關系為

構造關鍵環(huán)節(jié)為使兩小波函數濾波器滿足式（2），構造方法有線性正交法、共同因子法及Q平移方法（quarter shift）。Q平移法設計的濾波器為

式（3）精確滿足幅值條件，近似滿足相位條件，其相位表達式為

對稱中心在點n＝0．5（N－1）－0．25附近。利用Q平移方法設計的雙樹復小波變換整體為精確線性相位，由于φc（t）的解析性，DTCWT具有近似平移不變與抗頻帶混疊等優(yōu)良特性。前者有益于多重諧波信號檢測；后者對提取周期性沖擊特征非常有利。

2 雙樹復小波降噪閾值及閾值函數

利用小波降噪核心為閾值及閾值函數選擇。含噪信號表示為：s＝f＋se，e為服從標準正態(tài)分布的高斯白噪聲，其強度可用穩(wěn)健估計方法計算，即σ＝中值運算利用stein無偏似然估計，將某次層系數平方按由小到大排列，獲得某向量W＝［w1，w2，…，wn］，其中w1≤…wn，n為小波系數個數。由此估算向量R＝［r1，r2，…，rn］，其中：

式中：0≤α≤1，不同a值對應不同閾值函數，實用性強，缺點亦明顯，即具有間斷點且與真實系數有一定偏差。

為克服軟、硬閾值缺點，構造分段圓弧平滑閾值函數，設圓1半徑r1，圓2半徑r2，幾何關系見圖1。由圖1看出，圓弧1半徑為r1，圓弧2半徑為r2，圓1圓心坐標為O1（1－r1，r1），圓2圓心坐標為O2（1＋r2，r1），半徑滿足幾何關系r1＋r2＝βt。其中0≤β≤1。

圖1 幾何關系示意圖Fig．1 Geometric relation schematic

分段圓弧平滑閾值函數可表示為

改進的閾值函數在各點處一階可導，利于數學處理，避免信號重構時出現振蕩。調節(jié)參數b及半徑r1即可得不同折衷閾值函數。b＝1時即可平滑過渡至真實小波系數。a＝0．5，b＝1時示意見圖2。

圖2 閾值函數Fig．2 Threshold function

利用典型非平穩(wěn)信號Heavy sine（初始信噪比6．4 dB）、Doppler（初始信噪比6．4 dB）、Bump（初始信噪比7．3 dB）、Blocks（初始信噪比9．6 dB）驗證分段圓弧平滑閾值函數降噪的有效性。信噪比依次（每組SNR順序從左至右，先第一行后第二行）見表1。綜合比較可知，雙樹復小波去噪性能較好，從而驗證分段光滑圓弧平滑閾值函數降噪的有效性。

表1 Heavy sine，Doppler，Bump，Blocks降噪后SNR值Tab.1 SNR of denoised Heavy sine Dopp ler，Bum p，Blocks

3 單轉子仿真碰摩響應雙樹復小波分析

Jeffcott單轉子碰摩模型見文獻［2］，本文重設參數為

式中：M為轉子質量；c為軸阻尼；k為軸剛度；ks為非線性剛度系數；f為摩擦系數；kc為定子徑向剛度。

計算轉子動態(tài)位移響應其分叉圖見圖3。令s＝ω／ωc，ωc為臨界轉速。分析s＝5．65時動態(tài)響應，利用Wolf方法求得lyapunov指數為0．022，時間歷程及龐加萊映射見圖4（a），響應呈非線性；其頻譜見圖4（b），存在的高頻率成分未完全顯示，如s＝84、95．2、112，所取信號頻率成份較豐富，能較好反映轉子碰摩故障現象。

圖3 轉子系統(tǒng)分叉特性Fig．3 Rotor system bifurcation characteristic

圖4 混沌信號時域及頻域特征Fig．4 Chaos response time and frequency characteristic

圖5 不同分解方法所得頻譜Fig．5 Spectrum via different extraction methods

利用小波、小波包及雙樹復小波對混沌信號進行分析，所得頻譜見圖5。由圖5（a）看出，小波分解各層存在嚴重頻譜泄漏；而圖5（b）中除頻譜泄漏外亦存在較嚴重頻率混疊，如第2層存在5個虛假頻率，第3層節(jié)點在高頻處存在虛假頻率，第4層在中高頻處存在虛假頻率，小波及小波包分解的頻率泄漏由濾波器非理想截止特性造成，小波包分解產生的虛假頻率由分解過程中下抽樣造成采樣率降低，從而導致頻率折疊及在重構過程中隔點插零操作形成。由圖5（c）看出，第2層提取的響應高頻成份包括s＝84，95．2，112，得益于波重構過程中雙樹分解系數互補，減小信息丟失。各層幾乎不存在虛假頻率。此因雙樹復小波帶通頻率響應只在正半頻率軸上（除第一層外），能減小頻率混疊。而由于濾波器非理想截止特性，第4層存在頻譜泄漏。

4 轉子實驗故障信號雙樹復小波分析

選單轉子及三轉子碰摩位移響應信號進行分析，其中三轉子實驗裝置示意圖見圖6。單轉子臨界轉速4 216 r／min、2 754 r／min時的碰摩位移及頻譜見圖7。由圖7看出，基頻45．9 Hz及2倍頻91．8 Hz清晰可辨，而若干高次倍頻成份淹沒在噪聲中。利用雙樹復小波進行降噪處理，選構造分段光滑閾值函數對信號重構，頻譜見圖8。由圖8看出，已提取信號的高頻細節(jié)信息，其中第二層含505．1 Hz，593．2 Hz，643 Hz等；第三層含328．3 Hz（7倍頻），367．4 Hz等；第四層含183．8 Hz（4倍頻），204．1 Hz；第5層為2倍頻，頻帶間有輕微混疊，與相同分解層下小波包分解混疊頻譜成份相比（圖9，線框內），頻率混疊成份大幅降低。

圖6 轉子裝置示意圖Fig．6 Rotor unit schematic

圖7 轉子位移響應時間歷程及頻譜Fig．7 Rotor displacement time history and its spectrum

圖8 位移頻譜Fig．8 Displacement spectrum

圖9 位移頻譜Fig．9 Displacement spectrum

三轉子臨界轉速為1 867 r／min、轉速為5 700 r／min時轉子振動較大，軸與軸套發(fā)生碰摩。位移傳感器測得碰摩信號及頻譜見圖10，除基頻95．8 Hz外，亦產生2倍頻及若干高次倍頻。

圖10 轉子位移響應時間歷程及其頻譜Fig．10 Rotor displacement time history and its spectrum

利用雙樹復小波進行降噪處理，選構造分段光滑閾值函數對信號重構所得頻譜見圖11，已提取信號的高頻細節(jié)信息，其中第二層含669．5 Hz，765．3 Hz，860 Hz，分別為7、8、9倍頻；第三層含382 Hz，為4倍頻；第四層含191．2 Hz，為2倍頻。各頻帶間無混疊，有用信息獲得充分顯示。

圖11 位移頻譜Fig．11 Disp lacement spectrum

5 結 論

（1）構造的分段圓弧平滑閾值函數在整個區(qū)間一階導數連續(xù)，能避免重構信號振蕩，可有效消除非平穩(wěn)信號中噪聲。

（2）利用雙樹復小波變換分析碰摩非線性動態(tài)響應信號提取的轉子動態(tài)響應高頻信號與小波及小波包變換相比，雙樹復小波變換能抑制模式混疊，并能提取更多細節(jié)信息。

（3）由兩平行但不同的離散小波組成的雙樹復小波變換其雙樹分解系數信息互補，可減少由下采樣造成的信息丟失，能完整提取碰摩故障信號特征；其雙樹實小波解析特性可抗頻帶混疊并具有平移不變特性，減少信號信息特征丟失，在故碰摩故障信號特征提取中效果較好。

［1］Han Q，Zhang Z，Wen B．Periodic motions of a dual-disc rotor system with rub-impact at fixed limiter［J］．Proceedings of the Institution of Mechanical Engineers，Part C：Journal of Mechanical Engineering Science，2008，222（10）：1935－1946．

［2］Chu F，Zhang Z．Bifurcation and chaos in a rub-impact Jeffcott rotor system［J］．Journal of Sound and Vibration，1998，210（1）：1－18．

［3］Lin J，Qu L．Feature extraction based on Morletwavelet and its application for mechanical fault diagnosis［J］．Journal of Sound and Vibration，2000，234（1）：135－148．

［4］Peng Z，He Y，Lu Q，et al．Feature extraction of the rubimpact rotor system by means of wavelet analysis［J］．Journal of Sound and Vibration，2003，259（4）：1000－1010．

［5］Kingsbury N．Complex wavelets for shift invariant analysis and filtering of signals［J］．Applied and Computational Harmonic Analysis，2001，10（3）：234－253．

［6］Selesnick IW，Baraniuk R G，Kingsbury N C．The dualtree complex wavelet transform［J］．Signal Processing Magazine，2005，22（6）：123－151．

［7］Selesnick IW，Baraniuk R G，Kingsbury N C．The dual-tree complex wavelet transform［J］．Signal Processing Magazine，IEEE，2005，22（6）：123－151．

［8］Chang SG，Yu B，VetterliM．Adaptivewavelet thresholding for image denoising and compression［J］．Image Processing，IEEE Transactions on，2000，9（9）：1532－1546．

［9］Donoho D L．De-noising by soft-thresholding［J］．Information Theory，IEEE Transactions on，1995，41（3）：613－627．

［10］Donoho D L，Johnstone I M．Adapting to unknown smoothness via wavelet shrinkage［J］．Journal of the American Statistical Association，1995，432（90）：1200－1224．

［11］吳定海，張培林，任國全，等．基于雙樹復小波包的發(fā)動機振動信號特征提取研究［J］．振動與沖擊，2010，29（4）：160－163．

WU Ding-hai，ZHANG Pei-lin，REN Guo-quan，et al．Feature extraction of an engine vibration signal based on dualtree wavelet package transformation［J］．Journal of Vibration and Shock，2010，29（4）：160－163．

［12］Stein C M．Estimation of themean of a multivariate normal distribution［J］．The Annals of Statistics，1981，9（6）：1135－1151．

Rotor-stator rub-impact response analysis based on DTCWT

LI Ji-yong1，LIShun-ming1，CHEN Xiao-hong2，WANG Xing-xing1
（1．College of Energy and Power Engineering，Nanjing University of Aeronautics and Astronautics，Nanjing 210016，China；2．College of Science，Nanjing University of Aeronautics and Astronautics，Nanjing 210016，China））

The de-noising of rub-impact response signals of a rotor-stator system and its characteristics extraction are key problems of rub fault diagnosis．A piecewise smooth threshold function was constructed forwavelet de-noising，and non-stationary signals were de-noised separately by using Wavelet，Wavelet packets and dual-tree complex wavelet transform（DTCWT）．SNR（signal to noise ratio）was taken as an index to evaluate the de-noising effect．The results reveal that the constructed threshold function is highly active in de-noising．As an example，the dynamic response of a single rotor under rub-impactwas calculated and decomposed by DTCWT．The simulation and experiment results show that the dynamic response can be exactly derived due to the complementary property between approximate analyticity and sampling interpolation of DTCWT．It is also proved that the proposedmethod is effective．

rotor-stator rub；de-noising；dynamic response；threshold function；DTCWT

V233．1

：A

10．13465／j．cnki．jvs．2014．22．019

航空自然科學基金（2012ZD52054）；國家自然科學基金項目（61403193）；南航基本科研業(yè)務費科研項目（NS2014081）

2013－07－24 修改稿收到日期：2013－11－21

李紀永男，博士生，1984年生

李舜酩男，教授，1962年生郵箱：smli＠nuaa．com