资源描述:
《参数法功率谱估计.docx》由会员上传分享,免费在线阅读,更多相关内容在教育资源-天天文库。
1、参数法功率谱估计一、信号的产生(一)信号组成在本实验中,需要事先产生待估计的信号,为了使实验结果较为明显,我产生了由两个不同频率的正弦信号(频率差相对较大)和加性高斯白噪声组成的信号。(二)程序N=1024;n=0:N-1;xn=2*cos(2*pi*0.2*n)+cos(2*pi*0.213*n)+randn(1,1024);这样就产生了加有白噪声的两个正弦信号其波形如下二、参数模型法功率谱估计(一)算法原理简介1.参数模型法是现代谱估计的主要内容,思路如下:①假定所研究的过程是由一个白噪声序列激
2、励一个因果稳定的可逆线性系统的输出;②由已知的,或其自相关函数估计的参数;③由的参数来估计的功率谱。2.自回归模型,简称AR模型,它是一个全极点的模型。“自回归”的含义是:该模型现在的输出是现在的输入和过去p个输出的加权和。此模型可以表现为以下三式:①;②;③。3.AR模型的正则方程建立了参数和的自相关函数的关系,公式如下:时,时。(二)两种AR模型阶次的算法1.Yule-Walker算法(自相关法)(1)算法主要思想Yule-Walker算法通过解Yule-Walker方程获得AR模型参数。从低阶
3、开始递推,直到阶次p,给出了在每一个阶次时的所有参数。公式如下:①;②;③。(2)运算简要框图估计功率谱密度(2M-1)点FFTX(n)求自相关估计参数输出Yule-Walker法谱估计运算简要框图(3)程序示例clearall;closeall;N=512;n=0:N-1;xn=2*cos(2*pi*0.2*n)+cos(2*pi*0.213*n)+2*randn(1,512);Rx=zeros(1,N+1);%从课本上的公式来看,Rx(m)中的m属于(0,m),即共有m+1个,故在这里设Rx是一
4、个一行,N+1列的向量figure(1)plot(n,xn);title('(a)两正弦信号加白噪声波形')%下面用书中所讲自相关函数估计中的渐进无偏估计来估计自相关函数form=1:N+1;%由于在matlab中,下角标不能是0,m属于(0,m),在此只能从1到N+1sum=0;forn=1:(N+1-m);%同样道理,把书中公式里m换成m-1,N换成N+1,求和下限变为1sum=sum+xn(n).*xn(n+m-1);endRx(m)=sum/N;%切记,这里的Rx(1)才是自相关函数在0点的
5、取值。Rx(m)只是一个存储数据的代号,为了跟书中公式一致,才叫Rxend%下面估计各参数P=50;a=zeros(P,P);%a中有两个变量m,i,所以设a是P行P列的向量km=zeros(1,P);%因为阶次是P,故反射系数有P个p=zeros(1,P+1);%由于matlab中没有ρ,故用p来代替表示,ρ的范围是(0,P)共有P+1个%下面计算AR模型参数的各个初始化值p(1)=Rx(1);a(1,1)=-Rx(2)/Rx(1);km(1)=a(1,1);p(2)=Rx(1).*(1-abs(
6、a(1,1).^2));form=2:P%由于m=1时的各个值在上面已经给出,故从m=2开始求sum1=0;fori=1:m-1sum1=sum1+a(m-1,i).*Rx(m-i+1);enda(m,m)=-(Rx(m+1)+sum1)/p(m);%km(m)=a(m,m);求出kmfori=1:m-1a(m,i)=a(m-1,i)+a(m,m)*a(m-1,m-i);endp(m+1)=p(m)*(1-abs(a(m,m)).^2);endz=[1,a(P,:)];G=sqrt(p(P));[H
7、w]=freqz(G,z,512);%调用计算数字滤波器频响的函数figure(2)plot(w/(2*pi),10*log10(abs(H).^2));title('自相关法')ylabel('10log(PSD)')title('(b)yule-walker法估计功率谱密度')(4)结果分析从波形图中可以十分清晰的分辨出两个不同频率的正弦波2.Burg法(1)算法主要思想Burg法不是直接估计AR模型的参数,而是先估计反射系数。使用线性预测的方法来计算不同阶数下的反射系数,其同时使用前向和后向线
8、性预测,使前向和后向预测误差的平均功率相对各阶反射系数最小,将反射系数代入Levinson-Durbin公式即可求解。(2)运算简要框图用L-D公式求解估计反射系数X(n)输出(3)程序示例clearall;closeall;N=512;n=0:N-1;xn=2*cos(2*pi*0.2*n)+cos(2*pi*0.213*n)+randn(1,512);Rx=zeros(1,N+1);%从课本上的公式来看,Rx(m)中的m属于(0,m),即共有m+1个,故在这里设R