matlab经典、现代功率谱估计

1、上机作业：】、假设一平稳随机信号为x(n)08x(n-1)+w(n),其中是均值为0,方差为1的白噪声，数据长度为1024。、产生符合要求的w(n)和x(n)；、给出信号x(n)的理想功率谱；、编写周期图谱估计函数，估计数据长度N=1024及256时信号功率谱，分析估计效果。、编写Bartlett平均周期图函数，估计当数据长度N=1024及256时，分段数L分别为2和8时信号的功率谱，分析估计效果。一、解题思路w(n)可以通过随机序列randn(1,N)来产生，x(n)可以通过对w(n)滤波产生(由递推式可得系统的传递函数)，也可以直接由递推式迭代产生。由于线性系统的输出功率谱等于输入功率谱乘

2、以传递函数模的平方，X(n)可以看做w(n)通过一线性系统的输出，H(z)=1/(1-0.8z)。所以x(n)的理想功率谱P(ejw)=o2|H(ejw)|2。w周期图方法：直接对观测数据做FFT变换，变换的结果取模的平方再除以数据长度，作为估计的功率谱。256个观测点时可以对原观测数据以4为间隔提取得到。Bartlett法：将L组独立的观测数据分别求周期图，再将L个周期图求平均作为信号的功率谱估计。L组数据可以通过对原观测数据以L为间隔提取得到。二、MATLAB实现程序及注解clc;clear;closeall;Fs=500;%采样率N=1024;%观测数据w=sqrt(1)+randn(1

3、,N);%0均值,方差为1的白噪声,长度1024x=w(1)zeros(1,N-1);%初始化x(n),长度1024,x(1)=w(1)fori=2:Nx(i)=0.8*x(i-1)+w(i);%迭代产生观测数据x(n)end%理想功率谱h,w1=freqz(x);figure,plot(w1*500/(2*pi),10*log10(abs(h).人2);gridon;title(理想功率谱);xlabel(频率);ylabel(功率db);%周期图法%1024个观测点Pxx=abs(fft(x).A2/N;%周期图公式Pxx=10*log10(Pxx(index+1);%化为dbfigure

4、;plot(k,Pxx);gridon;title(周期图1024点)；xlabel(频率);ylabel(功率db);%周期图256个观测点x1=x(1:4:N);Pxxl=abs(fft(xl,1024).人2/N;Pxx1=10*log10(Pxx1(index+1);%化为dbfigure;plot(k,Pxx1);gridon;title(周期图256点)；xlabel(频率);ylabel(功率db);%Bartlett平均周期图,N=1024%分段L=2L=2;x_21=x(1:L:N);x_22=x(2:L:N);Pxx_21=abs(fft(x_21,1024).人2/len

5、gth(x_21);Pxx_22=abs(fft(x_22,1024).人2/length(x_22);Pxx_2=(Pxx_21+Pxx_22)/L;figure;subplot(2,2,1),plot(k,10*log10(Pxx_2(index+1);gridon;title(N=1024,L=2);xlabel(频率);ylabel(功率db);%分段L=8L1=8;x3=zeros(L1,N/L1);%产生L1行,N/L1列的矩阵用以存储分组的数据fori=1:L1x3(i,:)=x(i:L1:N);%将原始数据分为8组endPxx3=zeros(L1,1024);%产生L1行,10

6、24列矩阵用以存储分组的周期图fori=1:L1Pxx3(i,:)=abs(fft(x3(i,:),1024).人2/length(x3(i,:);%分别求周期图，结果保存在Pxx3中,FFT长度为1024endfori=1:1024Pxx3_m(i)=sum(Pxx3(:,i)/L1;%求平均endsubplot(2,2,2),plot(k,10*log10(Pxx3_m(index+1);gridon;title(N=1024,L=8);xlabel(频率);ylabel(功率db);%Bartlett平均周期图,N=256，求法同上%分段L=2,分别计算周期图，再取平均x=x(1:4:N

7、);L2=2;x_31=x(1:L2:length(x);x_32=x(2:L2:length(x);Pxx_31=abs(fft(x_31,1024).A2/length(x_31);Pxx_32=abs(fft(x_32,1024).A2/length(x_32);Pxx_3=(Pxx_31+Pxx_32)/L2;subplot(2,2,3),plot(k,10*log10(Pxx_3(index+1);gridon;title(N=256,L=2);xlabel(频率);ylabel(功率db);%分段L=8L3=8;x4=zeros(L3,length(x)/L3);fori=1:L3

8、x4(i,:)=x(i:L3:length(x);%将原始数据分为8组endPxx4=zeros(L3,1024);fori=1:L3Pxx4(i,:)=abs(fft(x4(i,:),1024).人2/length(x4(i,:);%分别求周期图,FFT长度为1024endfori=1:1024Pxx4_m(i)=sum(Pxx4(:,i)/L3;%求平均endsubplot(2,2,4),plot(k,10*log10(Pxx4_m(index+1);gridon;title(N=256,L=8);xlabel(频率);ylabel(功率db);三、结果及分析理想功率谱80I117。一一丄

9、一一一-ooO54320OrHJ图1理想功率谱图2周期图1024点及256点从上图可以看出，周期图法得到的功率谱估计，谱线的起伏较大，即估计所得的均方误差较大。当N增加时，摆动的频率加快，而摆动的幅度变化不大。且N=1024时，谱的分辨率较N=256时大。N=1024,L=2N=1D24=8频率频率频率频率图3Bartlett平均周期图法由上图可以看出，采用平均处理后，谱线上下摆动的幅度减小（即均方误差有所降低），曲线的平滑性也较周期图法好。N相同时，分段越多，方差越小，曲线越平滑。这是因为,N一定时，L加大，每一段的数据量就会相应减少，因此估计方差减小，偏移加大，从而分辨率降低。N一定时，L

10、与每段数据量相互矛盾，需择中选取。综上，Bartlett法相对于周期图法来说，较好地减小了估计误差。2、假设均值为0方差为1的白噪声w(n)中混有两个正弦信号，该正弦信号的频率分别为100Hz和110Hz，信噪比分别为10dB和30dB，初始相位都为0,采样频率为1000Hz。、采用自相关法、Burg法、协方差法、修正协方差法估计功率谱，分析数据长度和模型阶次对估计结果的影响(可采用MATLAB自带的功率谱分析函数)。、调整正弦信号信噪比，分析信噪比的降低对估计效果的影响。一、解题思路信噪比为信号平均功率与噪声的平均功率之比，即SNR=10lgS/N,单位dB。假设正弦信号的幅度为A，对其在一

11、个周期内的平方进行积分再除以周期长度可得，平均功率Ps=A2/2。s=白噪声的平均功率为1,所以又SNR的计算式可得题中两个正弦信号的幅值分别为sqrt(20)和sqrt(2000)。题中，所给观测信号为白噪声混有两个正弦信号，所以功率谱中应有两个明显的尖峰，所以选择AR模型来进行估计。AR模型的估计方法有自相关法，Burg法，协方差法和修正协方差法等。本次编程采用matlab自带的功率谱分析函数。自相关，Burg法，协方差法，修正协方差法分别为pyulear()，pburg()，pcov()，pmcov()。二、MATlAB实现程序及注解clc;clear;closeall;fs=1000;

12、%采样率1000N=l;%改变数据长度p=50;%AR模型阶数nfft=512;%fft长度t=0:1/fs:N;wn=sqrt(1)+randn(1,N*fs+1);%白噪声,均值0,方差1s1=sqrt(20)*sin(2*pi*100*t);%正弦信号1,信噪比10dbs2=sqrt(2000)*sin(2*pi*110*t);%正弦信号2,信噪比30dbx=s1+s2+wn;%观测数据%figure,plot(t,x);x1=xcorr(x,biased);Pxx,f=pyulear(x1,p,nfft,fs);%Yule-Walker方程figure,plot(f,10*log10(

13、Pxx);gridon;title(自相关法)；Pxx1,f1=pcov(x,p,nfft,fs);figure,plot(f1,10*log10(Pxx1);gridon;title(协方差法)；Pxx2,f2=pmcov(x,p,nfft,fs);figure,plot(f2,10*log10(Pxx2);gridon;title(修正协方差法)；Pxx3,f3=pburg(x,p,nfft,fs);figure,plot(f3,10*log10(Pxx3);gridon;title(Burg法);三、结果及分析图150阶AR模型谱估计从上图可以看出，采用参数建模的谱估计方法得到的功率谱曲

14、线平滑（方差小），分辨率高，可以明显地观察到两个谱峰。ED2D203QQ4DD日法2D2D-3D40D300白相关法检片经法陽正协方莹法dOO1DD200300图220阶AR模型谱估计降低模型阶次后，可以发现，谱的分辨率降低（两个谱峰几乎变成一个谱峰），但是曲线平滑性变好（估计误差降低）。图350阶AR模型（4倍观测数据）谱估计与图1相比，图3几乎没有什么变化，这是因为AR模型隐含着数据的外推。在经典谱估计中，对于观测数据以外的数据均假设为0（相当于时域加窗），从而造成谱的分辨率降低。AR谱估计中，由于隐含自相关函数的外推，而使分辨率大大提高。从而观测数据的多少对分辨率影响不大，对估计误差会有影响。数据越多，估计误差越小，谱线波动越小，这一点图中可以很明显地看出。