基于雙樹復(fù)小波變換的轉(zhuǎn)靜碰摩響應(yīng)分析

李紀(jì)永，李舜酩，陳曉紅，江星星

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

基于雙樹復(fù)小波變換的轉(zhuǎn)靜碰摩響應(yīng)分析

李紀(jì)永1，李舜酩1，陳曉紅2，江星星1

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

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

轉(zhuǎn)靜碰摩；降噪；動態(tài)響應(yīng)；閾值函數(shù)；雙樹復(fù)小波

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

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

雙樹復(fù)小波變換（Dual-tree Complex Wavelet Transform，DTCWT）由Kingsbury［6］提出。Salesnick等［7］對雙樹復(fù)小波濾波器構(gòu)造進(jìn)行深入研究。雙樹復(fù)小波變換具有的良好抗頻帶混疊能力及平移不變特性可保證雙樹復(fù)小波變換在濾波器特性及周期復(fù)合沖擊特征檢測方面優(yōu)于經(jīng)典離散小波變換、第二代小波變換及經(jīng)驗?zāi)Ｊ椒纸獾葧r頻域分析方法。

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

本文采用stein無偏估計計算閾值，利用構(gòu)造的閾值函數(shù)對雙樹復(fù)小波變換后非平穩(wěn)信號進(jìn)行降噪處理，計算單轉(zhuǎn)子碰摩響應(yīng)，并對碰摩位移響應(yīng)用雙樹復(fù)小波進(jìn)行特征提??；利用單轉(zhuǎn)子及多轉(zhuǎn)子碰摩實驗驗證本文所提方法。

1 雙樹復(fù)小波分解

DTCWT由兩實值小波φh（t）與φg（t）構(gòu)成復(fù)小波，即

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

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

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

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

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

2 雙樹復(fù)小波降噪閾值及閾值函數(shù)

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

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

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

圖1 幾何關(guān)系示意圖Fig．1 Geometric relation schematic

分段圓弧平滑閾值函數(shù)可表示為

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

圖2 閾值函數(shù)Fig．2 Threshold function

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

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

3 單轉(zhuǎn)子仿真碰摩響應(yīng)雙樹復(fù)小波分析

Jeffcott單轉(zhuǎn)子碰摩模型見文獻(xiàn)［2］，本文重設(shè)參數(shù)為

式中：M為轉(zhuǎn)子質(zhì)量；c為軸阻尼；k為軸剛度；ks為非線性剛度系數(shù)；f為摩擦系數(shù)；kc為定子徑向剛度。

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

圖3 轉(zhuǎn)子系統(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

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

4 轉(zhuǎn)子實驗故障信號雙樹復(fù)小波分析

選單轉(zhuǎn)子及三轉(zhuǎn)子碰摩位移響應(yīng)信號進(jìn)行分析，其中三轉(zhuǎn)子實驗裝置示意圖見圖6。單轉(zhuǎn)子臨界轉(zhuǎn)速4 216 r／min、2 754 r／min時的碰摩位移及頻譜見圖7。由圖7看出，基頻45．9 Hz及2倍頻91．8 Hz清晰可辨，而若干高次倍頻成份淹沒在噪聲中。利用雙樹復(fù)小波進(jìn)行降噪處理，選構(gòu)造分段光滑閾值函數(shù)對信號重構(gòu)，頻譜見圖8。由圖8看出，已提取信號的高頻細(xì)節(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，線框內(nèi)），頻率混疊成份大幅降低。

圖6 轉(zhuǎn)子裝置示意圖Fig．6 Rotor unit schematic

圖7 轉(zhuǎn)子位移響應(yīng)時間歷程及頻譜Fig．7 Rotor displacement time history and its spectrum

圖8 位移頻譜Fig．8 Displacement spectrum

圖9 位移頻譜Fig．9 Displacement spectrum

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

圖10 轉(zhuǎn)子位移響應(yīng)時間歷程及其頻譜Fig．10 Rotor displacement time history and its spectrum

利用雙樹復(fù)小波進(jìn)行降噪處理，選構(gòu)造分段光滑閾值函數(shù)對信號重構(gòu)所得頻譜見圖11，已提取信號的高頻細(xì)節(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 結(jié) 論

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

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

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

［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ù)小波包的發(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

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

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

李紀(jì)永男，博士生，1984年生

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