1、数字信号处理数字信号处理引言引言对各态遍历随机信号对各态遍历随机信号 X(n),自相关函数和功率谱密度均可用时,自相关函数和功率谱密度均可用时间平均来定义:间平均来定义:维纳维纳-辛钦定理:辛钦定理:22()1()lim()lim2121jNjj nxNNnNX eP eEx n eENN1()()()lim()()()21NXxMnNrmE X n X nmx n x nmr mN()()jj mxxmP er m e第十三章第十三章 经典功率谱估计经典功率谱估计1.周期图法(直接法)周期图法(直接法)2.间接法间接法3.直接法和间接法的关系直接法和间接法的关系4.直接法和间接法估计的质量、
2、直接法和间接法估计的质量、5.直接法的改进直接法的改进6.经典功率谱估计总结经典功率谱估计总结7.短时傅里叶变换短时傅里叶变换1.周期图法(直接法)周期图法(直接法)周期图法是把随机信号周期图法是把随机信号 X(n)的的 N 点观察数据点观察数据 xN(n)视为一能量视为一能量有限信号,直接取有限信号,直接取 xN(n)的的 DTFT 得到得到 XN(ej),然后再取其幅,然后再取其幅值的平方,并除以值的平方,并除以 N,作为对真实功率谱,作为对真实功率谱 P(ej)的估计:的估计:表示用周期图法估计出的功率谱。因为功率谱密度直表示用周期图法估计出的功率谱。因为功率谱密度直接由傅里叶变换得到,
3、所以周期图法又称直接法。自从接由傅里叶变换得到,所以周期图法又称直接法。自从 1965 年年 FFT 出现后,该方法就成了谱估计中的一个常用方法。将出现后,该方法就成了谱估计中的一个常用方法。将 在单位圆上等间隔取值,得在单位圆上等间隔取值,得由于由于 XN(k)可以用可以用 FFT 快速计算,所以可以方便地求出快速计算,所以可以方便地求出 。21()()jjPERNPeXeN()jPERPe21()()PERNPkXkN()jPERPe1.周期图法(直接法)周期图法(直接法)比较以下两种计算方法:比较以下两种计算方法:易知,直接法包含了下述假设及步骤:易知,直接法包含了下述假设及步骤:把平稳
4、随机信号把平稳随机信号 X(n)视为各态遍历的,用其一个样本视为各态遍历的,用其一个样本 x(n)来来代替代替 X(n),并且仅利用,并且仅利用 x(n)的的 N 个观察值个观察值 xN(n)来估计功率谱来估计功率谱P(ej)。从记录到一个连续信号从记录到一个连续信号 x(t)到估计出到估计出 ,还包括了对,还包括了对 x(t)的离散化、必要的预处理(如除去均值和趋势项、滤波等)。的离散化、必要的预处理(如除去均值和趋势项、滤波等)。22()1()lim()lim2121jNjj nxNNnNX eP eEx n eENN21()()jjPERNPeXeN()PERPk1.周期图法(直接法)周
5、期图法(直接法)一个实际的例子(一个实际的例子(fs=250Hz):):2.间接法间接法间接法的理论基础是维纳间接法的理论基础是维纳-辛钦定理。辛钦定理。1958年年Blackman和和Tukey给出了这一方法的具体实现,即先由给出了这一方法的具体实现,即先由 xN(n)估计出自相关函数,估计出自相关函数,然后求自相关函数的傅里叶变换得到的功率谱,记之为然后求自相关函数的傅里叶变换得到的功率谱,记之为 ,并以此作为对并以此作为对 P(ej)的估计,即的估计,即因为这种方法求出的功率谱是通过自相关函数间接得到的,所因为这种方法求出的功率谱是通过自相关函数间接得到的,所以称为间接法,又称自相关法或
6、以称为间接法,又称自相关法或BT法。当法。当 M 较小时,上式计较小时,上式计算量不是很大,因此该方法是算量不是很大,因此该方法是 FFT 问世之前常用的谱估计方法。问世之前常用的谱估计方法。与维纳与维纳-辛钦定理相比较:辛钦定理相比较:()jBTPe()(),1Mjj mBTmMPer m eMN()()jj mxxmP er m e2.间接法间接法如果如果 X(n)是各态遍历随机信号,是各态遍历随机信号,x(n)是其一个样本函数,则自是其一个样本函数,则自相关函数可定义如下:相关函数可定义如下:实际中的信号大多是因果信号,所以上式可以表示为:实际中的信号大多是因果信号,所以上式可以表示为:
7、本章所涉及的都是自相关函数,因此将本章所涉及的都是自相关函数,因此将 rx(m)简写为简写为 r(m)。如。如果观察值的个数为有限值,则求果观察值的个数为有限值,则求 r(m)的一种方法为:的一种方法为:由于由于x(n)只有只有N个观察值,因此对于每一个固定的延迟个观察值,因此对于每一个固定的延迟 m,可以,可以1()lim()()21NxNnNr mx n x nmN101()lim()()NNnr mx n x nmN101()()()NNNnr mxn xnmN2.间接法间接法利用的数据只有利用的数据只有 N-1-|m|个,且在个,且在 0N-1 的范围内,的范围内,xN(n)=x(n)
8、,所以实际计算时,上式变为:所以实际计算时,上式变为:的长度为的长度为 2N-1,它是以,它是以 m=0 为偶对称的。为偶对称的。由偏差的定义可知:由偏差的定义可知:1|01()()()Nmnr mx n x nmN()r m()()biaEE1|01|1|00()()()1()()()11 ()()()()()|()()()NmnNmNmnnbia r mE r mr mEx n x nmr mNE x n x nmr mr mr mNNNmmr mr mr mNN 2.间接法间接法可以看出:可以看出:对于一个固定的延迟对于一个固定的延迟|m|,当,当 N时,时,因此,因此 是对是对 r(m
9、)的渐进无偏估计;的渐进无偏估计;对于一个固定的对于一个固定的 N,只有当,只有当|m|N时,时,的均值才接近于的均值才接近于真值真值r(m),即当,即当|m|越接近于越接近于 N 时,估计的偏差越大;时,估计的偏差越大;的均值是真值的均值是真值 r(m)和一三角窗函数和一三角窗函数的乘积,的乘积,w(m)的长度是的长度是 2N-1。该窗函数对。该窗函数对r(m)加权,致使加权,致使 产产生了偏差。生了偏差。|()()mbia r mr mN ()0bia r m()r m()r m()r m|,01()0,NmmNNw mmN()r m|()()NmE r mr mN2.间接法间接法三角窗三
10、角窗w(m):当我们对一个信号做自然截短时,就不可避免地对该数据施加当我们对一个信号做自然截短时,就不可避免地对该数据施加了一个矩形窗,由此矩形窗就产生了加在自相关函数上的三角了一个矩形窗,由此矩形窗就产生了加在自相关函数上的三角窗,该三角窗影响自相关函数的估计质量。窗,该三角窗影响自相关函数的估计质量。2.间接法间接法由方差的定义可知:由方差的定义可知:当当 N 时,时,又因为,又因为 ,所以,对固,所以,对固定的延迟定的延迟|m|,是是 r(m)的渐进一致估计。的渐进一致估计。22212(1)var()()()()()11()()()NmiNmr mEr mE r mE r mE r mm
11、irir im r imNN 2var()EEvar()0r mlim()0Nbia r m()r m222mse()()var EEEbia2.间接法间接法计算计算 时,如果时,如果 N 和和 m 都比较大,则需要的乘法次数很多。都比较大,则需要的乘法次数很多。可以利用可以利用FFT实现对实现对 的快速计算。的快速计算。上式也可以写为:上式也可以写为:求求 的离散时间傅里叶变换,得:的离散时间傅里叶变换,得:1|01()()()Nmnr mx n x nmN()r m()r m101()()()NNNnr mxn xnmN()r m111(1)(1)0110(1)1()()()1()()NN
12、Nj mj mNNmNmNnNNj mNNnmNr m exn xnm eNxnxnm eN 2.间接法间接法把把 xN(n)补补 N 个零,得个零,得 x2N(n),即:,即:记记 x2N(n)的傅里叶变换为的傅里叶变换为 X2N(ej),则有,则有其中其中 X2N(ej)为有限长信号为有限长信号 x2N(n)的能量谱,除以的能量谱,除以 N 以后即为以后即为功率谱。这说明自相关函数的估计值功率谱。这说明自相关函数的估计值 和和 x2N(n)的功率谱是的功率谱是一对傅里叶变换。一对傅里叶变换。2(),0,1,1()0,1,21NNxnnNxnnNN1211()22(1)0(1)2121220
13、0221()()()1()()1()NNNj mj njm nNNmNnmNNNj nj lNNnljNr m exn exnm eNxn exl eNXeN()r m2.间接法间接法利用利用FFT计算自相关函数的步骤:计算自相关函数的步骤:对对 xN(n)补补 N 个零,得个零,得 x2N(n),对,对 x2N(n)做做 DFT 得得 X2N(k),k=0,1,2N-1;求求 X2N(k)的幅平方,然后除以的幅平方,然后除以N,得,得 ;对对 做逆变换,得做逆变换,得 。将将 中中 的部分向右平移的部分向右平移 2N点后形成的序列即点后形成的序列即为为 。221()NXkN221()NXkN
14、0()r m()r m(1)0Nm0()r m3.直接法和间接法的关系直接法和间接法的关系直接法:直接法:间接法:间接法:其中自相关函数其中自相关函数 与与 x2N(n)的功率谱是一对傅里叶变换:的功率谱是一对傅里叶变换:因此有因此有21()()PERNPkXkN2221()()NPERNPkXkN()()NNXkDFT xk22()()NNXkDFT xk()(),1Mjj mBTmMPer m eMN122(1)1()()Nj mNmNr m eXkN()r m221()|()()NNBTMNBTPERPkPkPk1)1(2)()(NNmnjNBTemrkP令令M=N-13.直接法和间接法
15、的关系直接法和间接法的关系由此可知,直接法可以看作是间接法的一个特例,即当间接法由此可知,直接法可以看作是间接法的一个特例,即当间接法中使用的自相关函数的最大延迟中使用的自相关函数的最大延迟 M=N-1 时,二者是相等的。前时,二者是相等的。前面已经指出:面已经指出:这就意味着,当这就意味着,当 M 较大,特别是接近于较大,特别是接近于 N-1 时,时,对对 r(m)的的估计偏差变大,此时估计出的功率谱的质量也必然下降。因此,估计偏差变大,此时估计出的功率谱的质量也必然下降。因此,在使用间接法时,都是取在使用间接法时,都是取MN-1,很显然此时,很显然此时令令MN-1,这意味着对最大长度为,这
16、意味着对最大长度为 2N-1 的自相关函数的自相关函数 的的截短,亦即施加了一个窗函数,记之为截短,亦即施加了一个窗函数,记之为 v(m),得:,得:|()()mbia r mr mN()r m()()BTPERPkPk()r m()()()Mrmr m v m3.直接法和间接法的关系直接法和间接法的关系 的均值等于真实的自相关函数的均值等于真实的自相关函数 r(m)乘以三角窗乘以三角窗 w(m),这,这是第一次加窗。该三角窗是由数据截短而产生的,其宽度为是第一次加窗。该三角窗是由数据截短而产生的,其宽度为 2N-1。v(m)是对自相关函数是对自相关函数 r(m)的第二次加窗,宽度为的第二次加
17、窗,宽度为 2M-1,MN-1。因为。因为 v(m)的宽度远小于的宽度远小于 w(m),所以,所以 v(m)的频谱的频谱V(ej)主瓣的宽度远大于主瓣的宽度远大于 w(m)的频谱的频谱 W(ej)主瓣的宽度。这时,主瓣的宽度。这时,对对 r(m)施加施加 v(m)的作用等效于在频域做的作用等效于在频域做 和和V(ej)的卷积,的卷积,这样就起到了对周期图平滑的作用。这样就起到了对周期图平滑的作用。直接法和间接法往往结合起来使用,步骤如下:直接法和间接法往往结合起来使用,步骤如下:对对 xN(n)补补N 个零,求个零,求 ;做做 的傅里叶逆变换得的傅里叶逆变换得 ,这时,这时|m|M=N-1;对
18、对 加窗函数加窗函数 v(m),这时,这时|m|MN-1,得,得 ;求求 的傅里叶变换即可得到的傅里叶变换即可得到 。()r m()jPERPe2()NPERPk2()NPERPk()r m()r m()Mrm()Mrm()jBTPe3.直接法和间接法的关系直接法和间接法的关系直接法和间接法的关系可归纳如下:直接法和间接法的关系可归纳如下:x(n)截短截短乘矩形窗乘矩形窗d0(n)求线性求线性相关函数相关函数DFTIDFTxN(n)221()NXkN()r m2()NPERPk1MN()()r m v m()BTPk2()()()NBTPERPkPkV k3.直接法和间接法的关系直接法和间接法
19、的关系一个实际的例子(一个实际的例子(fs=250Hz):):直接法直接法间接法间接法M=1004.直接法和间接法估计的质量直接法和间接法估计的质量当当 M=N-1 时,直接法和间接法估计出的结果时,直接法和间接法估计出的结果 和和 是是相同的,因此可以把这两个估计方法的质量一起讨论。为了书相同的,因此可以把这两个估计方法的质量一起讨论。为了书写的方便,下面把写的方便,下面把 和和 都简写为都简写为 。由偏差的定义可知:由偏差的定义可知:其中其中()jBTPe()jPERPe()jBTPe()jPERPe()jP e()()()jjjbia P eE P eP e|()()NmE r mr m
20、N1(1)1(1)1(1)()()()|()()()()()1()()2Njj mmNNj mmNNj mjjmNjjE P eEr m eNmr m eNr m w m eP eW eP eW ed 4.直接法和间接法估计的质量直接法和间接法估计的质量三角窗三角窗w(n)是由矩形窗是由矩形窗 d0(n)做自相关得到的。记做自相关得到的。记 W(ej)和和D0(ej)分别是分别是 w(n)和和 d0(n)的傅里叶变换,则有:的傅里叶变换,则有:因此有:因此有:当当 时,时,d0(n)趋于无限宽,趋于无限宽,W(ej)和和 D0(ej)都趋于都趋于 函数,函数,这时这时因此,对于固定的数据长度因
21、此,对于固定的数据长度 N,是有偏估计。当是有偏估计。当 时,时,它的期望值等于真值它的期望值等于真值P(ej),所以它又是渐进无偏的。,所以它又是渐进无偏的。222011()()sin(/2)/sin(/2)jjW eD eNNN()()()()()()jjjjjjbia P eE P eP eP eW eP e(1)/20()sin(/2)/sin(/2)jjND eeNN lim()0jNbia P e()jP eN lim()()jjNE P eP e4.直接法和间接法估计的质量直接法和间接法估计的质量由协方差的定义可知:由协方差的定义可知:当当=1=2时,依据时,依据12112212
22、1212122()()00()(00cov(),()()()()()()()()()1()()()21 ()()(2jjjjjjjjjjjjjjjjP eP eEP eE P eP eE P eE P eP eE P eE P eP eD eD edNP eD eD eN2)d()1()()()2jjjE P eP eW ed 222011()()sin(/2)/sin(/2)jjW eD eNNN4.直接法和间接法估计的质量直接法和间接法估计的质量可以得到估计的方差为:可以得到估计的方差为:当当 时,上式右边第一项趋近于零,第二项趋近于时,上式右边第一项趋近于零,第二项趋近于 。由此可知,由
23、此可知,是真实功率谱是真实功率谱 P(ej)的渐进无偏估计,但却不的渐进无偏估计,但却不是一致估计。不管是一致估计。不管 N 如何选取,估计值的方差总大于等于估计如何选取,估计值的方差总大于等于估计值均值的平方。值均值的平方。我们知道,我们知道,是是 r(m)的一致估计,但把的一致估计,但把 做傅里叶变换得做傅里叶变换得到的功率谱到的功率谱 却不是却不是 P(ej)的一致估计,所以功率谱的估计的一致估计,所以功率谱的估计要比相关函数的估计复杂的多。要比相关函数的估计复杂的多。22001var()()()()()2jjjP eP eDDdE P eNN 2()jP e()jP e()r m()r
24、 m()jP elim()()jjNE P eP e如果选择一个好的窗函数,使其频谱在主瓣以外的部分基本为如果选择一个好的窗函数,使其频谱在主瓣以外的部分基本为零(这样的窗函数是不存在的),如左下图所示,其中零(这样的窗函数是不存在的),如左下图所示,其中B1是主是主瓣的宽度。如果限定瓣的宽度。如果限定B1/2(-B1/2),则有,则有 D0(-)D0(+)=0,如右下图所示。如右下图所示。很显然,此时估计的方差为:很显然,此时估计的方差为:4.直接法和间接法估计的质量直接法和间接法估计的质量222001var()()()()()=()2jjjjP eP eDDdE P eE P eND0()
25、0 0 12B12BD0(+)D0(-)12B12B12B如果限定如果限定 1 和和 2 在在 0 B1,如,如下图所示。下图所示。此时估计的协方差为:此时估计的协方差为:可见可见 在在1和和 2上是不相关的,因此上是不相关的,因此 呈现较大的起伏。呈现较大的起伏。4.直接法和间接法估计的质量直接法和间接法估计的质量 0 0 12B12BD0(1-)D0(-2)12B12BD0(2+)D0(1-)12B12B1-2 2 1 1220102201021cov(),()()()()21 ()()()02jjjjP eP eP eDDdNP eDDdN()jP e()jP e4.直接法和间接法估计的
26、质量直接法和间接法估计的质量回忆平稳随机信号回忆平稳随机信号 X(n)的自相关函数和功率谱的定义:的自相关函数和功率谱的定义:通常求不出集总平均意义下的自相关函数和功率谱,因而假定通常求不出集总平均意义下的自相关函数和功率谱,因而假定X(n)是各态遍历的,取其一个样本是各态遍历的,取其一个样本 x(n),于是有,于是有:尽管自相关函数可以用时间平均来代替集总平均,但功率谱必尽管自相关函数可以用时间平均来代替集总平均,但功率谱必须保留集总平均。这是因为对随机过程须保留集总平均。这是因为对随机过程 X(n)的每一次实现的每一次实现 x(n),其傅里叶变换仍然是一个随机过程,在每一个频率处,它都是其
27、傅里叶变换仍然是一个随机过程,在每一个频率处,它都是一个随机变量,因此集总平均是必须的。这也说明,一个随机变量,因此集总平均是必须的。这也说明,P(ej)并不并不具备各态遍历性。具备各态遍历性。()()()r mE X n X nm()()jj mmP er m e1()lim()()21NNnNr mx n x nmN21()lim()21Njj nNnNP eEx n eN4.直接法和间接法估计的质量直接法和间接法估计的质量对功率谱的估计:对功率谱的估计:既无求均值运算,也无求极限运算,它只能看作是对真实功率既无求均值运算,也无求极限运算,它只能看作是对真实功率谱谱 P(ej)做均值运算时
28、的一个样本。缺少了统计平均,当然也做均值运算时的一个样本。缺少了统计平均,当然也就产生了较大的方差。就产生了较大的方差。根据根据也可以看出来,由单个样本也可以看出来,由单个样本 x(n)估计出的功率谱不会收敛到真估计出的功率谱不会收敛到真实功率谱。实功率谱。21()()jjNP eXeN2lim var()()0 jjNP eE P e4.直接法和间接法估计的质量直接法和间接法估计的质量 的频率分辨率是指的频率分辨率是指 保持真实功率谱保持真实功率谱 P(ej)中两个靠的中两个靠的很近的谱峰仍被分辨出来的能力。决定很近的谱峰仍被分辨出来的能力。决定 分辨率的主要因素分辨率的主要因素是所使用的数
29、据的长度,也即数据窗是所使用的数据的长度,也即数据窗 d0(n)的宽度。若数据的长的宽度。若数据的长度为度为 tp,采样频率为,采样频率为 fs,采样后的点数为,采样后的点数为 N,即,即 tp=N/fs,那么,那么,估计谱估计谱 的分辨率正比于的分辨率正比于 fs/N 或或 2/N。长度为。长度为 N 的各种窗的各种窗函数,其主瓣宽度为函数,其主瓣宽度为 2k/N,所以的分辨率也正比于,所以的分辨率也正比于 2k/N。若若P(ej)中有两个相距为中有两个相距为 BW 的谱峰,为了要区分他们,则要求的谱峰,为了要区分他们,则要求这样,数据的长度这样,数据的长度 N 应该满足:应该满足:为了保证
30、为了保证 的分辨率,希望的分辨率,希望 N 要大。但增大要大。但增大 N 时,又使时,又使起伏加剧。起伏加剧。()jP e()jP e()jP e()jP e2/kNBW2/NkBW()jP e()jP e4.直接法和间接法估计的质量直接法和间接法估计的质量当当 MN-1 时,直接法和间接法估计出的结果时,直接法和间接法估计出的结果 和和 不不相等,相等,是对是对 的平滑。此时:的平滑。此时:p 间接法也是一种有偏估计。由于在间接法也是一种有偏估计。由于在 上施加了一个较短的上施加了一个较短的窗口,使得间接法估计的偏差大于直接法,但方差小于直接法。窗口,使得间接法估计的偏差大于直接法,但方差小
31、于直接法。p 谱的平滑(方差的减小)是以频率分辨率为代价的。谱的平滑(方差的减小)是以频率分辨率为代价的。由于由于V(ej)的主瓣比的主瓣比W(ej)的主瓣宽,因而致使其分辨率下降。的主瓣宽,因而致使其分辨率下降。()jBTPe()jPERPe()jBTPe()jPERPe()()()|()()()()()()()()Mjj mBTmMMj mmMMj mmMjjjPer m v m eNmr m v m eNr m w m v m eP eW eV e()r m()jBTPe4.直接法和间接法估计的质量直接法和间接法估计的质量p对对 ,在,在 0 B1 时,时,和和 是不相关的,这时主瓣的宽
32、度是不相关的,这时主瓣的宽度B1=4/N。对,。对,,也可认为在上述频率范围内,当也可认为在上述频率范围内,当|2-1|B1时,时,和和 是不相关的。不过此时的是不相关的。不过此时的 B1=4 /M,因为,因为MN-1,所以所以 B1 增大。因此临近频率上的估计值变得较为相关。从这一增大。因此临近频率上的估计值变得较为相关。从这一角度也可以解释角度也可以解释 对对 平滑的原因。平滑的原因。()jBTPe()jPERPe1()jPERPe2()jPERPe1()jBTPe2()jBTPe()jBTPe()jPERPe5.直接法估计的改进直接法估计的改进直接法的估计结果直接法的估计结果 性能不好。
33、当数据长度性能不好。当数据长度N太大时,谱太大时,谱线起伏加剧,线起伏加剧,N 太小时,谱的分辨率又不好。因此需要改进。太小时,谱的分辨率又不好。因此需要改进。这里所说的改进,主要是改进其方差特性。这里所说的改进,主要是改进其方差特性。p 间接法是对直接法的一种改进,又称为周期图的平滑。间接法是对直接法的一种改进,又称为周期图的平滑。p 对直接法进行改进的另外一种方法是所谓的平均法,其指导对直接法进行改进的另外一种方法是所谓的平均法,其指导思想是把长度为思想是把长度为 N 的数据的数据 xN(n)分成分成 L 段,分别求每一段的功段,分别求每一段的功率谱,然后加以平均,以达到减小方差的目的。几
34、种主要的改率谱,然后加以平均,以达到减小方差的目的。几种主要的改进方法:进方法:()jPERPeBartlett 法法Welch 法法Nuttall 法法5.直接法估计的改进直接法估计的改进由概率论可知,对具有相同均值由概率论可知,对具有相同均值 和方差和方差 2 的独立随机变量的独立随机变量X1,X2,XL,新随机变量新随机变量 X=(X1+X2+XL)/L 的均值也是的均值也是,但方差为但方差为2/L。此即。此即 Bartlett 法的指导思想。法的指导思想。将数据分成将数据分成 L 段,每段的长度都是段,每段的长度都是 M,即,即 N=LM,第,第 i 段数据段数据为:为:然后分别计算每
35、一段的功率谱:然后分别计算每一段的功率谱:求平均,得到平均周期图:求平均,得到平均周期图:1()(1)(),1iNNxnxniM d niL 2101()()Mijij nPERNnPexn eM2111011()()()LLMjijij nPERPERNiinPePexn eLLM 12111()()()1()()()()LjijijPERPERPERijjjjE PeE PeE PeLP eD eP eW eL 的均值为:的均值为:式中式中 D1(ej)是矩形窗是矩形窗d1(n)的频谱,的频谱,W1(ej)是是 d1(n)做自相关所做自相关所得到的三角窗得到的三角窗 w1(n)的频谱,的频
36、谱,w1(n)的长度是的长度是 2M-1。因为。因为W1(ej)主瓣的宽度远大于主瓣的宽度远大于W(ej),所以平均后偏差增大,分辨率下降。,所以平均后偏差增大,分辨率下降。M 的选择主要取决于所需的分辨率。因为的选择主要取决于所需的分辨率。因为 W1(ej)的主瓣宽度是的主瓣宽度是4/M,若,若 P(ej)中有两个相距为中有两个相距为 BW 的谱峰,为了要分辨它们,的谱峰,为了要分辨它们,需要需要 4/M4/BW。如果数据长度。如果数据长度 N 已确定,根据已确定,根据所需的所需的 M,L 也就自然被确定。如果也就自然被确定。如果 N 可以变化,则应根据方可以变化,则应根据方差的要求确定差的
37、要求确定 L,然后再确定要记录的数据长度,然后再确定要记录的数据长度 N。M5.直接法估计的改进直接法估计的改进()jPERPe5.直接法估计的改进直接法估计的改进Welch 法是对法是对 Bartlett 法的改进。改进之一是在对法的改进。改进之一是在对 xN(n)分段时,分段时,允许每一段的数据有部分的交叠。例如,若每一段数据重合一允许每一段的数据有部分的交叠。例如,若每一段数据重合一半,这时的段数:半,这时的段数:式中式中 M 仍然是每段的长度。仍然是每段的长度。改进之二是,每一段的数据窗口可以不是矩形窗,例如汉宁窗改进之二是,每一段的数据窗口可以不是矩形窗,例如汉宁窗或汉明窗,记之为或
38、汉明窗,记之为 d2(n)。按按Bartlett 法求每一段的功率谱,记之为法求每一段的功率谱,记之为 ,即,即其中其中 U 是归一化因子,使用它是为了保证所得到的功率谱是渐是归一化因子,使用它是为了保证所得到的功率谱是渐进无偏估计。进无偏估计。/2/2NMLM()ijPERPe21201()()()Mijij nPERNnPexn dn eMU12201()MnUdnM5.直接法估计的改进直接法估计的改进如果如果 d2(n)是一矩形窗,则平均后的功率谱为:是一矩形窗,则平均后的功率谱为:其均值为:其均值为:记记 D2(ej)是是 d2(n)的频谱,及的频谱,及则有:则有:21211011()
39、()()()LLMjijij nPERPERNiinPePexn dn eLLMU 11()()()LjijijPERPERPERiE PeE PeE PeL 2221()()jjW eD eMU2()()()ijjjPERE PeP eW e(证明过程略)(证明过程略)5.直接法估计的改进直接法估计的改进Welch 法中各段允许交叠,因而段数增大,这样方差就可以得法中各段允许交叠,因而段数增大,这样方差就可以得到更大程度地改善。但是数据的交叠减小了每一段数据的不相到更大程度地改善。但是数据的交叠减小了每一段数据的不相关性,使方差的减小不会达到理想的程度。关性,使方差的减小不会达到理想的程度。
40、Welch 法允许分段时交叠,这样就增加了段数,当然也就增加法允许分段时交叠,这样就增加了段数,当然也就增加了做了做 FFT的次数,如果采用的数据窗非矩形窗,这又大大增加的次数,如果采用的数据窗非矩形窗,这又大大增加了做乘法的次数,因此了做乘法的次数,因此 Welch 法的计算量较大。法的计算量较大。Nuttall 等人提出的算法步骤如下:等人提出的算法步骤如下:与与Bartlett 法相同,即对法相同,即对 xN(n)自然分段(加矩形窗),且自然分段(加矩形窗),且不交叠,求得平均后的功率谱不交叠,求得平均后的功率谱 。由由 做反变换,得到做反变换,得到 对应的自相关函数对应的自相关函数 ,
41、其最大宽度是其最大宽度是 2M-1,M=N/L。()jPERPe()jPERPe()r m()jPERPe5.直接法估计的改进直接法估计的改进此步如同间接法,对此步如同间接法,对 加延迟窗加延迟窗w2(m),w2(m)的最大单边的最大单边宽度为宽度为M1,这样得到,这样得到 ,即:,即:由由 做正变换得到对做正变换得到对 x(n)功率谱的估计:功率谱的估计:很显然,该方法是把间接法和间接法结合了起来,同时也把平很显然,该方法是把间接法和间接法结合了起来,同时也把平滑和平均也结合了起来。这样就保持了平滑和平均减小方差的滑和平均也结合了起来。这样就保持了平滑和平均减小方差的优点,同时计算量也小于优
42、点,同时计算量也小于 Welch 法。法。()r m1()Mrm121()()(),Mrmr m w mmMM1()Mrm111()()Mjj mPBTMmMPerm e5.直接法的改进直接法的改进上述三种改进方法可以归纳为:上述三种改进方法可以归纳为:x(n)截短截短乘矩形窗乘矩形窗d0(n)不交叠分段不交叠分段d1(n)为矩形窗为矩形窗求平均求平均功率谱功率谱交叠交叠分段分段对每一段加权,对每一段加权,d2(n)可以不是矩形窗可以不是矩形窗平均平均做逆做逆变换变换加权加权w2(n)正正变换变换BartlettWelchNuttall5.直接法的改进直接法的改进仍然考虑前面的例子:仍然考虑前
43、面的例子:Bartlett,M=100Welch,M=100Welch,M=250M1=506.经典功率谱估计方法总结经典功率谱估计方法总结不论是直接法、间接法还是直接法的改进方法,都可以利用不论是直接法、间接法还是直接法的改进方法,都可以利用FFT 实现快速计算,且物理概念明确,仍是目前较常用的谱估实现快速计算,且物理概念明确,仍是目前较常用的谱估计方法。计方法。谱的分辨率较低,它正比于谱的分辨率较低,它正比于2/N,N 是所使用数据的长度。是所使用数据的长度。由于窗函数的影响不可避免,使得功率谱由于窗函数的影响不可避免,使得功率谱P(ej)在窗口主瓣内在窗口主瓣内的功率向边瓣部分泄露,降低
44、了分辨率。较大的边瓣有可能掩的功率向边瓣部分泄露,降低了分辨率。较大的边瓣有可能掩盖盖P(ej)中较弱的成分,或是产生假的峰值。中较弱的成分,或是产生假的峰值。方差性能不好,不是方差性能不好,不是P(ej)的一致估计,且的一致估计,且 N 增大时谱线起伏增大时谱线起伏加剧。加剧。周期图的平滑和平均与窗函数的使用密不可分。平滑和平均周期图的平滑和平均与窗函数的使用密不可分。平滑和平均主要是用来改善周期图的方差性能,但往往又减小了分辨率和主要是用来改善周期图的方差性能,但往往又减小了分辨率和增大了偏差。增大了偏差。7.短时傅里叶变换短时傅里叶变换对于平稳信号,前述的经典功率谱估计都是建立在传统的傅
45、里对于平稳信号,前述的经典功率谱估计都是建立在传统的傅里叶变换的基础之上。其实,傅里叶变换在信号的分析中自身就叶变换的基础之上。其实,傅里叶变换在信号的分析中自身就存在着不足。存在着不足。回忆福利叶变换的表达式:回忆福利叶变换的表达式:显然,如果我们想知道在某一个特定时间所对应的频率是多少,显然,如果我们想知道在某一个特定时间所对应的频率是多少,或对某一个特定频率所对应的时间是多少,那么傅里叶变换则或对某一个特定频率所对应的时间是多少,那么傅里叶变换则无能为力。傅里叶变换的这一缺点对统计特征随时间变化的非无能为力。傅里叶变换的这一缺点对统计特征随时间变化的非平稳随机信号来说,使用起来更加困难。
46、平稳随机信号来说,使用起来更加困难。dtetxjXtj)()(dejXtxtj)(21)(7.短时傅里叶变换短时傅里叶变换设信号设信号 x(n)由三个不同频率的正弦首尾相接所组成,即由三个不同频率的正弦首尾相接所组成,即其波形和幅频特性分别为:其波形和幅频特性分别为:1121232sin(),01()sin(),1sin(),1 nnNx nnNnNnNnN123nx(n)|()|jX e07.短时傅里叶变换短时傅里叶变换短时傅里叶变换(短时傅里叶变换(short-time Fourier transform,STFT)于)于 1946 年由年由 Gabor 提出,其定义为:提出,其定义为:上
47、式中窗函数上式中窗函数 w()的作用可以从不同的角度来解释:的作用可以从不同的角度来解释:当窗函数当窗函数 w()沿着沿着 t 轴移动时,相当于不断地截取一小段轴移动时,相当于不断地截取一小段又一小段的信号,然后对每一小段的信号做傅里叶变换。又一小段的信号,然后对每一小段的信号做傅里叶变换。尽管信号尽管信号 x(t)是非平稳的,但将它分成许多小段后,我们可是非平稳的,但将它分成许多小段后,我们可以假定它的每一小段是平稳的,窗函数以假定它的每一小段是平稳的,窗函数 w()的作用是尽可能的作用是尽可能地保证所截取的每一小段都是平稳的。地保证所截取的每一小段都是平稳的。窗函数的宽度越小,时域分辨率越
48、好,局部平稳性的假设窗函数的宽度越小,时域分辨率越好,局部平稳性的假设越成立。显然,与越成立。显然,与CFT相比,相比,STFT可以得到更多的时域信息。可以得到更多的时域信息。(,)()()jxSTFT t xwt eddtetxjXtj)()(7.短时傅里叶变换短时傅里叶变换设离散信号设离散信号 x(n)长度为长度为 L,窗函数,窗函数 w(n)宽度为宽度为 M,则离散,则离散 STFT定义为:定义为:上式中上式中 N 的大小决定了窗函数沿时间轴移动的间距。的大小决定了窗函数沿时间轴移动的间距。N 越小,越小,m 的取值越多,得到的时的取值越多,得到的时-频曲线就越密集。频曲线就越密集。可以
49、从两个方面对上述公式进行理解:可以从两个方面对上述公式进行理解:p当当 m 固定不变时,固定不变时,STFT 是序列是序列 x(n)w(n-mN)的离散傅里的离散傅里叶变换,此时叶变换,此时 STFT与离散傅里叶变换具有相同的性质。与离散傅里叶变换具有相同的性质。p当当 k 固定时,固定时,STFT是时间是时间 m 的函数,此时窗口的作用相当的函数,此时窗口的作用相当于一个滤波器。于一个滤波器。2(,)()()jnkMxnSTFT m kx n w nmN e7.短时傅里叶变换短时傅里叶变换设采样间隔为设采样间隔为 Ts,窗函数的宽度为,窗函数的宽度为 M,则时间分辨率和频率,则时间分辨率和频
50、率分辨率分别定义为:分辨率分别定义为:由此可见,若由此可见,若 M 较大,则时间分辨率低,但频率分辨率高,较大,则时间分辨率低,但频率分辨率高,可以看到频率的快变化;若可以看到频率的快变化;若 M 较小,则时间分辨率高,但频较小,则时间分辨率高,但频率分辨率低,看不到频率的快变化。时间分辨率和频率分辨率分辨率低,看不到频率的快变化。时间分辨率和频率分辨率相互矛盾的。率相互矛盾的。窗函数的宽度固定时,矩形窗主瓣宽度最窄,因此频率分辨窗函数的宽度固定时,矩形窗主瓣宽度最窄,因此频率分辨率高,但旁瓣衰减大,频谱泄露现象严重;汉宁窗主瓣宽度率高,但旁瓣衰减大,频谱泄露现象严重;汉宁窗主瓣宽度较宽,因此