《系统辨识作业2(11页).doc》由会员分享，可在线阅读，更多相关《系统辨识作业2(11页).doc（11页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、-第 1 页系统辨识作业 2-第 2 页系统辨识作业学 院：专 业：姓 名：学 号：日 期：-第-1-页系统辨识作业：系统辨识作业：以下图为仿真对象图中，v(k)为服从 N(0,1)正态分布的不相关随即噪声，输入信号采用循环周期 Np500 的逆 M 序列，幅值为 1，选择辨识模型为：加权因子1)(k，数据长度 L=500，初始条件取IP610)0(，001.0001.0001.0)0(要求：（1）采用一次完成最小二乘法对系统进行辨识，给出数据 u(k)和 z(k)，及LH，LZ和和)(J的值。（2）采用递推最小二乘法进行辨识，要给出参数收敛曲线以及新息)(kZ，残差)(k，准则函数)(kJ随
2、着递推次数 K 的变化曲线。（3）对仿真对象和辨识出的模型进行阶跃响应对比分析以检验辨识结果的实效。1 1、一次完成法对系统进行辨识：、一次完成法对系统进行辨识：估计值LTLLTLLSZHHH1)(，其中一次完成算法对系统辨识的 Matlab 程序见附录:部分输入、输出数据如下，全部的输入输出数据用图 1.1 所示输入数据 u(k)=Columns 1 through 160011110000100110输出数据 z(k)=Columns 1 through 9001.23723.93316.49877.99097.71326.59475.4523Columns 10 through 183.
3、22120.84190.62431.71100.71261.07122.80373.80474.6372图 1.1 输入数据 u(k)和输出数据 z(k)IH的值为(部分)：HL=-5.4523-6.594700-3.2212-5.452300-0.8419-3.22121.00000-3.5706-5.194400-0.6787-3.5706002.3020-0.678700IZ的值为(部分)：ZL=3.22120.8419-第-2-页0.62431.71100.71260.6787-2.3020-3.8270一次完成辨识后的值为：theta=-1.49180.69321.05410.469
4、1)(J的值为：J(theta)=493.18372 2、递推最小二乘法辨识系统：、递推最小二乘法辨识系统：初始条件：IP610)0(，001.0001.0001.0)0(递推最小二乘法辨识系统的 Matlab 程序见附录:其参数收敛曲线如图 2.1图 2.1 参数收敛曲线新息)(kZ，残差)(k，准则函数)(kJ随着递推次数 K 的变化曲线如图 2.2中依次所示:-第-3-页图 2.2 新息、残差、准则函数变化曲线3 3、仿真对象和辨识出的模型进行阶跃响应对比分析、仿真对象和辨识出的模型进行阶跃响应对比分析仿真对象和辨识模型阶跃响应对比 Matlab 程序见附录:图 3.1 分别给出了一次完
5、成算法辨识出来系统和辨识对象的阶跃响应对比图，递推算法辨识出来系统和辨识对象的阶跃响应对比图。附录：一次完成和递推法系统辨识附录：一次完成和递推法系统辨识 MatlabMatlab 程序程序%最小二乘法辨识系统;%叠加噪声为 1/(1-1.5z(-1)+0.7z(-2);%化为差分方程形式为;%e(k)=v(k)+1.5e(k-1)-0.7e(k-2);%仿真对象为(1z-1+0.5z-2)/(1-1.5z-1+0.7z-2);%化为差分方程形式为;%y(k)-1.5y(k-1)+0.7y(k-2)=u(k-1)+0.5u(k-2);%辨识模型为;%z(k)=-a1z(k-1)-a2z(k-2
6、)+b1u(k-1)+b2u(k-2)+v(k);%主函数;function mainclose all;clc;clear;%产生逆 M 序列;-第-4-页X=0,1,1,0,1,0,0,1;%寄存器初始值;F=0,1,1,1,0,0,0,1;%特征多项式;N=1000;%生成长度;M=;%存放产生的 M 序列;%产生逆 M 序列函数调用;out=IMxulie(X,F,N);%阶梯图输出逆 M 序列;figure(1);M=out;subplot(2,1,1);stairs(M);xlabel(k)ylabel(逆 M 序列)title(移位寄存器产生的逆 M 序列);grid on;%一
7、次完成最小二乘法辨识系统;%e(k)=v(k)+1.5e(k-1)-0.7e(k-2);%y(k)=1.5y(k-1)-0.7y(k-2)+u(k-1)+0.5u(k-2);%z(k)=y(k)+e(k);%即 z(k)=1.5y(k-1)-0.7y(k-2)+u(k-1)+0.5u(k-2)+v(k)+1.5e(k-1)-0.7e(k-2);%产生均值为 0，方差为 1 的白噪声;v=;v=randn(1,600);%产生系统辨识所需要数据;y(1)=0;y(2)=0;e(1)=0;e(2)=0;for i=3:600y(i)=1.5*y(i-1)-0.7*y(i-2)+M(i-1)+0.5
8、*M(i-2);e(i)=v(i)+1.5*e(i-1)-0.7*e(i-2);z(i)=y(i)+e(i);endsubplot(2,1,2);stem(z);xlabel(k);-第-5-页ylabel(幅度值);title(输出数据);grid on;%辨识模型为;%z(k)=-a1z(k-1)-a2z(k-2)+b1u(k-1)+b2u(k-2)+v(k);for i=10:509H(i-9,:)=-z(i-1),-z(i-2),M(i-1),M(i-2);endZL=(z(10:509);ES=inv(H*H)*H*ZL;%inv 用来求矩阵的逆矩阵;%输出辨识的参数;disp(输入
9、数据 u(k)=);disp(M);disp(输出数据 z(k)=);disp(z);disp(HL=);disp(H);disp(ZL=);disp(ZL);disp(theta=);disp(ES);%由估计值计算准则函数的值;J=(ZL-H*ES)*(ZL-H*ES);disp(J(theta)=);disp(J);%递推最小二乘法辨识系统;%一式 s(k)=s(k-1)+k(k)(z(k)-h(k)s(k-1);%二式 k(k)=p(k-1)h(k)h(k)p(k-1)h(k)+1/w(k);%三式 p(k)=I-k(k)h(k)p(k-1);%s(0)=0.01,0.01,.0.01
10、;%p(0)=106I;%新息 ZX=z(k)-H(k)S(k-1);%残差 E=z(k)-H(k)S(k);%准则函数 J(k)%递推求解;z1=z(10:509);-第-6-页M1=M(10:509);P=106*eye(4);S(:,1)=0.01;0.01;0.01;0.01;S(:,2)=0.01;0.01;0.01;0.01;J1(1)=0;J1(2)=0;for i=3:500h=-z1(i-1);-z1(i-2);M1(i-1);M1(i-2);K=P*h*inv(h*P*h+1);S(:,i)=S(:,i-1)+K*(z1(i)-h*S(:,i-1);%待估计参数变化;zx(
11、i)=z1(i)-h*S(:,i-1);%新息变化;E(i)=z1(i)-h*S(:,i);%残差变化;J1(i)=J1(i-1)+(zx(i)*zx(i)/(h*P*h+1);%准则函数变化;P=(eye(4,4)-K*h)*P;end%待估计参数过度曲线;figure(2);i=1:500;plot(i,S(1,:),i,S(2,:),i,S(3,:),i,S(4,:);title(待估计参数过度过程);legend(参数 a1 的过渡过程,参数 a2 的过渡过程,参数 b1 的过渡过程,参数 b2的过渡过程);grid on;%disp(S(:,500);%新息变化曲线;figure(3
12、);subplot(3,1,1);i=1:500;plot(i,zx(i);title(新息变化曲线);grid on;%残差变化曲线;subplot(3,1,2);i=1:500;plot(i,E(i);title(残差变化曲线);-第-7-页grid on;%准则函数变化曲线;subplot(3,1,3);i=1:500;plot(i,J1(i);title(准则函数变化曲线);grid on;%disp(J1(500);%分别对一次完成，递推和辨识对象阶跃响应对比;%一次完成法辨识阶跃响应对比;%一次辨识后对仿真对象和辨识模型进行阶跃响应对比;figure(4);subplot(2,1,
13、1);%辨识模型阶跃响应;%z(k)+a1z(k-1)+a2z(k-2)=b1u(k-1)+b2u(k-2);a=1,ES(1),ES(2);b=0,ES(3),ES(4);%求离散系统阶跃响应;gn=dstep(b,a,30);stem(gn,ro);grid on;hold on;%加阶跃响应后辨识对象的系统函数为：%(1+1z-1+0.5z-2)/(1-1.5z-1+0.7z-2);b1=0,1.0,0.5;a1=1,-1.5,0.7;gn1=dstep(b1,a1,30);stem(gn1,bx);hold off;xlabel(k);ylabel(幅值 h(k);title(辨识对象
14、和一次完成法辨识系统阶跃响应对比);legend(估计值,理论值);%递推法辨识阶跃响应对比;subplot(2,1,2);%辨识模型阶跃响应;-第-8-页%z(k)+a1z(k-1)+a2z(k-2)=b1u(k-1)+b2u(k-2);a=1,S(1,500),S(2,500);b=0,S(3,500),S(4,500);%求离散系统阶跃响应;gn=dstep(b,a,30);stem(gn,ro);grid on;hold on;%加阶跃响应后辨识对象的系统函数为：%(1+1z-1+0.5z-2)/(1-1.5z-1+0.7z-2);b1=0,1.0,0.5;a1=1,-1.5,0.7;
15、gn1=dstep(b1,a1,30);stem(gn1,bx);hold off;xlabel(k);ylabel(幅值 h(k);title(辨识对象和递推法辨识系统阶跃响应对比);legend(估计值,理论值);end%子函数，产生周期大于 500 的逆 M 序列函数;function out=IMxulie(X,F,N)n=N;X1=X;F1=F;s=1;for i=1:nY=X1;t=F1(8)*Y(8);for j=7:-1:1t=xor(t,F1(j)*Y(j);endX1(1)=t;for k=2:8X1(k)=Y(k-1);end-第-9-页U(i)=Y(8);V(i)=xor(U(i),s);s=not(s);endout=V;end