r_ess(k+1)=sum(y_ess)/1000; end figure(2); k=[0:99];
stem(k,r_ess,'.');
title('根据样本点估计出的前100自相关序列值'); xlabel('k');ylabel('r_ess[k]'); hold on;
realvalue=[1,zeros(1,99)]; stem(k,realvalue,'r','.');
legend('根据样本点估计出的前100自相关序列值','真实的自相关序列'); error1=r_ess-realvalue;
mean_error_b=mean(error1) var_error_b=var(error1) %%
for k=0:99 for m=0:9
for n=k+1:100
y_ess2(m+1,n)=x(n+100*m)*x(n-k+100*m); end end
r_ess2(k+1)=sum(sum(y_ess2))/1000; end
figure(3); k=0:99;
stem(k,r_ess2,'b.'); hold on;
realvalue2=[1,zeros(1,99)]; stem(k,realvalue2,'r.','.');
title('Bartlett法估计功率谱方法得出的前100个自相关序列值'); xlabel('k');ylabel('r_ess2[k]');
legend('Bartlett法估计功率谱方法得出的前100个自相关序列值','真实的自相关序列');
error2=r_ess2-realvalue2;
mean_error_c=mean(error2) var_error_c=var(error2) %%
y=zeros(1,1000);
B=[1]; A=[1,-0.9];
y=filter(B,A,x);
r_ess3=zeros(1,100); for k=0:99
for n=(k+1):1000
r_ess3(k+1)=r_ess3(k+1)+y(n)*y(n-k); end
r_ess3(k+1)=r_ess3(k+1)/1000; end
figure(4);
stem(r_ess3,'.');
title('y[n]前100个自相关序列估计值'); xlabel('k'),ylabel('r_ess3(k)'); hold on;
p=[1,zeros(1,99)]; h=filter(B,A,p); for i=1:100
h1(i)=h(101-i); end
rh=conv(h,h1); rh=rh(100:199);
realvalue3=conv(p,rh);
realvalue3=realvalue3(1:100); stem(realvalue3,'r.','.');
legend('y[n]前100个自相关序列估计值','y[n]的真实自相关序列');
2、计算机练习2:AR过程的线性建模与功率谱估计。考虑AR过程:
x(n)?a(1)x(n?1)?a(2)x(n?2)?a(3)x(n?3)?a(4)x(n?4)?b(0)v(n)
v(n)是单位方差白噪声。
(a) 取b(0)=1, a(1)=0.-7348, a(2)=-1.8820, a(3)=-0.7057, a(4)=-0.8851,产生x(n)的N=256个样点。
?x(k),并与真实的自相关序列值相比较。 (b)计算其自相关序列的估计r?x(k)的DTFT作为x(n)的功率谱估计,即: (c) 将rN?1?(ej?)?Pxk??N?1??x(k)e?jk??r1|X(ej?)|2。 N(d)利用所估计的自相关值和Yule-Walker法(自相关法),估计a(1), a(2), a(3), a(4)和b(0)的值,并讨论估计的精度。
?(e)?(e)用(d)中所估计的a(k)和b(0)来估计功率谱为:Pxj?|b(0)|21??a(k)e?jk?k?142。
(f)将(c)和(e)的两种功率谱估计与实际的功率谱进行比较,画出它们的重叠波形。 (g)重复上面的(d)~(f),只是估计AR参数分别采用如下方法:(1) 协方差法;(2) Burg方法;(3) 修正协方差法。试比较它们的功率谱估计精度。 仿真结果: (a)
图2.1x(n)的N=256个样点
(b)
图2.2自相关序列的估计值与真实的对比
图2.2中估计的自相关序列与真实的自相关序列差异较大;在k>100后,真实的自相关序列就波动得很小,而估计的自相关序列则仍有较大的波动,但总体上来言两者都随着k的增大而逐渐衰减,当k>150时,真实自相关序列逐渐趋于0,而估计出的自相关序列却仍有较大的波动,这是因为估计的点数较少,使得估计精度不够。 (c)
图2.3 估计的功率谱与真实功率谱对比
(d)Yule-Walker法(自相关法) 估计的参数值为: b(0)= 1.1537