葉杰凱，湯小明，林青云，梁登，吳華盛

(麗水市特種設(shè)備檢測(cè)院，浙江 麗水 323000)

0 引言

信號(hào)處理是一種利用數(shù)學(xué)算子從信號(hào)中提取特征信息的反向處理方法。傳統(tǒng)的傅里葉方法只能處理平穩(wěn)信號(hào)，并且是基于全局的頻譜分析，對(duì)于頻率隨時(shí)間變化的非平穩(wěn)信號(hào)往往無能為力。而時(shí)頻分析方法能夠有效地揭示時(shí)變性特征，已被廣泛應(yīng)用于從非平穩(wěn)信號(hào)中提取特征信息[1-2]。

傳統(tǒng)的時(shí)頻分析方法有連續(xù)小波變換(CWT)和短時(shí)傅里葉變換(STFT)，常被用來描述信號(hào)在時(shí)頻平面上的特性[3]。然而，它們都有同樣的局限性，就是所謂的“不確定性原理”，即不能在時(shí)間和頻率上任意精確地定位一個(gè)信號(hào)，也就是常說的時(shí)頻模糊問題。這種局限性導(dǎo)致后續(xù)的瞬時(shí)頻率提取和信號(hào)重構(gòu)等方面會(huì)存在不足，難以較為清晰地提取出故障信號(hào)及其分量。因此，針對(duì)時(shí)頻平面去模糊化的研究一直都是一個(gè)熱門方向，時(shí)頻重排和同步壓縮變換(synchrosqueezing transform, SST)就是其中一個(gè)較好的思路。

自從DAUBECGES I等[4]提出了同步壓縮變換來揭示非平穩(wěn)信號(hào)復(fù)雜的時(shí)頻特性，在醫(yī)學(xué)、地球物理學(xué)和音頻等多個(gè)領(lǐng)域?qū)ζ涠加兴鶓?yīng)用。SST是一種特殊的重分配算法，它對(duì)應(yīng)于一種提高時(shí)頻分辨率的非線性算子。在SST中，時(shí)頻系數(shù)僅在頻率軸上重新分配，使其既能夠簡(jiǎn)化分配過程又能與原始參量保持聯(lián)系[5-6]。然而，SST只能使瞬時(shí)頻率“慢變”信號(hào)的時(shí)頻表示銳化，并且當(dāng)處理強(qiáng)時(shí)變信號(hào)如非線性調(diào)頻信號(hào)時(shí)，SST不能產(chǎn)生集中的時(shí)頻表達(dá)。這是因?yàn)镾ST中的頻率重分配算子不能為強(qiáng)時(shí)變情況提供準(zhǔn)確的無偏估計(jì)。為此，PHAM D H等[7]提出了高階同步壓縮變換方法(high-order synchrosgueezing transform, HSST)。在高階振幅和相位近似的基礎(chǔ)上定義了新的同步壓縮算子，并對(duì)快速變化的瞬時(shí)頻率(IF)信號(hào)產(chǎn)生高度集中的時(shí)頻表達(dá)。但是，HSST是以傳統(tǒng)短時(shí)傅里葉變換STFT為基礎(chǔ)，在時(shí)頻變換窗寬上無法做到自適應(yīng)選擇。同時(shí)，BERRIAN A等[8]提出了基于縫式短時(shí)傅里葉變換(quilted STFT, QSTFT)的同步壓縮變換方法，即SST-QSTFT。該方法對(duì)復(fù)雜多組分信號(hào)的瞬時(shí)頻率無法做到準(zhǔn)確估計(jì)，再加上對(duì)噪聲極為敏感，因此在強(qiáng)背景噪聲下，很難獲得清晰的時(shí)頻表示[9-10]。但是，QSTFT重新定義了一組隨時(shí)頻變化的自適應(yīng)窗口函數(shù)集合，具有明顯的自適應(yīng)特性，其可以使用最合適的窗函數(shù)來適應(yīng)不同的時(shí)變信號(hào)，自動(dòng)匹配信號(hào)的局部變化，以此來提高STFT的時(shí)頻分辨率。

綜上，本文引入基于縫式短時(shí)傅里葉變換QSTFT的自適應(yīng)窗口計(jì)算方法，并將該方法應(yīng)用到HSST中，提出一種基于縫式短時(shí)傅里葉變換的自適應(yīng)高階同步壓縮變換(adaptive high-order synchrosqueezing transform based on a quilted short-time fourier transform, AHSST-QSTFT)。

1 理論描述

1.1 基于QSTFT的同步壓縮變換

定義一個(gè)調(diào)幅調(diào)頻(AM-FM)信號(hào)，令其表達(dá)式為

s(t)=A(t)eiφ(t)

(1)

式中A(t)和φ(t)分別是其瞬時(shí)幅值和瞬時(shí)相位。

