實(shí)驗(yàn)一利用相關(guān)函數(shù)辨識(shí)脈沖響應(yīng)

1、實(shí)驗(yàn)一利用相關(guān)函數(shù)辨識(shí)脈沖響應(yīng)北京工商大學(xué)系統(tǒng)辨識(shí)課程實(shí)驗(yàn)報(bào)告（2014-20151學(xué)期）課程名稱(chēng)：系統(tǒng)辨識(shí)題目：利用相關(guān)分析法辨識(shí)脈沖響應(yīng)專(zhuān)業(yè)班級(jí)：控制工程學(xué)生姓名：指導(dǎo)教師：劉劉成績(jī)：2015年1月18日北京工商大學(xué)計(jì)算機(jī)與信息工程學(xué)院2014-2015 1學(xué)期、實(shí)驗(yàn)?zāi)康耐ㄟ^(guò)仿真實(shí)驗(yàn)掌握利用相關(guān)分析法辨識(shí)脈沖響應(yīng)的原理和方法。、實(shí)驗(yàn)內(nèi)容圖1為本實(shí)驗(yàn)的原理框圖。過(guò)程傳遞函數(shù)為G(s),其中K120,丁8.3Sec,T26.2Sec；u(k)和z(k)分別為過(guò)程的輸入和輸出變量；v(k)為過(guò)程測(cè)量白噪聲，服從正態(tài)分布，均值為零，方差為2,記作v(k)N(0,；)；g0(k)為過(guò)程的脈沖響應(yīng)理論

2、值，?(k)為過(guò)程脈沖響應(yīng)估計(jì)值，(k)為過(guò)程脈沖響應(yīng)估計(jì)誤差。過(guò)程的輸入驅(qū)動(dòng)采用M序列，輸出受到白噪聲v(k)的污染。根據(jù)過(guò)程的輸?入和輸出數(shù)據(jù)u(k),z(k),利用相關(guān)分析算法根據(jù)輸出過(guò)程的脈沖響應(yīng)值g(k),并與過(guò)程脈沖響應(yīng)理論值g0(k)比較，得到過(guò)程脈沖響應(yīng)估計(jì)誤差值(k),當(dāng)k時(shí)，應(yīng)該有(k)0。,9相關(guān)分析L(NPn=-烈k)圖i相關(guān)分析法辨識(shí)脈沖響應(yīng)原理框圖三、實(shí)驗(yàn)要求進(jìn)行方案設(shè)計(jì)，模擬過(guò)程傳遞函數(shù)，獲得輸出數(shù)據(jù)，用M»列作為辨識(shí)的輸入信號(hào)，噪聲采用標(biāo)準(zhǔn)正態(tài)分布的白噪聲，計(jì)算互相關(guān)函數(shù)，不同值的脈沖響應(yīng)估計(jì)值、脈沖響應(yīng)理論值和脈沖響應(yīng)估計(jì)誤差，計(jì)算信噪比，畫(huà)出實(shí)驗(yàn)流

3、程圖，用MATLAB程實(shí)現(xiàn)。四、實(shí)驗(yàn)原理1、采用串聯(lián)傳遞函數(shù)G（s）仿真G(s)T1T2 s 1/T1s 1/T2令K1=則G（s）的表達(dá)框圖為:KTjTu(2、一個(gè)單輸入單輸出線性定常系統(tǒng)的動(dòng)態(tài)特性可用它的脈沖響應(yīng)函數(shù)g（7）來(lái)描述。y（t）則y(t)0g()x(t)d上式兩端同乘x(t),進(jìn)而取時(shí)間均值，有1T1TTimT0y(t)x(t)dt0g()TmT0x(t)x(t)dtd則Rxy()0g()Rx()d這就是著名的維納霍夫積分方程C如果輸入是白噪聲，這時(shí)x(t)的自相關(guān)函數(shù)為Rx()k(),Rx()k()則根據(jù)維納霍夫積分方程可得Rxy( )0 g( )Rx()dkg()或者g()

4、Rxy()這樣，只要記錄x(t)、y(t)的值，并計(jì)算它們的互相關(guān)函數(shù)，即可求得脈沖響應(yīng)函數(shù)g(T)o而在系統(tǒng)有正常輸入的情形下，辨識(shí)脈沖響應(yīng)的原理圖如下圖所示。正常輸入五、實(shí)驗(yàn)框圖產(chǎn)一生MJ產(chǎn)幼捫'屜下歹urk為加角入U(xiǎn)I(k)aR"£x(k)美曲V(X)產(chǎn)廣1N跳7*v(k)I門(mén)口條聲心/K)會(huì)：力口手中至ijw(k)”一野.他”不自工.£里論"CT*LgOz<r±i_i+i5i2Sdeltag作僮(岳玉；六、實(shí)驗(yàn)代碼function clc;clear allclose allex2Np=63;%循環(huán)周期delta_T=1;

5、%寸鐘節(jié)拍a=1;%國(guó)度M(1)=1;M(2)=0;M(3)=0;M(4)=1;M(5)=1;m(6)=o；金女a(chǎn)化忻列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;%周期數(shù)u=repmat(M_XuLie,1,r+1);%|各M序列賦給輸入，作為輸入信號(hào)癱一步,從u(k)得到x(k),y(k)K=120;T0=1;%采樣時(shí)

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%獲取沒(méi)有白噪聲時(shí)候輸

7、出完畢%乍圖figure(1);plot(u,'r');holdon;plot(x,'k');plot(y,'b');legend('u(k)','x(k)','y(k)');%第二步，將白噪聲添加入輸出信號(hào)%產(chǎn)生白噪聲信號(hào)vfangcha=0.5;%遁意指定的方差v=fangcha*randn(1,63*4);北京工商大學(xué)計(jì)算機(jī)與信息工程學(xué)院2014-2015 1學(xué)期%言號(hào)疊加,輸出實(shí)際信號(hào)z(k)z=y+v;figure(2);%丁印無(wú)白噪聲污染信號(hào)plot(y,'b');hol

8、don;%丁印白噪聲信號(hào)plot(v,'m');%丁印白噪聲污染后的信號(hào)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北京工商大學(xué)計(jì)算機(jī)與信息工程學(xué)院2014-2015 1學(xué)期%計(jì)算脈沖響應(yīng)估計(jì)值g1g1=Np*(Rmz+c)/(N

9、p+1)*aA2*delta_T);%計(jì)算理論脈沖g0fork=1:Npg0(k戶K/(T1-T2)*(exp(-k*delta_T/T1)-exp(-k*delta_T/T2);end%計(jì)算脈沖響應(yīng)估計(jì)誤差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('脈沖響應(yīng)理論值g0(k)','脈沖響應(yīng)估計(jì)值g1');七、實(shí)驗(yàn)結(jié)果1、軸入u(k),中間軸入x(k),無(wú)干擾軸入(k)2、白噪聲標(biāo)準(zhǔn)差為1.5時(shí)，理想輸出y(k),帶干擾的輸出z(k),干擾v(k)3、輸入白噪聲標(biāo)準(zhǔn)差為1.5,周期數(shù)r為3時(shí)，脈沖響應(yīng)理論值與估計(jì)值脈沖響應(yīng)估計(jì)誤差