非平稳信号的时频分析方法

非平稳信号的时频分析方法

0非平稳信号的一般分析方法非稳定信号分析是当前信号分析领域的一个重要课题。早期对非平稳信号几乎没有什么正规的分析方法。而发展最早、最成熟的信号分析方法——傅立叶变换是全局变换,用其分析平稳信号很有效,但用于分析非平稳信号缺乏物理意义。非平稳信号是统计量(相关函数、功率谱等)随时间变化的信号,因此时频分析方法是分析这类信号的有效手段。已有的时频分析方法很多,目前较典型的有短时Fourier变换(STFT),Wigner_Ville分布,小波变换(WT)等。尽管这些方法对非平稳信号的分析做出了较大的贡献,在工程实际中也获得了较广泛的应用,但它们大都还是以傅立叶变换为其最终的理论依据。傅立叶变换理论中表征信号交变的基本量是与时间无关的频率,基本时域信号是平稳的简谐波信号。这些概念是全局性的,因而用它们分析非平稳信号容易产生虚假信号和假频等矛盾现象。对非平稳信号比较直观的分析方法是使用具有局域性的基本量和基本函数。瞬时频率是容易想到的具有局域性的基本量,也是很早就已提出的概念。瞬时频率的比较直观的定义是解析信号相位的导数,但以往这一定义会产生一些佯缪的结果,导致基于瞬时频率的时频分析方法和理论始终未真正建立和发展起来。1996年,美籍华人NordenE.Huang等人在对瞬时频率的概念进行了深入研究后,创立了Hilbert_Huang变换(HHT)的新方法。这一方法创造性地提出了固有模态信号的新概念以及将任意信号分解为固有模态信号组成的方法——经验模态分解法,从而赋予了瞬时频率合理的定义、物理意义和求法,初步建立了以瞬时频率为表征信号交变的基本量,以固有模态信号为基本时域信号的新时频分析方法体系。这一方法体系从根本上摆脱了傅立叶变换理论的束缚,能很好地解释以往将瞬时频率定义为解析信号相位的导数时容易产生的一些所谓“悖论”,在实际应用中也业已表现出了一些独特的优点。但是这一新的方法还处在发展阶段,在建立严密的理论和方法的完善方面还有许多事要做。本文主要对这一方法做一些理论上的探索和完善。1经验模态分解法Hilbert_Huang变换首先假设:任一信号都是由若干固有模态信号(IntrinsicModeSignal,简称IMS)或固有模态函数(IntrinsicModeFunction,简称TMF)组成的,任何时候,一个信号都可以包含许多固有模态信号,如果固有模态信号之间相互重叠,便形成复合信号。其中固有模态信号是满足以下两个条件的信号:(1)整个数据中,零点数与极点数相等或至多相差1;(2)信号上任意一点,由局部极大值点确定的包络线和由局部极小值点确定的包络线的均值均为0,即信号关于时间轴局部对称。在Hilert_Huang变换中表征信号交变的基本量不是频率,而是瞬时频率(InstantaneousFrequency,简记IF)。瞬时频率可以通过Hilbert变换获得,即先对信号s(t)作Hilbert变换得解析信号z(t)=s(t)+jH[s(t)]=a(t)ejΦ(t)(1)幅值函数a(t)=√s2(t)+Η2[s(t)](2)a(t)=s2(t)+H2[s(t)]−−−−−−−−−−−−√(2)和相位函数Φ(t)=arctanΗ[s(t)]s(t)(3)Φ(t)=arctanH[s(t)]s(t)(3)再对相位函数求导即得瞬时角频率和瞬时频率ω(t)=dΦ(t)dt(4)ω(t)=dΦ(t)dt(4)f(t)=12πdΦ(t)dt(5)f(t)=12πdΦ(t)dt(5)但NordenE.Huang等人分析认为,通过以上过程求得的瞬时频率只对固有模态信号才具物理意义。实际信号常常是复合信号,因此对实际信号进行时频分析时,需要先将信号分解成IMS的和。为此Hilbert_Huang变换的创立者又提出了一种经验模态分解方法(theEmpiricalModeDecomposition,简称EMD),这是一种经验筛法,其过程介绍如下:对任一信号s(t),首先确定出s(t)上的所有极值点,然后将所有极大值点和所有极小值点分别用一条曲线连接起来,使两条曲线间包含所有的信号数据。将这两条曲线分别作为s(t)的上、下包络线。若上、下包络线的平均值记作m,s(t)与m的差记作h,则s(t)-m=h(6)将h视为新的s(t),重复以上操作,直到当h满足一定的条件(如h变化足够小)时,记c1=h(7)将c1视为一个IMF,再作s(t)-c1=r(8)将r视为新的s(t),重复以上过程,依次得第二个IMFc2,第三个IMFc3,…。当cn或r满足给定的终止条件(如分解出的IMF或残余函数r足够小或r成为单调函数)时,筛选过程终止,得分解式s(t)=n∑t=1ci+r(9)其中,r称为残余函数,代表信号的平均趋势。对式(9)中的每个IMF分别作Hilbert变换后得s(t)=Ren∑i=1ai(t)ejΦi(t)=Ren∑i=1ai(t)ej∫ωi(t)dt(10)这里省略了残余函数r,Re表示取实部。称展开式(10)为Hilbert幅值谱,简称Hilbert谱,记作H(ω,t)=Ren∑i=1ai(t)ej∫ωi(t)dt(11)进一步可以定义边际谱h(ω)=∫∞-∞H(ω,t)dt(12)以上的EMD和与之相应的Hilbert谱信号分析方法统称为Hilbert_Huang变换。展开式(10)中,每个组成成分的幅值和相位是随时间可变的,而同样信号s(t)的Fourier变换展开式为s(t)=Re∞∑i=1aiejωit(13)其中ai,ωi为常数。因此Hilbert_Huang变换可以看作是Fourier变换的一般化,其创新主要体现在两个方面:创造性地提出了固有模态信号(IMS)的概念;创造性地提出了经验模态分解方法(EMD)。目前Hilbert_Huang变换已在海洋、地震、生物工程、桥梁健康监测等实际应用中显示出了独特的优点。但这一方法提出的时间还不长,在理论的建立和方法的完善等方面还有待继续深入研究。例如,固有模态信号的描述性定义有待提出具体的数学模型;通过极值点拟合固有模态信号包络线的合理性有待进一步论证;将任意信号分解成IMS组合的EMD方法有待理论论证和进一步完善或寻找新的分解方法。另外,在算法中,边界处理也是需要解决的一个重要问题。以下介绍作者在Hilbert_Huang变换的研究中所做的一些理论探索工作。2角速度t和st转换规则NordenE.Huang等人创造性地提出了固有模态信号的概念,但对这一概念的定义迄今仍是描述性的。为建立固有模态信号的数学模型,本文作者作了如下探索。由式(1)—(3)容易得到,IMS可用以下数学形式表示c(t)=a(t)cosΦ(t)(14)但并非任何形如a(t)cosΦ(t)的信号都是IMS。本文指出,形如a(t)cosΦ(t)的信号成为IMS必须满足一定的数学条件。以下用振动模型对此进行说明。在振动理论中,频率是一个表征信号交变的基本变量,它是指单位时间内物体往复振动(或位移交变)的次数。如图1(a)所示,当物体以角速度ω沿半径为a0的圆周运动时,物体在直径上的投影P的运动是一简谐振动,其位移为s(t)=a0cosϑ(t)=a0cosωt(15)当时间t经过2π/ω个单位以后,P往复振动一次,于是P的频率f=1/(2π/ω)=ω/2π(Hz)。显然,当角速度ω越大,P往复振动越快,频率f越大。反之,当f越大时,也必然要求ω越大。因此“频率”本质上可以理解为一个表示往复振动快慢的物理量,即瞬时频率。以上振动中,角速度ω和半径a0均为常数,这样形成的信号是平稳信号,其瞬时频率处处相等。但实际中物体绕原点运动的半径往往不为常数,运动的角速度也不均匀,如图1(b)所示。于是投影P的表达式变为s(t)=a(t)cosΦ(t)(16)这样形成的信号是非平稳信号,其中决定投影P振动快慢的量仍然是角速度,因此式(5)表示了信号s(t)=a(t)cosΦ(t)的瞬时频率,即信号交变快慢的局部特性。显然这时的瞬时频率是时变的。式(16)与式(5)体现了非平稳信号随时间变化的根本特征,因而是分析非平稳信号的基础。但由图1可以直观得出,欲使式(5)能始终表示信号s(t)=a(t)cosΦ(t)的交变特征,调幅信号a(t)相对于被调制信号cosΦ(t)应为缓变信号,其取值不应影响cosΦ(t)的单调性,即s(t)与cosΦ(t)应具有同向单调性,亦即s′(t)(cosΦ(t))′≥0(17)否则s(t)的交变性将与cosΦ(t)的交变性不一致,这样就不能再用式(5)定义s(t)的瞬时频率,即不能再用cosΦ(t)的交变性表征s(t)的交变性。如图2所示,(a)为信号cosΦ(t)的振动模型,(b)和(c)均为经信号a(t)调制后的信号s(t)=a(t)cosΦ(t)的振动模型,图2表示,当a(t)为满足式(17)的缓变信号时,s(t)与cosΦ(t)具有同向单调性,从而s(t)与cosΦ(t)的交变性一致,当a(t)为不满足式(17)的非缓变信号时,s(t)与cosΦ(t)不具有同向单调性,从而s(t)与cosΦ(t)的交变性不一致,因此,欲用式(5)能计算信号s(t)=a(t)cosΦ(t)的瞬时频率,s(t)与cosΦ(t)应满足式(17)的关系式。式(16)可用式(14)代替,因而式(17)给出了固有模态信号应满足的一个重要数学条件。称式(14)中的cosΦ(t)为固有模态信号的本征信号。由条件(17)可以导出一个重要结论:固有模态信号与其本征信号具有相同的极点和零点。对固有模态信号,其瞬时频率显然应是非负的,因此由式(5)得Φ′(t)≥0,而调制信号a(t)>0也是必须的,由此本文给出固有模态信号的数学模型:当式(14)满足条件(1)a(t)>0;(2)Φ′(t)≥0,即Φ(t)单调上升;(3)s(t)′(cosΦ(t))′≥0,即c(t)与cosΦ(t)具有同向单调性时,称信号c(t)为固有模态信号。以后称条件(3)为固有模态信号的本征条件。3不同频率的at对信号相干合成的影响基于以上提出的固有模态信号的数学模型,可以研究固有模态信号局部对称性要求的必要性和用极值点拟合固有模态信号包络线的合理性。首先研究是否存在不为常值的正值函数a(t)满足固有模态信号的本征条件。考察信号C(t)=a(t)cost(18)易得C′(t)(cost)′=a(t)sin2t-12a′(t)sin2t(19)因为a(t)sin2t≥0‚12sin2t是有界的,因此,对∀t,存在足够小但不一定为零的a′(t)使得C′(t)(cost)′≥0。即对式(19)存在a(t)≠常数,使C′(t)(cost)′≥0(20)虽然欲使式(18)满足条件式(20),不要求a(t)为常数,但a′(t)仍需受到一定限制。为研究方便,不妨令a(t)=cos(t/d)+4,则a′(t)=-1dsin(t/d),显然|a′(t)|越小,即a(t)频率(或瞬时频率)12πd越低,式(20)越易得到满足。另一方面,图3示出了不同频率的a(t)对信号C(t)=a(t)cost对称性的影响,其中图3(c)、(f)、(i)分别是对图3(b)、(e)、(h)圆圈内图形的放大,以观察C(t)与其Hilbert变换包络线的交点和相邻极值点间的位置差。由图3(a)、(d)、(g)可以看出:随着d的增大,信号C(t)=a(t)cos(t)=[cos(t/d)+4]cost关于时间轴的对称性越来越明显,同时,由图3(b)、(c)、(e)、(f)、(h)、(i)可知,由Hilbert变换求得的C(t)的包络线(即a(t))与C(t)的交点越来越接近。事实上,由C(t)=a(t)cosΦ(t)知,当cosΦ(t)=1时,C(t)=a(t),当cosΦ(t)=-1时,C(t)=-a(t),因此C(t)的局部对称性主要体现在cosΦ(t)的极值点时刻,C(t)的上包络线可以通过其上使cosΦ(t)=1的时刻对应的点拟合得出,同理C(t)的下包络线可以通过其上使cosΦ(t)=-1的时刻对应的点拟合得出(根据上一节的结论,若C(t)满足条件固有模态信号的本征条件,则cosΦ(t)=1时刻对应的点正是C(t)的极大值点,cosΦ(t)=-1时刻对应的点正是C(t)的极小值点)。显然,当a(t)相对于cosΦ(t)的频率较小时,这种拟合性较好。虽然由图3(b)、(e)、(h)的Hilbert变换包络线观察,似乎也可以认为这三个信号关于时间轴呈对称性,但当d=2.0[图3(b)]时C′(t)(cost)′=-[(cost2+4)cost]′sint=(12sint2cost+(cost2+4)sint)sint=12sint2[2cos2t2-1+4(cost2+4)cost2]sint=sin2t2[2cos2t2-1+4(cost2+4)cost2]cost2=sin2t2(6cos2t2+16cost2-1)cost2令cost2=0.04,则6cos2t2+16cost2-1=-0.350,于是由式(21)得C′(t)(cost)′<0,因而式(20)不成立。而由式(19),当d越大时,(20)越易成立。综上所述,固有模态信号的局部对称性要求一方面使得其包络线a(t)与其本身的交点较接近,因而可以通过固有模态信号的极值点拟合其包络线,另一方面使得其满足固有模态信号的本征条件,再一方面也保证用Hilbert变换能求出其相位函数,进而求出瞬时频率,因为存在下述结论:当A(f)=F{a(t)}完全位于区间|f|<f0内,而F{cosΦ(t)}位于区间|f|<f0之外时,H[a(t)cosΦ(t)]=a(t)sinΦ(t)。其中f0的含义是使[-f0,f0]包含a(t)频率范围的一个数,H=[a(t)cosΦ(t)]是a(t)cosΦ(t)的Hilbert变换。由此论证了固有模态信号局部对称性要求和用极值点拟合固有模态信号包络线的合理性。4经验模态分解emd的基本原理设c(t)=a(t)cosΦ(t)为一IMF,则cosΦ(t)=c(t)a(t)(22)-Φ′(t)sinΦ(t)=(c(t)a(t))´(23)[Φ′(t)]2[1-cos2Φ(t)]=[(c(t)a(t))´]2(24)[Φ′(t)]2(1-(c(t)a(t))2)=[(c(t)a(t))´]2(25)Φ′(t)=√[(c(t)a(t))´]2/(1-(c(t)a(t))2)(26)由以上式子可以说明经验模态分解(EMD)的理论依据。因为,由前面的讨论,式(26)中的a(t)应为缓变信号,因此[c(t)/a(t)]′主要由c′(t)决定,即固有模态信号的瞬时频率与该处信号的导数有关,且与|c′(t)|成近似正比关系。对任意信号s(t),令s(t)=c(t)+r(t)(27)其中c(t)为IMF。两边求导得s′(t)=c′(t)+r′(t)(28)如果对几乎所有的t,|c′(t)|相对|r′(t)|较大,则|s′(t)|主要由|c′(t)|决定,即s(t)的极值点主要由c(t)的极值点决定(根据极值点的判定定理)。因此,通过信号s(t)的极值点拟合其包络线,求出包络线均值m(t)(理想情况下应为r(t)),即可求得包含c(t)几乎所有信息的信号s(t)-m(t)(这一过程称为筛选)。实际中,m(t)≠r(t),即s(t)-m(t)≠c(t)。因此还需继续对s(t)