www.pudn.com > RELS11111.rar > RELS.m

u0=idinput(15,'prbs'); 
u0=u0'; 
u=[]; 
for i=1:4 
    u=[u,u0]; 
end                        %取长度为60的M序列 
figure(1); 
subplot(211); 
stem(u),grid on  %画出M序列输入信号 
v=randn(1,60);%产生一组60个正态分布的随机噪声 
subplot(212); 
plot(v),grid on %画出随机噪声信号 
u,v  %显示输入信号和噪声信号 
%建模 
a1=-1.5;a2=0.7;b1=1;b2=0.5;d1=-1;d2=0.2; 
z=zeros(1,60); 
z(1)=v(1); 
z(2)=-a1*z(1)+b1*u(1)+v(2)+d1*v(1); 
for k=3:60 
    z(k)=-a1*z(k-1)-a2*z(k-2)+b1*u(k-1)+b2*u(k-2)+v(k)+d1*v(k-1)+d2*v(k-2); 
end 
%辨识 
theta0=[0.001 0.001 0.001 0.001 0.001 0.001]'; 
p0=10^6*eye(6,6); 
%E=0.00000000005;%相对误差E=0.00000000005 
theta=[theta0,zeros(6,14)]; 
%e=zeros(6,15); 
for k=3:60 
    h1=[-z(k-1),-z(k-2),u(k-1),u(k-2),v(k-1),v(k-2)]'; 
    x=h1'*p0*h1+1; 
    x1=inv(x); 
    k1=p0*h1*x1;%求出K值 
    m1=z(k)-h1'*theta0; 
    theta1=theta0+k1*m1;%辨识参数theta 
    v(k)=z(k)-h1'*theta1; 
   % e1=theta1-theta0; 
    %e2=e1./theta0; 
   % e(:,k)=e2; 
    theta0=theta1; 
    theta(:,k)=theta1;%把递推出的辨识参数C的列向量加入辨识参数矩阵 
    p1=p0-k1*h1'*p0;%find p(k) 
    p0=p1;%给下次使用 
   % if e2