递推的极大似然法辨识程序&
递推的极大似然法辨识程序
,,?例5.2系统模型如图5.5所示。试用递推的极大似然法对系统辨识的参数集
,1,210.2,z，z λv(k) 1
+ e(k) ,1,2 z，0.5zy(k)+ u(k) u ,1,21,1.2z，0.6z z(k)
图5.5 例5.2系统模型
,1v(k)随机信号，输入信号为幅值为的M序列或随机信号，要求画出程序流程图，打印出程序(程序中带有注释)和辨识中的参数、误差曲线。
解: 首先解释编程所用的部分字母:由于在MATLAB语言中无法用希腊字母描述、无
ˆP(k),Pi,(i,0,1)K(k),K,(k)法用上标及下标，故用 ‘o’和‘o1’
表
关于同志近三年现实表现材料材料类招标技术评分表图表与交易pdf视力表打印pdf用图表说话 pdf
示;令;;产生M序列时，a(i),b(i),c(i),d(i)表示四级移位寄存器的第1，2，3，4级寄存器的输出;
? 编程如下(光盘上该程序:FLch5RMLeg2.m，可在MATLAB6.I下直接运行):
编程如下:
clear %清零
a(1)=1;b(1)=0;c(1)=1;d(1)=0;u(1)=d(1);z(1)=0;z(2)=0; %初始化 for i=2:1200 %产生m序列u(i)
a(i)=xor(c(i-1),d(i-1));
b(i)=a(i-1);
c(i)=b(i-1);
d(i)=c(i-1);
u(i)=d(i);
end
u; %若取去‘;’可以在程序运行中观测到m序列
v=randn(1200,1); %产生正态分布随机数
V=0; %计算噪声方差
for i=1:1200
V=V+v(i)*v(i);
end
V1=V/1200;
for k=3:1200 %根据v和u计算z
z(k)=1.2*z(k-1)-0.6*z(k-2)+u(k-1)+0.5*u(k-2)+v(k)-v(k-1)+0.2*v(k-2);
end
o1=0.001*ones(6,1);p0=eye(6,6); %赋初值
zf(1)=0.1;zf(2)=0.1;vf(2)=0.1;vf(1)=0.1;uf(2)=0.1;uf(1)=0.1;
%迭代计算参数值和误差值
for k=3:1200
h=[-z(k-1);-z(k-2);u(k-1);u(k-2);v(k-1);v(k-2)];
hf=h;
K=p0*hf*inv(hf'*p0*hf+1);
p=[eye(6,6)-K*hf']*p0;
v(k)=z(k)-h'*o1;
o=o1+K*v(k) ;
p0=p;
o1=o;
a1(k)=o(1);
a2(k)=o(2);
b1(k)=o(3);
b2(k)=o(4);
d1(k)=o(5);
d2(k)=o(6);
e1(k)=abs(a1(k)+1.2);
e2(k)=abs(a2(k)-0.6);
e3(k)=abs(b1(k)-1.0);
e4(k)=abs(b2(k)-0.5);
e5(k)=abs(d1(k)+1.0);
e6(k)=abs(d2(k)-0.2);
zf(k)=z(k)-d1(k)*zf(k-1)-d2(k)*zf(k-2);
uf(k)=u(k)-d1(k)*uf(k-1)-d2(k)*uf(k-2);
vf(k)=v(k)-d1(k)*vf(k-1)-d2(k)*vf(k-2);
hf=[-zf(k-1);-zf(k-2);uf(k-1);uf(k-2);vf(k-1);vf(k-2)];
end
o1 %若取去‘;’可以在程序运行中观测到参数
V1
%绘图
subplot(4,1,1)
k=1:1200;
plot(k,a1,'k:',k,a2,'b',k,b1,'r',k,b2,'m:',k,d1,'g',k,d2,'k');
xlabel('k')
ylabel('parameter')
legend('a1=-1.2,','a2=0.6','b1=1.0','b2=0.5','d1=-1.0','d2=0.2'); %图标炷
title('The parameter idendification of the RML'); end
subplot(4,1,2)
k=1:1200;
plot(k,e1,'k',k,e2,'b',k,e3,'r',k,e4,'m',k,e5,'g',k,e6,'k'); xlabel('k')
ylabel('error')
%title('误差曲线')
end
subplot(4,1,3)
k=1:1200;
plot(k,u);
xlabel('k')
ylabel('input')
%title('系统输入信号')
end
subplot(4,1,4)
k=1:1200;
plot(k,v);
xlabel('k')
ylabel('random noise')
%title('系统所加的随机噪声') end
? 程序运行结果如图5.7 所示
The parameter idendification of the RML2
a1=-1.2,a2=0.60b1=1.0b2=0.5parameterd1=-1.0-2d2=0.2020040060080010001200k1.5
1error0.5
00200400600800100012001k
0.5input
00200400600800100012005k
0
random noise-5020040060080010001200k
图5.7 RML辨识参数曲线
浙公网安备 33021202002483