试验2相关分析法辨识脉冲响应

1、实验2 相关分析法辨识脉冲响应实验报告哈尔滨工业大学航天学院控制科学与工程系专业：自动化班级：1004102姓名：籍洋日期： 2013年10月10日1实验题目:相关分析法辨识系统脉冲响应2.实验目的通过仿真实验掌握利用相关分析法辨识脉冲响应的原理和方法。3实验主要原理一个单入单出线性定常系统的动态特性可用它的脉冲响应函数g( R)描述。X(t)线性系统旳# g何r tO0则 y(t)二 o g(匚)x(t上式两端同乘x(t-.),进而取时间均值，有1 T:1 TTim t 0 y(t)x(t -)dt = ° g(匚)1帆 丁 ° x(t - ；)x(t - .)dtd匚

2、则R/ )° g(；)Rx( -；)d匚这就是著名的维纳-霍夫积分方程。如果输入是 白噪声，这时x(t)的自相关函数为Rx( ) =k、()，RxC -二)( -匚)则根据维纳-霍夫积分方程可得Rxy( )° g(；)Rx( -Gd ：； -kg( )或者g()二Rxy()k这样，只要记录x(t)、y(t)的值，并计算它们的互相关函数，即可求得脉冲响应函数g( t)而在系统有正常输入的情形下，辨识脉冲响应的原理图如下图所示线牲系统-yod)延迟丁乘法器积分器4.实验对象或参数下图为本实验的原理框图。系统的传递函数为G (s),其中K =120,=8.3Sec,T2 =6.2

3、Sec ； u(k)和z(k)分别为系统的输入和输出变量；v(k)为测量白噪声，服从正态分布，均值为零，方差为:，记作 v(k) N (0,石2 ) ； g0 (k)为系统的脉冲响应理论值，g(k)为系统脉冲响应估计 值，(k)为系统脉冲响应估计误差。系统的输入采用M序列(采用实验1中的M序列即可)，输出受到白噪声v(k)的污染。根据过程的输入和输出数据 占(k), z(k)，利用相关分析法计算?出系统的脉冲响应值g(k)，并与系统的脉冲响应理论值go(k)比较，得到系统 脉冲响应估计误差值(k)，当kT闵时，应该有(k)t o。v(k)go (k)l( k ) = go (k)-g (k)1

4、、模拟过程传递函数G(s)，获得过程的输入和输出数据 (k), z( k)?(采 样时间取1秒)。(1)惯性环节u(k)Ky(k)s +1/T其中，T为惯性环节的时间常数，K为惯性环节的静态放大倍数。若采样时间 记作To，则惯性环节的输出可写成：y(k) = e°0/Ty(k -1) TK(1 -e0/T )u(k -1)+ TKT(e。0T 一1)+ To)】U(k)一U(J"To传递函数G(s)仿真(串联)K 11G (s):T1T2 s+1/s + 10T2K令K1=，则G (s)的表达框图为：T1T2'u(k)K1x(k)1y(k)s +1 / T1s+1

5、/ T22、互相关函数的计算1(r 1)NpRmz (k)=' u(i -k)z(i)rN p Tp 出其中，r为周期数，i二Np T表示计算互相关函数所用的数据是从第二个周期 开始的，目的是等过程仿真数据进入平稳状态。(可分别令r =1、3,对比仿真 结果)3、c的补偿补偿量C应取Rmz ( Np -1)，不能取Rmz (Np )。因为Rmz (k)是周期函数,则有 Rmz (Np )= Rmz (0)，故不能取Rmz ( Np )。4、计算脉冲响应估计值?Ngo(k)- ?(k)脉冲响应估计误差;g_脉冲响应估计值g（ k八（Np.） a2-JRMz（k） C1Np二(g。(k)2

6、k5程序框图开始初始化参数T0,T1,T2,K产生M序列和白噪声v 求出系统在M序列作用下的输出yZ=Y+V1!求Rmz!求估计脉冲响应G1F求理论脉冲响应g0JL计算估计误差绘图结束6 .程序代码fun cti onsigma=resp on se(r)x=0,1,0,1,1,0;%初始化Np=2A6-1;%游列长度a=1;%振幅t=1;for i=1:Np* (叶1)y(i)=x(6);temp=xor(x(5),x(6);for j=5:-1:1x(j+1)=x(j);endx(1)=temp;endfor i=1:Np* (叶1)if (y(i)=0)u(i)=a;elseu(i)=-

7、a;endendK=120;T1=8.3;T2=6.2;T0=1;K1=K/T1/T2;x(1)=0;y(1)=0;for k=2:Np*(r+1)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)=exp(-T0/T2)*y(k-1)+T2*(1-exp(-T0/T2)*x(k-1)+T1* (T 2*(exp(-T0/T2 )-1)+T0)*(x(k)-x(k-1)/T0;%未经白噪声污染的输岀endv=white noise(1,le

8、 ngth(y);%产生白噪声z=y+v;%系统实际输岀for k=1:Npsum=0;for i=Np+1:(r+1)*Npsum=u(i-k)*z(i)+sum;endRmz(k)=1/(r*Np)*sum;endc=-Rmz(Np-1);%补偿量 c% 计算脉冲响应估计值for k=1:NpG(k)=Np/(Np+1)*aA2*t)*(Rmz(k)+c);gO(k)=K/(T1-T2)*(exp(-k*t/T1)-exp(-k*t/T2);end%计算脉冲响应估计误差SUM1=0;SUM2=0;for k=1:Npe(k)=g0(k)-G(k);SUM1=e(k)A2+SUM1;SUM2

9、=g0(k)A2+SUM2;endsigma=sqrt(SUM1/SUM2);step=0:Np-1;plot(step,Rmz(63),Rmz(1:62);hold on;plot(step,G(63),G(1:62),'r')；plot(step,g0(63),g0(1:62),'g')；legend('互相关函数','脉冲响应估计值','脉冲响应理论值')end产生白噪声的函数：lenfunctionsig=whitenoise(sigma,len)%白噪声产生函数，sigma 为均方差,为白噪声序列数据长度

10、a=65539;M=2147483647;b=100;x(1)=12345;r(1)=x(1)/M;%第一部分为产生0-1的均匀分布随机数for i=1:12*lenx(i+1)=mod(a*x(i)+b,M);r(i+1)=x(i+1)/M;end%n = 12;for i = 1 : lensig(i)=0;for j = 1 : nsig(i) = sig(i) + r(n * (i - 1) + j);%第二部分产生正态分布，方差为sigma的随机序列endendsig = (sig - 12*0.5)*sigma;end7 实验结果及分析 另r = 1,白噪声均方差 sigma=0.

11、5,运行命令： sigma=response(1)运行结果如下：得 sigma = 0.0416543210-1-2010 2030405060互相关函数脉冲响应估计值脉冲响应理论值另 r = 3,白噪声均方差 sigma=0.5,运行命令：sigma=response(3,0.5)运行结果如下：得 sigma = 0.0373图像如下：76543210-1互相关函数脉冲响应估计值脉冲响应理论值-2 t010 203040506070比较发现r=1和r=3时产生的曲线基本相似，但是脉冲响应误差在r=3时更小。另r = 1,白噪声均方差 sigma=1,运行命令：sigma= response(1,1)运行结果如下：得 sigma = 0.0594互相关函数脉冲响应估计值脉冲响应理论值76543210-1-2-3 010203040506070比较第一个与第三个结果，发现脉冲响应的估计误差是随着输入白噪声标准