实验一利用相关函数辨识脉冲响应

1、实验一利用相关函数辨识脉冲响应北京工商大学系统辨识课程实验报告（2014-20151学期）课程名称：系统辨识题目：利用相关分析法辨识脉冲响应专业班级：控制工程学生姓名：指导教师：刘刘成绩：2015年1月18日北京工商大学计算机与信息工程学院2014-2015 1学期、实验目的通过仿真实验掌握利用相关分析法辨识脉冲响应的原理和方法。、实验内容图1为本实验的原理框图。过程传递函数为G(s),其中K120,丁8.3Sec,T26.2Sec；u(k)和z(k)分别为过程的输入和输出变量；v(k)为过程测量白噪声，服从正态分布，均值为零，方差为2,记作v(k)N(0,；)；g0(k)为过程的脉冲响应理论

2、值，?(k)为过程脉冲响应估计值，(k)为过程脉冲响应估计误差。过程的输入驱动采用M序列，输出受到白噪声v(k)的污染。根据过程的输?入和输出数据u(k),z(k),利用相关分析算法根据输出过程的脉冲响应值g(k),并与过程脉冲响应理论值g0(k)比较，得到过程脉冲响应估计误差值(k),当k时，应该有(k)0。,9相关分析L(NPn=-烈k)图i相关分析法辨识脉冲响应原理框图三、实验要求进行方案设计，模拟过程传递函数，获得输出数据，用M»列作为辨识的输入信号，噪声采用标准正态分布的白噪声，计算互相关函数，不同值的脉冲响应估计值、脉冲响应理论值和脉冲响应估计误差，计算信噪比，画出实验流

3、程图，用MATLAB程实现。四、实验原理1、采用串联传递函数G（s）仿真G(s)T1T2 s 1/T1s 1/T2令K1=则G（s）的表达框图为:KTjTu(2、一个单输入单输出线性定常系统的动态特性可用它的脉冲响应函数g（7）来描述。y（t）则y(t)0g()x(t)d上式两端同乘x(t),进而取时间均值，有1T1TTimT0y(t)x(t)dt0g()TmT0x(t)x(t)dtd则Rxy()0g()Rx()d这就是著名的维纳霍夫积分方程C如果输入是白噪声，这时x(t)的自相关函数为Rx()k(),Rx()k()则根据维纳霍夫积分方程可得Rxy( )0 g( )Rx()dkg()或者g()

4、Rxy()这样，只要记录x(t)、y(t)的值，并计算它们的互相关函数，即可求得脉冲响应函数g(T)o而在系统有正常输入的情形下，辨识脉冲响应的原理图如下图所示。正常输入五、实验框图产一生MJ产幼扪'屉下歹urk为加角入UI(k)aR"£x(k)美曲V(X)产广1N跳7*v(k)I门口条声心/K)会：力口手中至ijw(k)”一野.他”不自工.£里论"CT*LgOz<r±i_i+i5i2Sdeltag作僮(岳玉；六、实验代码function clc;clear allclose allex2Np=63;%循环周期delta_T=1;

5、%寸钟节拍a=1;%国度M(1)=1;M(2)=0;M(3)=0;M(4)=1;M(5)=1;m(6)=o；金女a化忻列M_XuLie(Np)=0;forn=1:Nptemp=xor(M(6),M(5);if(temp=0)M_XuLie(n)=a;elseM_XuLie(n)=-a;endM(6)=M(5);M(5)=M(4);M(4)=M(3);M(3)=M(2);M(2)=M(1);M(1)=temp;end%生成游列完毕r=3;%周期数u=repmat(M_XuLie,1,r+1);%|各M序列赋给输入，作为输入信号瘫一步,从u(k)得到x(k),y(k)K=120;T0=1;%采样时

6、间T1=8.3;T2=6.2;K1=K/(T1*T2);%初始化X(k),Y(k)为0K2=1x(63)=0;y(63)=0fork=2:63*4%取得x(k)序列x(k)=exp(-T0/T1)*x(k-1)+T1*K1*(1-exp(-T0/T1)*u(k-1)+T1*K1.*(T1*(exp(-T0/T1)-1)+T0)*(u(k)-u(k-1)/T0%取得y(k)序列y(k尸exp(-T0/T2)*y(k-1)+T2*K2*(1-exp(-T0/T2)*x(k-1)+T2*K2.*(T2*(exp(-T0/T2)-1)+T0)*(x(k)-x(k-1)/T0end%获取没有白噪声时候输

7、出完毕%乍图figure(1);plot(u,'r');holdon;plot(x,'k');plot(y,'b');legend('u(k)','x(k)','y(k)');%第二步，将白噪声添加入输出信号%产生白噪声信号vfangcha=0.5;%遁意指定的方差v=fangcha*randn(1,63*4);北京工商大学计算机与信息工程学院2014-2015 1学期%言号叠加,输出实际信号z(k)z=y+v;figure(2);%丁印无白噪声污染信号plot(y,'b');hol

8、don;%丁印白噪声信号plot(v,'m');%丁印白噪声污染后的信号plot(z,'k');legend('y(k)','v(k)','z(k)');%#算Rmz(k)fork=1:NpRmz(k)=0;%初始化为0fori=(Np+1):(r+1)*Np)Rmz(k)=Rmz(k)+u(i-k)*z(i);endRmz(k)=Rmz(k)/(r*Np);end%十算cc=-Rmz(Np-1);10北京工商大学计算机与信息工程学院2014-2015 1学期%计算脉冲响应估计值g1g1=Np*(Rmz+c)/(N

9、p+1)*aA2*delta_T);%计算理论脉冲g0fork=1:Npg0(k户K/(T1-T2)*(exp(-k*delta_T/T1)-exp(-k*delta_T/T2);end%计算脉冲响应估计误差delta_gdelta_g=sqrt(sum(g0-g1).A2)/sum(g0.A2);figure(3);plot(g0,'k');holdon;plot(g1,'r');%axis(0,100,0,10);legend('脉冲响应理论值g0(k)','脉冲响应估计值g1');七、实验结果1、轴入u(k),中间轴入x(k),无干扰轴入(k)2、白噪声标准差为1.5时，理想输出y(k),带干扰的输出z(k),干扰v(k)3、输入白噪声标准差为1.5,周期数r为3时，脉冲响应理论值与估计值脉冲响应估计误差