繼而信號(hào)的STFT變換可以寫成下面的形式：

(2)

式中：g(t)是窗函數(shù)；g*(t)是g(t)的復(fù)共軛。

那么，SST-STFT可以用下面的公式表示：

(3)

(4)

將時(shí)頻面任意一個(gè)局部區(qū)域定義為Ω?R2,且該區(qū)域?qū)?yīng)一個(gè)縫窗口函數(shù)集合hΩ。令ht,ω=hΩ，對(duì)其中每一個(gè)(t,ω)∈Ω。則函數(shù)集合h的表達(dá)式為

h(τ,t,ω)=ht,ω(τ)

(5)

式中t,τ∈R且ω∈R+。

對(duì)信號(hào)s(t)，本文研究的QSTFT的表達(dá)式為

(6)

因此，由式(3)和式(6)可得SST-QSTFT表達(dá)式為

(7)

1.2 高階同步壓縮變換

SST采用零階估計(jì)來對(duì)信號(hào)的頻率進(jìn)行估計(jì)，僅在處理慢時(shí)變的中頻信號(hào)時(shí)效果較好。而高階SST方法則是基于信號(hào)振幅和相位的高階Taylor展開，對(duì)信號(hào)的IF進(jìn)行更加精確的估計(jì)，以此來提高在特征頻率上面的能量，達(dá)到去模糊化的目的。

AM-FM信號(hào)的Taylor展開式如下：

(8)

式中φ(k)(t)表示φ(t)的第k階導(dǎo)數(shù)。

那么，式(2)可以改寫成如下形式：

(9)

(10)

(11)

此時(shí)，HSST的表達(dá)式可以寫成

(12)

1.3 基于QSTFT的自適應(yīng)高階同步壓縮變換

(13)

最后，根據(jù)式(7)、式(11)和式(13)可知，AHSST-QSTFT的表達(dá)式為

(14)

因此，本文提出的基于QSTFT的自適應(yīng)高階同步壓縮變換算法AHSST-QSTFT的具體實(shí)現(xiàn)步驟如下：

1)定義縫窗口函數(shù)，對(duì)STFT的窗函數(shù)進(jìn)行改進(jìn)，得到自適應(yīng)窗函數(shù)集合；

3)按照式(14)的方法重置高階瞬時(shí)頻率估計(jì)，從而獲得改進(jìn)的自適應(yīng)高階時(shí)頻表達(dá)QHTs(t,ω)，并畫出時(shí)頻圖。

2 仿真分析

為了驗(yàn)證該方法的可行性，構(gòu)造兩個(gè)具有AM-FM性質(zhì)的分量信號(hào)：

(15)

其中：信號(hào)的采樣頻率設(shè)為4 096 Hz；信號(hào)長(zhǎng)度為4 096；信號(hào)的理想瞬時(shí)頻率曲線如圖1所示(本刊為黑白印刷，如有疑問請(qǐng)咨詢作者)。為了更加貼近實(shí)際情況，使模擬信號(hào)更有代表性，將兩個(gè)信號(hào)進(jìn)行了疊加，并且添加了SNR=-2的高斯白噪聲(信噪比公式SNR=10lg(Ps/Pn)，得到信號(hào)x(t)的時(shí)域波形如圖2所示，頻域波形如圖3所示。由FFT計(jì)算的結(jié)果可以看出，在噪聲干擾較大的情況下，無法從頻譜圖中區(qū)分和辨識(shí)時(shí)變的信號(hào)特征。

圖1 信號(hào)x1(t)和x2(t)的頻率曲線

圖2 信號(hào)x1(t)、x2(t)和x(t)的時(shí)域波形圖

圖3 信號(hào)x(t)的頻域波形圖

因此，為了說明本文所提出算法的優(yōu)勢(shì)，分別畫出了仿真信號(hào)的4階HSST和AHSST-QSTFT的時(shí)頻圖，并且采用相同的脊曲線提取算法對(duì)分量進(jìn)行了提取，結(jié)果見圖4。其中虛線為真實(shí)的瞬時(shí)頻率，實(shí)線為提取的脊線。

圖4 時(shí)頻表示及其脊線提取結(jié)果

通過兩組時(shí)頻圖的對(duì)比，可以得出AHSST-QSTFT方法能夠獲得更好的時(shí)頻表達(dá)效果，并且提取到的瞬時(shí)頻率更加接近于真實(shí)情況。

3 實(shí)驗(yàn)分析

故障診斷對(duì)于機(jī)械設(shè)備的預(yù)知維修具有重要意義。根據(jù)RANDALL R B等提出的滾動(dòng)軸承模型[11]，為了進(jìn)一步驗(yàn)證基于AHSST-QSTFT方法在機(jī)械設(shè)備故障診斷中的有效性，本文設(shè)計(jì)了如下實(shí)驗(yàn)。

圖5為試驗(yàn)裝置的故障模擬器。設(shè)備由550 W的交流電機(jī)驅(qū)動(dòng)。故障軸承安裝在最右側(cè)的皮帶輪上，利用電火花加工在目標(biāo)軸承加工4個(gè)點(diǎn)蝕坑，坑的深度為0.8 mm，并盡可能靠近支撐軸承的機(jī)殼上方垂直固定振動(dòng)加速度傳感器，用以進(jìn)行數(shù)據(jù)采集，傳感器位置如圖中箭頭所示。實(shí)驗(yàn)采用B&K3560C數(shù)據(jù)采集分析儀進(jìn)行數(shù)據(jù)采集，采樣頻率為16 384 Hz，轉(zhuǎn)速1 450 r/min；滾動(dòng)體個(gè)數(shù)Z=9，試驗(yàn)中的滾動(dòng)軸承型號(hào)為6207，壓力角為α=0°，故障類型為軸承外圈故障，可以根據(jù)下面的計(jì)算得到故障頻率：

圖5 軸承故障模擬實(shí)驗(yàn)臺(tái)

(16)

式中：fr為轉(zhuǎn)頻，Hz；fo為軸承外圈故障特征頻率，Hz。

選擇8 192個(gè)采樣數(shù)據(jù)分析，實(shí)際測(cè)得的軸承振動(dòng)信號(hào)如圖6所示，在時(shí)域波形上能夠明顯看出瞬態(tài)沖擊故障。圖7是該信號(hào)的頻譜分析，可以觀察到峰值區(qū)域主要集中在1 000 Hz、4 000 Hz附近。這是因?yàn)檩S承發(fā)生故障時(shí)，當(dāng)轉(zhuǎn)子經(jīng)過缺陷位置時(shí)，相當(dāng)于是在軸承上加上了一個(gè)沖擊力所激發(fā)出來的信號(hào)，反映的是該系統(tǒng)的固有頻率。

