数字信号处理_第四章

第四章功率谱估计4.1引言4.2经典谱估计4.3现代谱估计中的参数建模4.4AR模型谱估计的性质4.5AR谱估计的方法4.6最大熵谱估计与最大似然谱估计4.7特征分解法谱估计NO.24.1引言对信号和系统进行分析研究、处理有两类方法：一类是在时域进行，我们前面学习的维纳滤波、卡尔曼滤波和自适应滤波都属于这种方法；本章则是在频率域进行研究的一类方法。这两类方法都是信号处理的重要方法。对确定性信号傅里叶变换是在频率域分析研究的理论基础，但对于随机信号，其傅里叶变换不存在，因此转向研究它的功率谱。谱，就是信号的某些特征在频域随频率的分布。功率谱反映了随机信号功率的分布特性，有着很广泛的应用：在雷达信号处理中，回波信号的频谱提供了运动目标的位置、强度、速度等信息；在声纳系统中，为了寻找舰艇或潜艇也要对混有噪声的信号进行谱分析；在语音处理中，谱分析用来探测语音语调共振峰；在电子战中，还利用频谱对目标进行分类。mjmxxjxxmrPe)()e(按照Weiner-Khintchine定理，信号的功率谱和其自相关函数服从傅里叶变换关系：实际所能得到的随机信号的长度总是有限的，用有限长度的信号计算得到的功率谱只是真实功率谱的一个估计。功率谱估计分为经典谱估计和现代谱估计两大类。dPmrnjjxxxxeeπ)(21)()]()([)(*mnxnxEmrxx上式称做功率谱的定义，对于平稳随机信号，服从各态历经定理，集合平均可以用时间平均代替，由上式还可以推出功率谱的另一个定义：)()(121lim)(*mnxnxNmrNNnNxx将计算自相关函数中的集合平均用时间平均代替：mjNNnNmjxxmnxnxNPee)()(121lim)(*代入功率谱定义式，得mmnjnjNNnNmnxnxN)(*)()(121limee令l=n+m,则ljlNNnnjNxxlxenxNPeej)()(121lim)(*mmnjnjNNnNjxxmnxnxNeP)(*)()(121lim)(ee2121jnNNnNnxNe)(lim上式中x(n)是观测数据，Pxx(ejω)是随机变量，必须对Pxx(ejω)取统计平均值，得到2jω)(121lim)e(NNnnjNxxenxNEP该式被认为是功率谱的另一定义，周期图法谱估计Weiner-khintchine定理表明功率谱是无限多个自相关函数的函数，但观测数据只有有限个,只能得到有限个自相关函数。按照上面的定义公式求功率谱,也需要无限多个观测数据。因此根据有限个样本数据，计算随机信号的功率谱，是一个功率谱的估计问题现代谱估计是以信号模型为基础，下图表示x(n)的信号模型，输入白噪声w(n)的均值为0，方差为σ2w，x(n)的功率谱为：2j2j|)e(|)e(HPwxxH(z)w(n)x(n)如果由观测数据能够估计出信号模型的参数，则信号的功率谱可以按照上式计算出来，这样,估计功率谱的问题就变成了估计信号模型参数的问题。信号模型有很多种，如AR模型、MA模型等等，针对不同的情况，需要选择不同的模型。现代谱估计的质量比经典谱估计的质量有很大的提高。但遗憾的是，尚无任何理论能指导选择一个合适的模型，只能根据功率谱的一些先验知识，或者说一些重要的谱特性，来选择模型。4.2经典谱估计BT法是先估计自相关函数,然后进行傅里叶变换得到功率谱。设对随机信号x(n)，只观测到一段样本数据，n=0,1,2,…,N-1。根据这一段样本数据估计自相关函数,有两种估计方法，即有偏自相关函数估计和无偏自相关函数估计。有偏自相关函数估计的误差相对较小，这种估计是一种渐近一致估计：1||0*)()(1)(ˆmNnxxmnxnxNmr4.2.1BT法对上式进行傅里叶变换，得到BT法的功率估计:mnjxxjmrP-BTe)(ˆ)e(ˆ为了减少谱估计的方差，经常用窗函数w(m)对自相关函数进行加权，此时谱估计公式为1)1(-)()(ˆ)(ˆMMmnjxxjmwmrPeeBT式中0)()(mwmw-(M-1)≤m≤(M-1)其它为了采用FFT计算傅里叶变换，必须将求和域(-M+1,M-1)移到(0~L-1)，功率谱的计算公式为：10-)()(ˆLmωmxxmSPjjωBTee上式也被称为加权协方差谱估计。它要求加窗后的功率谱仍是非负的，这样窗函数w(m)的选择必须满足一个原则，即它的傅里叶变换必须是非负的，例如巴特利特窗就满足这一条件。k=0,1,2,…,L-1)()(ˆmSFFTkPxxBT)()(ˆ0)()(ˆ)(LmwLmrmwmrmSxxxxxx0≤m≤M-1M≤m≤L-ML-M+1≤m≤L-1按照有偏自相关函数公式估计自相关函数，已经证明这是渐近一致估计，但经过傅里叶变换后得到功率谱的估计，功率谱估计却不一定仍是渐近一致估计，可以证明它是非一致估计，是一种不好的估计方法。BT法中用有偏自相关函数进行估计时，它和用周期图法估计功率谱是等价的，因此BT法的估计质量和周期图法的估计质量是一样的。将功率谱的另一定义式重写如下:2jje)(121lim)e(NNnnNxxnxNEP如果忽略上式中求统计平均的运算，假设观测数据为:x(n)0≤n≤N-1，便得到周期图法的定义:210j-je)(1)e(ˆNnnxxnxNP4.2.2周期图法用周期图法计算功率谱框图观测数据x(n)FFT取模的平方1/NPxx(ej)^由周期图法功率谱估计公式推导它与BT法的等价关系210)(1)(ˆNnnjjxxnxNPee令m=k-n,即k=m+nmjmNnNNmxxnmxnxNPee10*1)1(j)()(1)(ˆ1.周期图与BT10*10)()(1NnnjNkkjnxkxNee1010)(*)()(1NkNnnkjnxkxNe方括号中的部分是有偏自相关函数的计算公式：mjxxNNmxxmrP)e(ˆ)e(ˆ1)1(j)()(1)(ˆ*||10nmxnxNmrmNnxx)e(ˆ)e(ˆBTjjxxPP利用有偏自相关函数的BT法和周期图法是等价的关系已知自相关函数的估计值，m=-(N-1),…-1,0,1,…,N-1，按照BT法求功率谱的统计平均值：)(ˆmrxxmxxNNmxxmrEPEj-1)1(je)(ˆ)]e(ˆ[有偏自相关函数的统计平均值在第一章中已确定，将结果代入上式,得到mxxNNmxxmrNmNPEj-1)1(j)e(||)]e(ˆ[2.周期图法谱估计质量分析1)mmxxBxxmrmwPEj-j)e()()]e(ˆ[式中其它0||||)(NmNmNmwB两序列乘积的傅里叶变换，其频域服从卷积关系)]([)e(mrFTPxxjxx)()()]e(ˆ[jxxjBjxxePeWPE2)2/sin()2/sin(1)]([)e(NNmwFTWBjBWB(ejω)称为三角窗的谱函数。上式表明，周期图的统计平均值等于它的真值与三角窗函数频谱的卷积，因此周期图是有偏估计，但当N→∞时，wB(m)→1，三角窗函数的频谱趋近于δ函数，周期图的统计平均值趋于它的真值，因此周期图属于渐近无偏估计。周期图的方差的精确表达式很繁冗，为分析简单：假设x(n)是实的零均值的正态白噪声信号，方差是σx2，即功率谱是常数σx2，其周期图用IN(ω)表示，N表示观测数据的长度。按照周期图的定义:2|)(|1)(jeXNIN推导得周期图的方差：2)周期图的方差24)sin()sin(1)](var[NNIxNnNnNkknxkxNj-jee1010)()(1当N趋于无限大时，周期图的方差并不趋于0，而趋于功率谱真值的平方，即4)](var[xNNI所以无论怎样选择N，周期图的方差总是和σ4x同一个数量级。信号功率谱的真值是σ2x，这说明周期图的方差很大。用这种方法估计的功率谱在σ2x附近起伏很大，所以周期图是非一致估计，是一种很差的功率谱估计方法。图4.2.2白噪声的周期图01235432100123543210Pxx(ej)Pxx(ej)54321001230123543210Pxx(ej)Pxx(ej)//(a)(b)//(c)(d)N＝32N＝64N＝256N＝128可以看到，随着N的增大，功率谱曲线并没有趋向真值，而是起伏变的越来越剧烈。周期图法估计功率谱不是一致估计，均方误差很大，使估计出的功率谱不可靠。其频率分辨力低是根本缺点，原因在于观测数据只有一段。由于BT法和周期图法具有等效的功率谱估计质量，因此BT法也不是一致估计，分辨率低，估计误差大。4.2.3经典谱估计方法的改进基本思想：对随机变量进行观测，得到L组独立的记录数据，用每一组数据求其周期图，然后将L个周期图加起来求平均。这样得到的周期图，其方差将是用一组数据得到的周期图的方差的1/L。210)(1)(MnnjiienxMI1.平均周期图法假设随机信号x(n)的观测区间为：0≤n≤M-1，共进行了L次独立观测，得到L组记录数据，某一组记录数据用xi(n),i=1,2,…,L表示，第i组的周期图为：将得到的L个周期图进行平均，作为信号x(n)的功率谱估计：LiixxILP1j)(1)e(ˆ为了分析偏移，对上式求统计平均:其中，窗函数的频谱为2j)2/sin()2/sin(1)]([)e(MMmwFTWBB)]([)(1)(ˆ1IEIELPELiixxje)()()]([jxxjBePeWIE表明平均周期图仍然是有偏估计，偏移和每一段的数据个数M有关；由于MN,平均周期图的偏移比周期图的偏移大，表现在三角谱窗主瓣的宽度比周期图主瓣的宽度宽。所以其分辨率更加降低，因此也可以说，偏移的大小反映分辨率的低与高。对平均周期图求方差，由于是L次独立观测，L个周期图相互独立，因此平均周期图的方差为)](var[1)]e(ˆvar[ijxxILP即平均周期图的估计方差是周期图的方差的1/L。即：是以分辨率的降低为代价换取了估计方差的减少平均周期图法Pxx(ej)/032101230321012303210123Pxx(ej)/Pxx(ej)/(a)(b)(c)N＝256，L＝2N＝256，L＝4N＝256，L＝8这种方法是用一适当的窗函数W(ejω)与周期图进行卷积，来达到使周期图平滑的目的。deπeπ)()(21)(ˆ)(jNlwxxWIP式中mxxNNmNmrIj-1)1(e)(ˆ)(是有偏自相关函数)(ˆmrxxde)e(π21)(jjππnWnw-(M-1)≤n≤M-12.mxxMMmjxxemwmrPje)()(ˆ)(ˆ1)1(周期图的窗函数法和前面提到的BT法的加权协方差谱估计是类似的。窗函数法中，周期图和窗函数的频谱卷积得到功率谱，等效于在频域对周期图进行修正，使周期图通过一个线性系统，滤除掉周期图中的快变成分，谱窗函数需具有低通特性。对上式求统计平均,得到mxxMMmxxmwmrEPEj-1)1(je)()]