四種時頻分析方法的頻率分辨率研究

摘 要：時頻分析是分析時變譜的有力工具，其頻率分辨率是值得研究的關鍵問題之一。通過對基于短時傅里葉變換、小波變換、Choi-Williams分布和Hilbert-Huang變換的四種時頻分析方法的頻率分辨率的實驗比較， 說明Hilbert-Huang變換法具有最高的頻率分辨率，其次是Choi-Williams分布，小波緊隨其后，最差是短時傅里葉變換。



關鍵詞：短時傅里葉變換；小波變換；Choi-Williams分布；Hilbert-Huang變換；時頻分析；頻率分辨率

中圖分類號：TN911.72 文獻標識碼：A

文章編號：1004-373X(2008)11-037-03

Experimental Research on Frequency Resolution of Time-Frequency

Distribution Using Four Methods



YANG Shenshen1，LI Ming1，GU Xuekang2

(1.School of Information Science and Technology，East China Normal University，Shanghai，200241，China;

2.China Ship Scientific Research Center，Wuxi，214082，China)



Abstract：Time-Frequency (TF) is a powerful tool for the analysis of time-varying spectrum.One of key issues in this regard is frequency resolution.Experimentally compared the frequency resolution of four types of TF analysis methods，namely，Short Time Fourier Transform (STFT)，Wavelet Transform (WT)，Choi-William Distribution (CWD) and Hilbert-Huang Transform (HHT).The results suggest that HHT has the highest frequency resolution，CWD is better than WT and STFT is the worst.



Keywords：short-time Fourier transform;wavelet transform;Choi-William distribution;Hilbert-Huang transform;time-frequency analysis;frequency resolution

時頻分析(TF)的基本思想是設計時間和頻率的聯合函數，用它同時描述信號在不同時間和頻率的能量密度或強度是分析非平穩信號的有力工具。信號特征的提取與TF分布的頻率分辨率有很大關系。在頻率分辨率的研究中，取頻率分量較少且頻率之差盡量小的正弦疊加信號是較簡便的實驗做法^[1]。本文用4個頻率分量(f1，f2，f3，f4)的正弦疊加信號作為待分析信號，把STFT，WT，CWD和HHT四種時頻分析拿出來做比較，展現四者在頻率分辨率上的差異。

1 算法原理

1.1 Short Time Fourier變換(STFT)

為了研究信號在時間t的特性，一個自然的想法是加強在那個時間的信號，壓縮在其他時間的信號。這可以通過用中心在t的窗函數h(t)乘信號來實現，產生的被改變的信號為:



st(τ)=s(τ)h(τ－t)

(1)



這里，改變的信號是兩個時間即t(固定時間)和τ(執行時間)的函數。窗函數決定留下的信號圍繞著時間t大體上不變，而在遠離t的地方信號被壓縮了許多倍。也就是:



st(τ)s(τ) 對于τ接近t

0 對于τ遠離t

(2)



因為改變的信號加強了圍繞著時間t的信號，所以，它的Fourier變換反映的是在時刻t的一個小的時間間隔內的頻率分布，即:



St(ω)=12π∫e－jωτst(τ)dτ

=12π∫e－jωτst(τ)h(τ－t)dτ

(3)



因此，在時間t的能量密度頻譜是:



Psp(t，ω)=|St(ω)|2

=12π∫e－jωτ s(τ)h(τ－t)2

(4)



對于每一個不同的時間，都可以得到一個不同頻譜，這些頻譜的總體就構成了一個時頻分布，即頻譜圖^[2]。

窗函數對STFT的性能有很大的影響：為了提高時間分辨率，需要選擇的窗函數h(t)盡可能短；但為了頻率分辨率，則要求選擇的窗函數h(t)時間寬度盡可能長，因此與時間分辨率的提高相矛盾。另一方面，對于要分析的非平穩信號而言，也許某一小時段上是以高頻信息為主，我們希望用短時間窗進行分析；而在某一場時間段上是一些低頻信息，我們希望用一個長時間窗進行分析。因此對一個時變的非平穩信號，利用STFT難以找到一個合適的時間窗口來適應于不同的時間段，這也是它的最大不足之處^[3]。

1.2 Wavelet變換(WT)

小波分析^[4，5]本質上是一種時間尺度分析。它克服了STFT的時頻窗不能改變形狀的缺陷，而可以有效地聚焦信號的瞬時結構。

小波的一般定義如下:



W(a，b;X，Ψ)=|a|－12∫∞－∞X(t)Ψ*(t－ba)dt

(5)



