韋堯兵,張學(xué)程
(蘭州理工大學(xué)機(jī)電工程學(xué)院,甘肅 蘭州 730050)
在旋轉(zhuǎn)機(jī)械變工況運(yùn)行過(guò)程中存在著豐富的故障信息,一些在平穩(wěn)工況運(yùn)行過(guò)程中不易反映出來(lái)的故障信息可能會(huì)在此時(shí)完全地顯現(xiàn)出來(lái),因此變工況運(yùn)行的旋轉(zhuǎn)機(jī)械振動(dòng)信息對(duì)旋轉(zhuǎn)機(jī)械故障的診斷具有極高的利用價(jià)值[1]。為了使旋轉(zhuǎn)機(jī)械安全穩(wěn)定地運(yùn)行,對(duì)旋轉(zhuǎn)機(jī)械進(jìn)行狀態(tài)監(jiān)測(cè)與故障診斷非常重要。目前,基于振動(dòng)信號(hào)的狀態(tài)監(jiān)測(cè)是主要的方式之一。旋轉(zhuǎn)機(jī)械升速、降速階段中的振動(dòng)信號(hào)為典型的非平穩(wěn)信號(hào),常規(guī)的分析方法已不適用。階次跟蹤法(order tracking technique,OT)是一種處理時(shí)變信號(hào)的方法,通過(guò)等角度重采樣將時(shí)域中非平穩(wěn)信號(hào)轉(zhuǎn)化為角域中平穩(wěn)信號(hào),然后再對(duì)處理后的平穩(wěn)信號(hào)進(jìn)行常規(guī)分析處理。然而,基于時(shí)-頻分布的階次跟蹤法[2]的重點(diǎn)和難點(diǎn)在于需要精準(zhǔn)地獲得參考軸轉(zhuǎn)速所對(duì)應(yīng)的瞬時(shí)頻率(instantaneous frequency,IF),它是準(zhǔn)確進(jìn)行階次跟蹤法中等角度重采樣的前提及關(guān)鍵。
在許多實(shí)際的信號(hào)分析處理過(guò)程中,估算一個(gè)時(shí)變過(guò)程的瞬時(shí)頻率是一項(xiàng)非常重要的工作。瞬時(shí)頻率的準(zhǔn)確估計(jì)是時(shí)-頻分布的關(guān)鍵,其估算方法很多,有相位法、峰值搜索法、過(guò)零點(diǎn)法、求根估算法等。相位估計(jì)法主要是用振動(dòng)信號(hào)的相位信息來(lái)得到瞬時(shí)頻率,分為相位差分法、相位建模法以及鎖相環(huán)法。相位差分法雖然計(jì)算簡(jiǎn)便,但是抗噪聲性能不高;相位建模法雖可以較精準(zhǔn)地估算出信號(hào)的瞬時(shí)頻率,但是其計(jì)算量大,對(duì)初始值要求嚴(yán)格,對(duì)適用場(chǎng)合要求較高;鎖相環(huán)法具有較高的抗噪能力,自適應(yīng)能力強(qiáng),但是不能追蹤變化過(guò)快的頻率。過(guò)零點(diǎn)法是用信號(hào)過(guò)零點(diǎn)的次數(shù)來(lái)求得信號(hào)的瞬時(shí)頻率,一般使用于周期性信號(hào)頻率的檢測(cè),在時(shí)變信號(hào)分析中其精準(zhǔn)度不高。求根估算法是在解方程組過(guò)程中利用零點(diǎn)與頻率的關(guān)系來(lái)獲得瞬時(shí)頻率,計(jì)算速度慢, 難以選取初始值[3]。
目前常用的瞬時(shí)頻率的估計(jì)方法是峰值搜索法[4],其基本思想是通過(guò)時(shí)-頻分布獲得待分析信號(hào)時(shí)間頻率的聯(lián)合分布,然后保留聯(lián)合分布中每個(gè)時(shí)刻的最大值所對(duì)應(yīng)的頻率,以此來(lái)獲得對(duì)信號(hào)瞬時(shí)頻率的估計(jì)。在轉(zhuǎn)速波動(dòng)工況下,旋轉(zhuǎn)機(jī)械的運(yùn)行可以看作是多次重復(fù)的升降速過(guò)程,因此振動(dòng)信號(hào)往往不表現(xiàn)出線性調(diào)頻特征。如果用線性調(diào)頻小波變換對(duì)其瞬時(shí)頻率進(jìn)行一階逼近,并不能夠準(zhǔn)確、完整地反映信號(hào)的頻率構(gòu)成。在線性調(diào)頻小波變換基礎(chǔ)上,Angrisani和D′Arco[5]對(duì)核函數(shù)進(jìn)行改進(jìn),引入頻率彎曲的概念,使線性調(diào)頻小波變換推廣為非線性調(diào)頻小波變換(nonlinear chirplet transform,NCT),增加了時(shí)-頻分布(time-frequency distribution,TFD)的自由度,提高了對(duì)非線 性調(diào)頻信號(hào)的刻畫(huà)能力。Peng、Yang 和方楊[6-8]等在此基礎(chǔ)上定義了更具普遍性的多項(xiàng)式核函數(shù)并進(jìn)行了深入的探討,使得NCT在分析信號(hào)時(shí)具有更高的時(shí)頻能量聚集性,因而抗噪性能更好,可以精確地估計(jì)信號(hào)的瞬時(shí)頻率。本文在NCT的基礎(chǔ)上提出了連續(xù)非線性調(diào)頻小波變換(CNCT),并將其應(yīng)用到瞬時(shí)頻率估計(jì)中,此方法首先通過(guò)基于NCT的峰值搜索法粗略估計(jì)出滾動(dòng)軸承的瞬時(shí)頻率,然后用最小二乘擬合方法對(duì)該頻率進(jìn)行連續(xù)擬合,當(dāng)達(dá)到設(shè)定的閾值時(shí)停止擬合,此方法有效地避免了交叉項(xiàng)的存在,提高了瞬時(shí)頻率的估計(jì)精度。
瞬時(shí)頻率的數(shù)學(xué)定義首先是基于Gabor解析信號(hào)提出的,在此基礎(chǔ)上又由Ville進(jìn)行了改進(jìn)[9-10]。現(xiàn)有一時(shí)域內(nèi)的連續(xù)信號(hào):
v(t)=β(t)cos[φ(t)]
(1)
式中:β(t)為信號(hào)的幅度;φ(t)為信號(hào)的相位。連續(xù)周期信號(hào)s(t)經(jīng)過(guò)Hilbert變換可得:
(2)
則其解析信號(hào)可表示為:
(3)
將瞬時(shí)頻率定義為解析信號(hào)的相位對(duì)時(shí)間的導(dǎo)數(shù),即:
(4)
此瞬時(shí)頻率定義的物理意義明確,解析信號(hào)為復(fù)數(shù)平面的一個(gè)向量,該向量幅角的轉(zhuǎn)速即為瞬時(shí)頻率。基于時(shí)-頻分布的瞬時(shí)頻率的另一種定義如下:
(5)
式中:w(t,f)為信號(hào)的時(shí)-頻分布, 即某時(shí)刻的瞬時(shí)頻率為該信號(hào)在此時(shí)的頻率加權(quán)平均,是時(shí)-頻分布的一階矩。
通過(guò)時(shí)-頻分布獲得的時(shí)變信號(hào)的瞬時(shí)頻率估計(jì)與工程實(shí)際更相符,因而得到了許多研究者的關(guān)注。
基于時(shí)-頻分布的瞬時(shí)頻率估計(jì)精度取決于時(shí)-頻平面內(nèi)的能量聚集程度,但凡信號(hào)具有局部能量的部位,它的時(shí)-頻分布也聚集在這些部位。高能量聚集的時(shí)-頻分布方法可提供更具抗噪能力、更加可靠的瞬時(shí)頻率估計(jì)。
引入時(shí)間平移、頻率平移等合成算子來(lái)表示頻率隨時(shí)間的非線性變化,時(shí)間平移、頻率平移、時(shí)間拉伸、時(shí)間傾斜、頻率傾斜可分別表示為:
(6)
式中:tc為時(shí)間中心;fc為頻率中心;Δt為尺度參數(shù);lgΔt為對(duì)數(shù)時(shí)寬;d為時(shí)間傾斜參數(shù);c為線性調(diào)頻斜率;q(t)為采用高斯函數(shù)作為窗函數(shù);*為卷積運(yùn)算符。合成算子中若某個(gè)參數(shù)等于零,則表明與之對(duì)應(yīng)的仿射變換不發(fā)生。
有了合成算子之后,將時(shí)頻分析中的基信號(hào)表示成:
qtc,fc,lgΔt,d,c(t)=Ntc,fc,lgΔt,d,cq(t)
(7)
式(7)稱(chēng)為分析函數(shù)。則時(shí)頻分布可表示成信號(hào)s(t)與分析函數(shù)的內(nèi)積形式,如短時(shí)傅立葉變換:
STFT(t,f)=
(8)
五維連續(xù)線性調(diào)頻小波變換:
CT(t,f)=
(9)
對(duì)于瞬時(shí)頻率呈非線性變化的信號(hào),擴(kuò)大其頻率傾斜算子,擴(kuò)大后的頻率傾斜算子稱(chēng)為頻率彎曲算子,如式(10)所示:
(10)
式中:cl(l=2,3,…,k)為頻率彎曲參數(shù),也即多項(xiàng)式系數(shù)。如果頻率彎曲算子的瞬時(shí)頻率與信號(hào)的瞬時(shí)頻率相一致,時(shí)頻分析就具有最佳時(shí)頻聚集性,所以只要選擇合適的頻率彎曲參數(shù),就可以實(shí)現(xiàn)對(duì)信號(hào)瞬時(shí)頻率的精確描述。
頻率彎曲算子的推廣令時(shí)頻分析描述信號(hào)頻率的變化更加自由,此時(shí)線性調(diào)頻小波變換可稱(chēng)為非線性調(diào)頻小波變換,是cn=0的特例(n=3,…,k)。在數(shù)學(xué)上,Weierstrass證明了多項(xiàng)式級(jí)數(shù)可對(duì)閉區(qū)間上的任意連續(xù)函數(shù)一致逼近,因此頻率彎曲算子可對(duì)任意頻率隨時(shí)間連續(xù)變化的信號(hào)進(jìn)行描述。
峰值檢測(cè)的目的是為了獲得各時(shí)刻的瞬時(shí)頻率,但是時(shí)頻面上的峰值檢測(cè)需具有良好時(shí)頻能量聚集[11]。根據(jù)時(shí)頻譜圖在參考軸轉(zhuǎn)頻附近進(jìn)行譜峰搜索,將峰值對(duì)應(yīng)的坐標(biāo)看做該時(shí)間點(diǎn)的瞬時(shí)頻率,算法如下:
ni=ni±1,ni±2,…,ni∈(0,N-1)
(11)
式中:arg max表示參數(shù)取最大值(在此為k);D為搜尋范圍;(ni,ki)為某坐標(biāo)對(duì)應(yīng)的峰值,也即ni時(shí)刻對(duì)應(yīng)的瞬時(shí)頻率為ki。為了選取最優(yōu)的頻率彎曲參數(shù),可先采用某一初始參數(shù)對(duì)信號(hào)進(jìn)行粗略NCT分析;然后利用最小二乘法對(duì)時(shí)頻譜峰進(jìn)行擬合,將擬合得到的多項(xiàng)式系數(shù)作為頻率彎曲參數(shù),再次對(duì)信號(hào)進(jìn)行NCT分析。NCT分析有限長(zhǎng)度的信號(hào)時(shí),由于能量泄漏產(chǎn)生端點(diǎn)效應(yīng),可通過(guò)自適應(yīng)波形匹配延拓改善端點(diǎn)效應(yīng),將失真部分隔離在待分析信號(hào)的外部[12]。
用最小二乘擬合法處理峰值檢測(cè)后得到的瞬時(shí)頻率數(shù)據(jù)。設(shè)該數(shù)據(jù)時(shí)-頻面上坐標(biāo)為(tk,fk),k=1,2,…,以二階方程擬合為例,其擬合方程為:
ui(t)=at2+bt+m
(12)
擬合誤差ε的平方和為:
(13)
(14)
式(14)對(duì)應(yīng)的數(shù)字算法為:
(15)
式中:ΔL=nΔt,Δt=1/us,us為瞬時(shí)頻率,由STFT時(shí)窗每次滑動(dòng)點(diǎn)數(shù)來(lái)確定n。
本文采用如下仿真信號(hào)進(jìn)行分析來(lái)證明連續(xù)非線性調(diào)頻小波變換在時(shí)頻分析方面的優(yōu)異性能。
對(duì)上述仿真信號(hào)運(yùn)用連續(xù)非線性調(diào)頻小波變換方法進(jìn)行處理,圖1、圖2為第一次NCT處理后的時(shí)-頻分布結(jié)果。從圖中可以看出,非線性調(diào)頻成分能量分布的頻帶相當(dāng)寬,所以該時(shí)-頻分布的能量聚集性較差,所對(duì)應(yīng)的瞬時(shí)頻率曲線波動(dòng)較大。當(dāng)進(jìn)行首次最小二乘法逼近后,其時(shí)-頻分布的能量聚集性有較大改善,瞬時(shí)頻率曲線波動(dòng)幅度 降低,如圖3、圖4所示。當(dāng)再次逼近該多項(xiàng)式系數(shù)時(shí),時(shí)-頻分布能量聚集性再次得到改善,如圖5、圖6所示,在時(shí)-頻分布圖上能量聚集性表現(xiàn)非常突出,由此估計(jì)的時(shí)變信號(hào)的瞬時(shí)頻率曲線更加靠近非線性調(diào)頻信號(hào)的準(zhǔn)確頻率。圖7、圖8是判定條件達(dá)到閾值時(shí)的最后一次逼近,此時(shí)瞬時(shí)頻率波動(dòng)已趨于穩(wěn)定。
圖1 第1次NCT后的時(shí)-頻分布
圖2 第1次NCT后的瞬時(shí)頻率
圖3 第2次NCT后的時(shí)-頻分布
圖4 第2次NCT后的瞬時(shí)頻率
圖5 第3次NCT后的時(shí)-頻分布
圖6 第3次NCT后的瞬時(shí)頻率
圖7 第4次NCT后的時(shí)-頻分布
圖8 第4次NCT后的瞬時(shí)頻率
基于CNCT的瞬時(shí)頻率估計(jì)法在低信噪比的情況下優(yōu)于傳統(tǒng)的瞬時(shí)頻率估計(jì)方法。采用最小二乘擬合方法可以進(jìn)一步有效地提高瞬時(shí)頻率估計(jì)精度。對(duì)于一般階次分析方法,本文的方法是對(duì)其的有力補(bǔ)充,尤其在旋轉(zhuǎn)機(jī)械處于變轉(zhuǎn)速工況的場(chǎng)合有突出的優(yōu)勢(shì)。仿真和實(shí)驗(yàn)測(cè)試驗(yàn)證了該方法的有效性。