圖6 實(shí)測(cè)信號(hào)的時(shí)域波形

圖7 實(shí)測(cè)信號(hào)的頻譜圖

考慮到實(shí)際信號(hào)的能量都集中在共振區(qū)間，不利于最終故障的識(shí)別。因此在本文中首先對(duì)實(shí)際信號(hào)進(jìn)行希爾伯特變換，得到對(duì)應(yīng)的包絡(luò)信號(hào)，然后再采用AHSST-QSTFT方法對(duì)信號(hào)進(jìn)行時(shí)頻分析。

圖8(a)和圖8(b)分別是在相同的參數(shù)設(shè)置下采用HSST和AHSST-QSTFT方法得到的時(shí)頻表示?？紤]到軸承外圈的故障頻率為87.01 Hz，一般地，如果能夠找到故障特征頻率對(duì)應(yīng)的多個(gè)倍頻就能夠確定故障。圖8(a)雖然能夠較為模糊地看到某些對(duì)應(yīng)的特征頻率，但受到背景噪聲的影響，它們?cè)跁r(shí)頻平面上表示得并不是很突出。圖8(b)為采用提出的AHSST-QSTFT方法得到的時(shí)頻表示，相比于圖8(a)，在圖中可以清晰地看到故障特征頻率fo(87.01 Hz附近)及其2倍頻(174 Hz附近)、3倍頻(262 Hz附近)、4倍頻(349 Hz附近)、5倍頻(436 Hz附近)、6倍頻(523 Hz附近)、7倍頻(610 Hz附近)瞬時(shí)頻率曲線，并且它們基本上凸顯在時(shí)頻平面上，因此，可以準(zhǔn)確地判斷故障類型為軸承外圈故障。

圖8 不同方法得到的時(shí)頻表示

實(shí)驗(yàn)結(jié)果表明，AHSST-QSTFT方法相較HSST具有更好的抗噪聲干擾能力，其最大程度地凸顯了故障特征信息，增強(qiáng)了時(shí)頻平面故障特征的辨識(shí)度。

4 結(jié)語(yǔ)

本文提出了一種基于QSTFT的自適應(yīng)高階同步壓縮變換算法，充分利用了QSTFT的窗口自適應(yīng)性和HSST的時(shí)頻壓縮特性。該方法提高了經(jīng)典STFT框架下振動(dòng)信號(hào)的窗口選擇適應(yīng)性，使重排后頻率估計(jì)最大程度地接近信號(hào)的真實(shí)瞬時(shí)頻率，提高了時(shí)頻平面的清晰度，進(jìn)一步增強(qiáng)了HSST方法在處理時(shí)變信號(hào)時(shí)的分辨能力。同時(shí)，數(shù)值仿真和軸承外圈故障實(shí)驗(yàn)表明，本文提出的方法能夠有效降低噪聲的干擾，并且在提取更為準(zhǔn)確的故障特征信息方面更具優(yōu)勢(shì)。