1、具有迭代结构的时频分布计算方法第 l1 卷第 3 期2006 年 6 月电路与系统J0URNAL0FCIRCUITSANDSYSTEMSVo1.11,No.3June,2006文章编号:10070249(2006)03 009605具有迭代结构的时频分布计算方法李滔,汤建龙,杨绍全(西安电子科技大学电了工程学院,陕西西安 7l0071)摘要 t 本文讨论了时频分布中的迭代算法的问题.通过选择特殊的计算窗口 ,时频分布的迭代运算形式能够有效地利用前面数据段的分析结果,从而避免了重复性的运算,使得计算效率得到提高.本文对原有的利用单边窗口的时频分布迭代算法的性能进行了分析,提出了采用对称型窗口的迭
2、代计算形式.与单边窗口相比,双边形式的计算窗口不但可以有效地提高时频表示精度,同时还能够更为准确的表示信号的瞬时频率.文章对各项理论分析结果提供了相应的仿真实验结果.关键词:迭代算法;时频分布;指数遗忘分布中图分类号:TN971.1 文献标识码:A1 导言作为对传统频率分析方法的扩展,时频分析方法通过时间和频率的二维函数反映出信号频率分量随时间的变化趋势,从而直观的刻画了非平稳信号的时变特征.虽然目前存在的各种时频表示方法的表现形式各异,但是通过考察可以发现,当前时频研究的重点主要集中在如何提高计算效率,消除交叉项以及改善时频分辨率等几个方向上.现有时频分布方法中的大多数都是基于傅立叶变换的,
3、这些时频分布形式通过使用滑动窗口对信号进行取样和加权,并对窗函数内部的信号进行频率分析.通常在计算不同时间点上的频率分布时,应当相应的移动窗口的位置,对各个时间点的邻域分别进行采样和加权.因此在这样的时频分布计算中,运算的复杂程度通常与窗口的宽度有关.简而言之,窗口越宽,则获得的时频分布聚集性越高,而相应的计算也越复杂.有鉴于此,人们总是希望在时频分布计算中引入更为高效的算法.其中一个较为实际的解决措施是采用迭代算法,使之能够在计算中有效的利用前面数据段的处理结果,从而避免重复性的计算.在文献【l3中提出了短时傅立叶变换的特殊形式 ,分析表明通过采用特殊的计算窗口可以有效地提高计算效率.因此在
4、讨论高效的时频分布计算形式时,可以把窗函数的选择作为一个出发点.文献【4】介绍了计算时间序列的时变傅立叶变换算法,当计算不同的时间点上的频率分布时,这种算法并不是简单的舍弃旧的数据,而是通过全极点计算窗口来不断的减少它们在计算中的影响.与通过 FFT 的实现算法相比,这种方法的计算量远小于后者.在文献【5,6】中是通过建立系统方程的方法来获得信号的傅立叶变换结果.通过迭代算法对傅氏系数的估计可以转化为状态方程的初始条件求解问题.而文献【7】提出的谱分析方法的特点在于:当新的数据出现时,频率分析所需要的计算数量与白相关计算中的时间滞后没有关系,因此采用这种迭代算法的谱估计方法被称为具有时移不变性
5、.本文讨论了时频分布的迭代运算形式.对能够利用迭代算法的时频分布形式进行了讨论.分析表明,与其他的时频分布计算形式相比,利用对称型窗口的时频分布能够更为准确的表示信号的瞬时频率的分布,同时也可以利用迭代算法使得运算效率得到提高.2 时频分布的迭代计算形式收稿日期:200312-26 修订日期:2004042基金项目 t 国防预研资助项日(41101030103)第 3 期李滔等:具有迭代结构的时频分布计算方法 97短时傅立叶变换通常也被称为时变傅立叶分析(timevaryingequencytransform,TVFT).设非平稳的时间序列为 x(m)(m=0,1,2),而 m)是时域计算窗
6、I51.当窗 I51 滑动时,可以表示为 w(nm),其中的参数表示窗 I51 的位置.计算窗 I51 对时间序列的加权结果为:(,)=一 m)x(m)m=0,1,2(1)如果设滑动窗 I51 的中心位置为:可以令原点的位置为 ,因此的时变的傅立叶变换为:TVFTz(co,一 ,)=+(71 一 f)】?P 一(2)这里代表采样的时间间隔.如果考虑(,)的镜像(,), 则有下面的关系:(,)=(一 m,)=m)x(nm)(3)与(2)式类似,(,)的傅氏变换为 :FAco,71)=(,71)?e= (一 m,71)?)?P(4)(,)与 TVFTz(co,一,)之间的关系可以用下面的公式表示:
7、(co,一,)=(-co,)-e(5)从上面的公式(4)可以推导出一类时域计算窗 I51,使得在(co,)的求解过程能够用到迭代算法.现在令:h(co,m)=m)?e-J因此公式(4)可以表示为 :(,71)= x(n 一)?h(oJ,)=x()?h(oJ,)(6)在上式中木表示的是卷积运算.因此可以获得公式(6)的 Z 变换形式:(co,Z)=X(z)?n(ro,z)(7)在式中(,z),(z),H(co,z)分别表示 (缈,),()与 h(oJ,)的 Z 变换.当时域计算窗 I51m)的 Z 变换具有如下的形式时:W(z)=01P=1,23-(8)rlye)由于在 W(z)中只存在极点,因
8、此这一类的窗 I51 可以称为全极点窗 I51,在公式中p 为分式的阶数.根据窗 I51 的 z 变换形式 ,可以求得其时域表示为:)=“f)从定义式(9)可以看到全极点窗口都是单边的,因此根据这些计算窗口得到的信号在某一个时间点上的频率分布仅仅依赖于信号在该时间点之前的历史表现.根据前面的分析可以从公式(9)定义的全极点窗口的 z 变换中推导出下面的结论:()南因此可以相应的得剑信号时频分布的 z 变换结果:,(11)从公式(11)可以得到利用全极点计算窗口的时频分布的迭代形式:Fx(co,71)=()一(一-e-jeJat),5(co,一)3 指数遗忘分布在具有迭代结构的时频分布中,最具有
9、代表性的是指数遗忘分布(exponentiallyforgettingtransform,EFT),该分布对应于公式中的特殊情况.正如其名称所表示的,EFT 使用单边指数窗口,98 电路与系统第 11 卷在 EFT 计算中信号的历史按照指数规律进行加权, 其作用相当于对历史进行逐渐的遗忘.因为在这种分布中,总是假定信号中距离需要分析的时问点较近的部分与其他部分相比,在频率分析中会起到更大的作用.离散形式的指数遗忘分布有下面的定义:X=()C-)P 州哪(12)在上式中0 被称为遗忘系数,该参数决定了单边指数窗口的形状.根据指数遗忘分布的定义式,可以得到如下的迭代结构:n+lX 川,:()?C-
10、)堋?Pn卜“=P ()?e-J 堋?m+川?P斛:P X+?P(13)=埘 =一公式(13)说明在 EFT 计算中 ,运算的复杂程度与所采用的单边指数窗口的长度无关.例如在时间点 n+1 处信号的频率分布可以通过前一个时间点 F/上的频率分布与信号在,z+1 处取值的线性组合来得到.因此与利用 FFT 算法的频率分析方法相比 ,EFT 算法具有更高的计算效率.在信号的时频分布中,人们总是希望时频分布的峰值能够相应的反映出瞬时频率的时间变化规律,这是目前广泛采用时频分布的瞬时频率估计方法的理论依据.因此从某种意义上来说,瞬时频率表示的准确程度也可以作为对时频分布进行评价的一个标准.通过分析可以
11、发现,尽管EFT 具有运算效率高的优点,它却不能够对信号的瞬时频率获得好的描述,EFT 的峰值分布仅仅能够无偏的表示平稳信号的瞬时频率.作为瞬时频率的估计方法,通过 EFT 的峰值选择进行 IF 估计的偏差与均方差(MSE)8可以分别表示为:bias:(14)E 缈 I21:一一(15).l6?SNR公式(14)表明 :利用 EFT 分布的峰值提取方法对 IF 的估计的偏差与两个参数有关,其一是信号的相位结构,其二是遗忘系数.当信号相位的二次微分不为零的时候,瞬时频率的估计值与真值相比将会产生漂移.而且相位的二次微分变化越剧烈,遗忘系数越小则漂移越大.因此限制了 EFT 的实际应用,因为在大多
12、数场合下信号的相位结构往往未知,所以瞬时频率估计的偏差也无法预料.上面的分析也说明,尽管全极点的计算窗口带来的迭代运算非常简单,但是由此产生的时频分布并不理想.以得到广泛应用的 WignerVille 分布为例,该分布对于线性调频信号(LFMsigna1)的瞬时频率是无偏的估计子,因此通常可以选取 LFM 信号作为时频分布表示性能的个评价标准.使用单边窗口的时频分布中,以 EFT 为例,其分布的峰值向瞬时频率逼近 ,但如果瞬时频率是时变的,则 EFT 表示的估计值只能够朝着瞬时频率的真值逼近 ,而永远不会到达瞬时频率.这是因为对频率分析而言,单边的计算窗口使得频率分析结果仅仅依赖于信号的历史信
13、息,因此估计的结果总是滞后于信号的瞬时频率变化.有鉴于此在时频分布中采用对称的计算窗口将是更为适当的选择.4 采用对称计算窗口的时频分布迭代形式本节将讨论具有对称计算窗口的时频分布的迭代形式.考虑对称的指数窗口:h(k)=ea,a0,该窗口的 Z 变换形式为:/(z):eat.z(16)方程(16)可以写为 :一l(z):P=za.ktCaZ.1Z+lZlCZz(e.一 C)z 一(C 一+P)z+1因此滑动窗口的傅立叶变换为:(17)(18)第 3 期李滔等:具有迭代结构的时频分布计算方法由此可以导出采用对称的指数窗口的时频分布迭代算法为:(,?)=(P 一+P).rFx(o),? 1).P
14、 一一 rFAo),?2)?P 一+(P 一 ea)?x(n1)-P 一(19)由于在这里使用的计算窗口是双边指数形的,因此可以利用类似于指数遗忘分布的命名方法将这种分布命名为双边指数遗忘分布(doublesidedexponentiallyforgettingtransform,DSEFT).利用 DSEFT的瞬时频率估计的偏差为:而相应的均方差为:bias-E】_(20)(21)公式(20),(21)表明利用双边指数遗忘分布的瞬时频率表示对于线性调频信号是无偏的,同时估计的均方差远小于指数遗忘分布.在前面所述的双边指数遗忘分布中使用的窗口长度是无限长的,这里还将讨论具有有限长度计算窗口时频
15、分布迭代算法.以余弦窗口为例,设窗口的宽度为 2N+1,则窗口可以表示为:()=COS(O)1)-NkN(22)其中=/2 .该窗口的 z 变换为 :(Z1:cos(N+1)O)1(zu+l-zu-,)(23)z一 2zCOS(O)1)+1而相应的可以得到对应于迭代计算的傅立叶变换为:=(24,?)=2X(O),?一 1)?cos1)?P 一础一,l 一 2)?ej2 础 r,c,一cos(N+1)1.x(nN1)?e-J+x(n+N1)?e 掰一】公式(25)表明 ,利用迭代算法的时频分布计算窗口不局限于全极点的计算窗口,也可以扩展到其他类型的窗口;此外在 EFT 与 DSEFT 变换中使用
16、了具有无限长度的窗口 ,但是对于有限宽度的窗口入余弦窗口,仍然可以使用迭代算法,其计算形式并不唯一.5 仿真结果为了验证提出的算法,这里使用了两种信号,第一种是相位的二次微分线性增长的信号“,而另一种信号的相位的三次微分线性增加 ej 枷.图 l(a)说明对于第一种信号利用 EFT 得到的瞬时频率估计的偏差是线性增加的,而利用 DSEFT 算法的瞬时频霉絮.簪 g嚣通遗忘子图 2 瞬时频率估计的 MSE与遗忘系数的关系率估计的偏差为零.图 l(b)婆 86萋.时间(a)信号 Pm.时间(b)信 e/2图 1 两种信号的瞬时频率估计偏差表明对于第二种信号,利用 DSEFT 的瞬时频率估计的偏差是
17、线性增加的,而利用 EFT 的瞬时频率估计的偏差则为二次函数, 而且偏差比前者大很多,因此与 EFT 相比 DSEFT 是更有效的频率估计方法.图 2 表示了在 EFT 与 DSEFT 中瞬时频率估计 MSE(均方差)与遗忘系数之间的关系,信噪比为 10dB.图中从上至下的三条曲线分别表示 EFT,DSEFT 与理论分析的结果的自然对数,从图中可以看到与 EFT 相比,DSEFT 的瞬时频率估计的均方差仅为前者的 1/8,这一点与理论分析的结果相吻合.图 3 用来比较分别采用全极点计算窗口与对称计算窗口的时频分布对正弦调频信号瞬时频率的估计性能,图中 3 条曲线分别代表信号的瞬时频0.450.4船 0.350.30u
