系统辨识作业_梯度法最小二乘法.doc

系统辨识作业 第 6 页系统辨识1、 最小二乘算法仿真对象为：，为服从N(0,1)分布的白噪声，输入信号采用4阶M序列。辨识模型为：，观测数据长度取，加权矩阵取。则matlab程序如下:% M array constructionNp=15;r=4;X1=1;X2=1;X3=1;X4=1;m_length = r*Np;a=1;for i=1:1:m_length Y4=X4;Y3=X3;Y2=X2;Y1=X1; X4=Y3;X3=Y2;X2=Y1; X1=xor(Y3,Y4); if Y4=0 M(i)=-a; else M(i)=a; endendfigure;i=1:1:m_length;plot(i,M);% 白噪声noise = zeros(1,m_length);for i=1:1:m_length temp = noise + 0.5*rands(1,m_length); noise = temp;endnoise = noise/12;%noise = temp;% parameter of systemn=2;d=1;a1=-1;a2=0.5;b1=1;b2=0.5;S_U0=0.2;S_Y0=0.2;% generate u,yu0=ones(1,m_length)*S_U0;U = M + u0 + noise;figure;i=1:1:m_length;plot(i,U);%formulationy(1)=0;y(2)=0;y(3)=0;Y(1)=S_Y0+y(1)+noise(1);Y(2)=S_Y0+y(2)+noise(2);Y(3)=S_Y0+y(3) +noise(3);for k=4:m_length y(k) = b1*U(k-1-d)+b2*U(k-2-d)-a1*y(k-1)-a2*y(k-2); Y(k)=S_Y0+y(k)+noise(k);endfigure;i=1:1:m_length;plot(i,Y);%辨识% premitive valuec=2;P = (c3)*eye(3*n+d);sita(:,3) = 0,0,0,0,0,0,0;alf = 0.995;% M%compute U0,Y0sum_U = 0;sum_Y = 0;for k=1:1:m_length sum_U = sum_U + U(k); sum_Y = sum_Y + Y(k);endt_U0 = sum_U/m_length;t_Y0 = sum_Y/m_length;for k=1:1:m_length t_u(k) = U(k) - t_U0; t_y(k) = Y(k) - t_Y0;endfigure;i=1:1:m_length;plot(i,t_u);figure;i=1:1:m_length;plot(i,t_y);v(1)=0;v(2)=0;v(3)=0;for k=4:1:m_length fai = -t_y(k-1),-t_y(k-2),t_u(k-1-d),t_u(k-2-d),v(k-1),v(k-2),v(k-3); v(k) = t_y(k) - fai*sita(:,k-1); W = P*fai; dan = 1/(alf + fai*W); sita(:,k) = sita(:,k-1) + dan*W*v(k); P = ( P - dan * W*W)/alf;end仿真结果如下图所示：则各个参数的估计值为：2、 梯度校正辨识方法仿真对象为，其中为服从N(0,1)分布的白噪声，取初值，。所选的辨识模型为，观测数据长度取。则matlab程序如下:%确定性系统的递推梯度校正参数估计（RGC）clear all; close all;a=1 -1.5 0.7; b=1 0.5; d=3; %对象参数na=length(a)-1; nb=length(b)-1; %na、nb为A、B阶次L=400; %仿真长度uk=zeros(d+nb,1); %输入初值：uk(i)表示u(k-i)yk=zeros(na,1); %输出初值u=randn(L,1); %输入采用白噪声序列theta=a(2:na+1);b; %对象参数真值thetae_1=zeros(na+nb+1,1); %参数估计初值alpha=1; %范围(0,2)c=0.1; %修正因子for k=1:L phi=-yk;uk(d:d+nb); y(k)=phi*theta; %采集输出数据 thetae(:,k)=thetae_1+alpha*phi*(y(k)-phi*thetae_1)/(phi*phi+c); %递推梯度校正算法 %更新数据 thetae_1=thetae(:,k); for i=d+nb:-1:2 uk(i)=uk(i-1); end uk(1)=u(k); for i=na:-1:2 yk(i)=yk(i-1); end yk(1)=y(k)