其中，Ψ*(#8226;)為小波基函數，a為尺度因子，b表示與坐標原點的平移。方程(5)的一個直觀的物理解釋是：W(a，b;X，Ψ)是X在時間t = b，尺度為a的能量。

在本文，我們使用Morlet小波計算信號x(t)的時頻分布。其計算表達式如下：



SCx(t，a;h)=|Tx(t，a;h)|2

=1|a|∫∞－∞x(s)h*s－tads2

(6)



其中a=v0v，v0是母小波函數h(t)的中心頻率^[6]。

1.3 Chio-Williamas分布(CWD)

20世紀60年代中期，Cohen發現眾多的時頻分布只是Wigner-Ville分布(WVD)的變形，可以用統一的形式表示，習慣稱之為Cohen類時頻分布，那么其表達式為^[2，7]：



Wx(t，ω)=∫∞－∞∫∞－∞∫∞－∞x(u+τ2)x(u－τ2)

φ(θ，τ)e－j(tτ+τω－uθ)dudθdτ

(7)



式中x(t)為解析信號，φ(θ，τ)為核函數。當取φ(θ，τ)=1時，即表示Wigner-Ville分布：



Wx(t，ω)=∫∞－∞x(t+τ2)x(t－τ2)e－jωπdτ

(8)



Choi和Williams^[7]提出了時頻分析核，其表現形式為：



φ(θ，τ)=e－θ2τ2/σ

(9)



Choi和Williams核函數與時間和頻率無關，是時延θ和頻延τ的函數，具有時頻移不變性。σ為衰減系數，它與交叉項的幅值成比例關系。

相應的分布為:



Wx(t，ω)=∫∞－∞∫∞－∞∫∞－∞x(u+τ2)x(u+τ2)#8226;

e－θ2τ2/σ2e－j(tτ+τω－ωθ)dudθdτ

(10)



此分布稱為Choi-Williams分布。注意到當σ→∞時，我們便得到WVD，因為此時核趨近于1。相反地，σ越小，相干項的衰減就越大。此分布具有以下數學性質^[6]：

P1 y(t)=x(t－t0)Wy(t，ω)=Wx(t－t0，ω)

y(t)=x(t)ej2πω0tWy(t，ω)=Wx(t，ω－ω0)



時移和頻移不變性:

P2 y(t)=kx(kt)

k>0Wy(t，ω)=Wx(kt，ωk) 時頻伸縮性。

P3 fx(t)=∫+∞－∞ωWxa(t，ω)dω∫+∞－∞Wxa(t，ω)dω(xa是信號x對應的解析信號)瞬時頻率。

P4 tx(ω)=∫+∞－∞tWxa(t，ω)dω∫+∞－∞Wxa(t，ω)dω群延遲。

我們在Matlab軟件中用時頻工具箱中tfrcw.m來計算。

1.4 Hilbert-Huang變換

非平穩信號的主要特征是其頻率是時間的函數。但研究表明：STFT，WT，CWD這些時頻分析法的最終理論依據都是傅里葉變換，因而不可避免地會受到傅里葉變換分析非平穩信號的缺陷的影響，如出現虛假頻率和多余信號分量等。N.E.Huang等人提出的HHT從根本上突破了傅里葉變換理論的限制，首次建立了一種基于瞬時頻率的、較系統和完整的信號分析方法^[8]。

文獻［9］提出了內蘊模式函數的概念和EMD篩選算法。固有模態函數(IMF)必須滿足兩個條件：信號極值點的數量與過零點的數量必須相等，或最多相差一個；在任一時間點上，信號極大值定義的上包絡和極小值定義的下包絡的局部均值為零。

EMD算法：對任一實信號s(t)，首先確定出s(t)上的所有極大值點和極小值點，然后將所有的極大值點和所有的極小值點分別用一條光滑的曲線聯接起來，使兩條曲線間包含所有的信號數據。將這兩條曲線分別作為s(t)的上下包絡線，計算出它們的平均值曲線m(t)，s(t)與m(t)的差記作h(t)，則^[8]:



s(t)-m(t)=h(t)



把h(t)作為原信號，重復上面介紹的步驟，直到當h(t)滿足一定的條件時，記:



c1(t)=h(t)



將c1(t)視為一個IMF，再作:



s(t)-c1(t)=r(t)



將r(t)視為新的s(t)，重復以上過程，一次得第二個IMF c2(t)，第三個IMF c3(t)，…。當cn(t)或r(t)滿足給定的終止條件時，分解過程終止，得分解式:



s(t)=∑Nn=1cn(t)+rn(t)

(11)



其中，r(t)稱為殘余函數，代表信號的平均趨勢。

對式(11)中的每個IMF分別作Hilbert變換后得:



s(t)=Re∑ni=1ai(t)ejΦi(t)=Re∑ni=1ai(t)ej∫ωi(t)dt

(12)



這里省略了殘余函數r(t)，Re表示取實部。稱展開式(12)為Hilbert幅值譜，簡稱Hilbert譜，記作:



H(ω，t)=Re∑ni=1ai(t)ej∫ωi(t)dt

(13)



2 實驗結果

我們選用四個頻率的余弦信號進行分析，如下：



x(t)=［cos(2πf1t)+cos(2πf2t)］［u(t)－u(t－100)］+

［cos(2πf3t)+cos(2πf4t)］［u(t－200)－u(t－300)］



其中f1，f2，f3，f4皆是頻率。令f1=0.05，f2=0.20，f3=0.25，f4=0.45，則:Δf1=0.15，Δf2=0.20。由此可見，兩組頻率之間的差是很微小的，通過觀察4種時頻方法的等高圖和立體圖，可以比較出它們頻率分辨率的優劣。

首先來看HHT中各IMF分量及其余量的Hilbert譜：

圖1 信號x(t)HHT中各IMF分量及其余量的Hilbert譜

STFT的窗函數選的是長度為11的hamming窗，進行WT分析時，利用時頻工具箱中的函數［tfr，t，f，wt］=tfrscalo(x，t，wave，fmin，fmax，N，trace)，其中wave取3；進行CWD分析時，利用時頻工具箱中的函數［tfr，t，f］=tfrcw(x，t，N，G，H，signal，trace)，其中時間窗G取的是長度為15的hamming窗，時間窗h取的是長度為11的hamming窗，衰減系數σ即函數中的signal取為1。

由圖2可以看到STFT的頻率分辨率很糟糕，兩組頻率交織在一起；WT的情況好一些，第一組頻率f1，f2尚能分辨，但第二組頻率f3，f4就分不清楚了；CWD顯然比上面兩種情況要改善了許多，兩組頻率很清晰地顯示出來，只是Cohen類的分布都屬于二次型的，因而產生的相干項不可避免的干擾了時頻表示的可讀性；HHT頻率分辨率最優。

圖2 信號x(t)的時頻等高圖

由立體圖圖3可以進一步證實上面的結論。另一方面我們研究的信號幅值都是相等的，而WT這種方法幅值的變化過于凸顯，與事實不符，其他三種方法都表現良好。

圖3 信號的立體圖比較

3 結 語

本文對基于STFT，WT，CWD和HHT的四種時頻分布的頻率分辨率進行了實驗比較。結果表明，STFT的頻率分辨率最低，WT較前者有了明顯改善，CWD表現更好， HHT的頻率分辨率最優。在這里只對頻率分辨率做了實驗比較，時間分辨率也是一個很好的研究方向，但本文還未涉及，作者將在以后的工作中繼續對這方面的工作進行深入研究。

參 考 文 獻

［1］Li Ming.Experimental Research for Processing of Choi-Williams Distribution and Bessel Distribution［A］.IEEE Conf.Proc.，1st Int.Conf.on Information，Communication Signal Processing，Singapore，1997，3:1 593-1 597.



［2］［美］ L.科恩.時-頻分析：理論與應用［M］.白居憲，譯.西安:西安交通大學出版社，1998.

［3］葛哲學，陳仲生.Matlab時頻分析技術及其應用［M］.北京:人民郵電出版社，2006.



［4］Mallat S.A Theory for Multi-Resolution Signal Decomposition Wavelet Representation ［J］.IEEE Trans.on Pattern Analysis and Machine Intell，1989，11(7):674-693.



［5］Coilman R R，Wickerhauser M W.Entropy-based Algorithms for Basis-Selection ［J］.IEEE Tran.Information Theory，1992，38(2):713-718.



［6］胡昌華，周濤，夏啟兵，等.基于MATLAB的系統分析與設計——時頻分析［M］.西安:西安電子科技大學出版社，2001.



［7］劉希強，沈萍，山長侖，等.數字化地震波形資料的時頻分析方法及應用［J］.西北地震學報，2004，26(2):118-125.



［8］鐘佑明，秦樹人.HHT的理論依據探討——Hilbert變換的局部乘積定理［J］.振動與沖擊.2006，25(2):12-19.



［9］Huang N E，Zheng Shen，Long S R，et al.The Empirical Mode Decomposition and the Hilbert Spectrum for Nonlinear Non-Stationary Time Series Analysis ［J］.Proc RSOC London，1998，454:903-995.



［10］Li Ming，Jia W，Zhao W.Analysis of Pattern Similarity between Choi-Williams Kernel and Born-Jordan Kernel［A］.Proc.of IEEE ICICS99，Singapore，1999.



作者簡介 楊哂哂 女，1984年出生，碩士研究生，從事信號處理研究。

李 明 男，1955年出生，PHD，華東師范大學教授，博士生導師。主要研究信號處理，計算機科學等。

顧學康 男，博士，中國船舶科學研究中心研究員，博士生導師。

注：本文中所涉及到的圖表、注解、公式等內容請以PDF格式閱讀原文。