Chapter+4中国科学院大学现代数字信号处理课程课件

第四章功率谱估计4.1引言4.2经典谱估计4.3现代谱估计中的参数建模4.4AR模型谱估计方法4.5最大熵谱估计方法4.6特征分析法谱估计功率谱定义功率谱估计的方法功率谱估计的应用信号的功率谱和其自相关函数服从一对傅里叶变换关系mjmxxjxxmrPe)()e(dPmrnjjxxxxe)e(π21)(经典谱估计方法间接方法：BT法直接方法：周期图法现代谱估计方法参数法：ARMA模型法（AR模型、MA模型、ARMA模型）非参数法：Pisarenko谐波分解法、MUSIC方法在信号处理的许多场所，要求预先知道信号的功率谱密度(或自相关函数)；常常利用功率谱估计来得到线性系统的参数估计；从宽带噪声中检测窄带信号。22)()(HPyyBT法周期图法改进的周期图法BT法是先估计自相关函数,然后进行傅里叶变换得到功率谱。有偏自相关函数估计的误差相对较小，是一种渐近一致估计：1||0*)()(1)(ˆmNnxxmnxnxNmr-BTˆˆ(e)()ejjωmxxmPrm周期图法的定义如下:212j-j011ˆ(e)()eNnjxxnPxnXeNN观测数据x(n)FFT取模的平方1/NPxx(ej)^1.周期图与BT法的等价关系2101111**()00001ˆ(e)()e11()e()e()()eNjjnxxnNNNNjkjnjknknknPxnNxkxnxkxnNN令m=k-n,即k=m+n，则1*01j(1)1(1)BTˆ(e)eeˆ(e1()()ˆ())NjmxxmNmnxNNjmmxNjPPxnxmnNrm利用有偏自相关函数的BT法和周期图法是等价的。2.周期图法谱估计质量分析1)1j-j(1)1-j(1)-jˆ()||(ˆ[(e)]ee()()e)NmxxmNNmmNmxxBxmxxxErmNmrEPwmrmmN其它0||||)(NmNmNmwB式中上式在频域表示为：j-1ˆ[(e)](e)ed2πjjxxBxxEPWPΩ式中)]([)e(mrFTPxxjxx2)2/sin()2/sin(1)]([)e(NNmwFTWBjB周期图的统计平均值等于它的真值卷积三角谱窗函数，因此周期图是有偏估计，但当N→∞时，wB(m)→1，三角谱窗函数趋近于δ函数，周期图的统计平均值趋于它的真值，因此周期图属于渐近无偏估计。2)用这种方法估计的功率谱在σ2x附近起伏很大，故周期图是非一致估计。白噪声的周期图01235432100123543210Pxx(ej)Pxx(ej)54321001230123543210Pxx(ej)Pxx(ej)//(a)(b)//(c)(d)N＝32N＝64N＝256N＝128Bartlett平均周期图法窗口处理法平均周期图Welch法（修正的周期图求平均法）主要思想：对序列x(n)进行L次独立观测或将其分成L段，计算每组观测数据的周期图，再将L个周期图加和后求平均。LiixxILP1j)(1)e(ˆ假设随机信号x(n)的观测数据区间为：0≤n≤N-1，共进行了L次独立观测，得到L组记录数据，每一组记录数据用xi(n),i=1,2,3,…,L表示；或对长为N的数据x(n)分成L段，每段有M个数据，N=LM，第i段数据表示为xi(n)=x(n+iM-M)。第i组的周期图用下式表示：210)(1)(MnnjiienxMI估计方法：将得到的L个周期图进行平均，作为信号x(n)的功率谱估计，公式如下：LiixxILP1j)(1)e(ˆ估计效果分析：平均周期图的估计方差是周期图的方差的1/L，L越大方差越小，功率谱越平滑；相应的，M越小，偏移越大，分辨率越低；估计的均方误差也减少；以分辨率的降低换取了估计方差的减少，估计量的方差和分辨率是一对矛盾。平均周期图法Pxx(ej)/032101230321012303210123Pxx(ej)/Pxx(ej)/(a)(b)(c)N＝256，L＝2N＝256，L＝4N＝256，L＝8主要思想：用一适当的功率谱窗函数W(ejω)与周期图进行卷积，来达到使周期图平滑的目的。ˆ(e)()(e)jjxxNPIW()0jWe其中，π()1ˆ(e)()(e)()(e)d2πjjjxxNNPIWIW式中mxxNNmNmrIj-1)1(e)(ˆ)(de)e(π21)(jjππnWnw-(M-1)≤n≤M-1估计方法：那么mxxMMmjxxemwmrPj1)1()()(ˆ)e(ˆmxxMMmxxmwmrEPEj-1)1(je)()](ˆ[)]e(ˆ[又1j-j(1)ˆ[()]()()()eMmxxxxBmMEPermwmwm偏移分析：估计效果分析：||ˆ[()]()()()xxxxxxBNmErmrmrmwmN可得周期图的窗函数法仍然是有偏估计,其偏移和wB(m)、w(m)两个窗函数有关。如果w(m)窗的宽度比较窄，M比N小得多，这样|m|N，则wB(m)~1,d)e()e(π21e)()()]e(ˆ[)-j(jππj-1)1(jWPmwmrPExxmxxMMmxx由于w(m)比wB(m)窄，W(ejw)的主瓣比WB(ejw)宽，故可以利用窗函数法进一步平滑周期图，减小估计方差；但相应的会增加偏移，降低频率分辨率。主要思想：对Bartlett法进行修正，使之更适合FFT计算。选择适当的窗函数w(n)，并在周期图计算前直接加进去；在分段时，可使各段之间有重迭，这样将会使方差减小。估计方法：首先把数据长度为N的信号x(n)分成L段，每一段数据长度为M，N=LM；然后把窗函数w(n)加到每一个数据段上，求出每一段的周期图，形成修正的周期图；再对每一个修正的周期图进行平均。第i段的修正周期图为210j-e)()(1)(MnωniinwnxUIi=1,2,3,…，L式中102)(1MnnwMU同样,将每一段的修正的周期图之间近似看成互不相关,最后功率谱估计为LiixxILP1j)(1)e(ˆ对上式求统计平均,得到d)e()e(π21)]e(ˆ[)-j(jπjWPPExxxx式中210jje)(1)e(MnnwMUW估计效果分析：估计是渐近无偏的；这种方法对窗函数没有限制，无论什么样的窗函数均可使谱估计非负；分段时，为了减少因分段数增加给分辨率带来的影响，可使各段间有重叠，例如重叠50%。存在问题：经典谱估计中，采用观测到的N个样本值估计功率谱，认为在此观察到的N个数据以外的x(n)=0，这与实际情况是不符合的，从而造成分辨率的降低。解决方法：根据已观察到的数据，选择一个正确的模型，认为x(n)是白噪声通过此模型产生的，那么就不必认为N个以外的数据全为零了。方法步骤：（1）选择合适的信号模型；（2）根据x(n)有限的观测数据，或者它的有限个自相关函数估计值，估计模型参数;（3）基于此模型，计算输出功率谱。2j2j|)e(|)e(HPwxx对于具有尖峰的谱，应该选用具有极点的模型，如AR和ARMA模型；对于具有平坦的谱峰和深谷的信号，可以选用MA模型；既有极点也有零点的谱应选用ARMA模型，相对地说,ARMA模型适用范围较宽。在选择模型合适的基础上，应尽量减少模型的参数。对AR(2)信号的模型选择1050－5－101050－5－101050－5－101050－5－101520250.0－0.1－0.2－0.3－0.4－0.50.10.20.30.40.50.0－0.1－0.2－0.3－0.4－0.50.10.20.30.40.50.0－0.1－0.2－0.3－0.4－0.50.10.20.30.40.50.0－0.1－0.2－0.3－0.4－0.50.10.20.30.40.5////(a)(b)(c)(d)dBdBdBdB真实AR(2)的PSDMA(2)模型真实AR(2)的PSDMA(5)模型真实AR(2)的PSDMA(2)模型真实AR(2)的PSDMA(10)模型1050－5－101520251050－5－101520250.0－0.1－0.2－0.3－0.4－0.50.10.20.30.40.50.0－0.1－0.2－0.3－0.4－0.50.10.20.30.40.5//(e)(f)dBdB真实AR(2)的PSDMA(10)模型真实AR(2)的PSDMA(5)模型对AR(2)信号的模型选择在计算得到时间序列数据的样本自相关系数之后,可以通过求解Yule—Walker方程得到自回归模型的系数估计量.假设模型的差分方程和系统函数分别用下式表示：qiiPiiinwbinxanx01)()()()()(1)(10zAzBzazbzHPiiiqiiiAR模型的系数和信号自相关函数之间的关系21(0)(1)()1(1)(0)(1)00()(1)(0)xxxxxxwxxxxxxpxxxxxxrrrparrrparprpr通过求解Yule-Walker方程求模型参数：PlwxxAPlxxAxxlrlhmrlhmr121)()()1()()(m≥1m=0估计功率谱的方法首先根据信号观测数据估计信号自相关函数；再按照所选择信号模型，解上面相应的方程，求出模型参数；最后按照下式求出信号的功率谱:21022jω2jωe1e|)e(|)e(PiiiqiiiwwxxabHP4.3.3AR模型谱估计的性质1、AR模型的线性预测piipiezazH11)(AR模型的系统函数为PiiizazH111)(线性一步预测误差滤波器的系统函数为当api=ai(i=1,2,3,…,p)时，He(z)和H(z)互为逆滤波器，He(z)=1/H(z)，因此He(z)也称为白化滤波器。He(z)x(n)e(n)H(z)w(n)()()()()()()()()()()eexnwnhnenxnhnwnhnhnwn利用上述AR模型与线性预测之间的关系，可以实现预测解卷积。最小相位系统(minimum-phasesystem)在一定的幅频特性情况下，其相移为最小的系统，也称最小相移系统。这种系统的系统函数（亦称网络函数或传递函数）与非最小相位系统相比，二者的幅频响应特性是相同的，但前者的相位绝对值则较后者为小。在保持系统函数的幅频响应特性不变的情况下，使其相位最小的充分必要条件是：对于模拟信号系统，要求其零点(即使系统函数为零的复频率值)仅位于S平面（即复频域平面）的左半平面或虚轴上；对于离散信号系统，则要求其零点仅位于Z平面（即离散信号复频域平面)的单位圆内或单位圆上。常可用于进行相位校正。2、预测误差滤波器的最小相位特性AR模型H(z)必须因果稳定，即极点均在单位圆内，才能保证信号x(n)是平稳随机信号，于是He(z)应为最小相位系统。当最佳P阶线性预测系数与AR模型参数相同时，由此得到的极点保证在单位圆内，AR滤波器稳定，预测误差滤波器He(z)或者A(z)是最小相位系统。3、AR模型