系统辨识作业及答案

1、一.问答题1 .介绍系统辨识的步骤。答：(1)先验知识和建模目的的依据；(2)实验设计；(3)结构辨识；(4)参数估计；(5)模型适用性检验。2 .考虑单输入单输出随机系统，状态空间模型101x(k1)=20x(k)1u(k)y(k)=11X(k)v(k)转换成ARMA模型。答：ARMA模型的特点是u(k)=0,10x(k1)x(k)20y(k)=11x(k)v(k)3.设有一个五级移位寄存器，反馈取自第2级和第3级输出的模2加法和。试说明:(1) 其输出序列是什么？(2) 是否是M序列？(3) 它与反馈取自第4级与第3级输出模2加法和所得的序列有何不同？(4) 其逆M序列是什么？答：(1)设

2、设输入序列11111CD1111190111017001112510011(2)01111100011118100112601001(3)00111111001119010012710100(4)1001112010012010100281101001001(13)1010021110102900111(6)10100141101022111013001110110101511101(23)01110310011111101160111024001113210011其输出序列为：1111100101不是M序列第4级与第3级模2相加结果CD1111191100117011112501100011

3、11100110018001112610110(300111111011019000112701011(4)0001112010112010001281010110001(13)1010121010002911010(6)01000141101022001003011101001001511101(23)10010311111010010161111024110013201111不同点：第2级和第3级模二相加产生的序列，是从第4时刻开始，每隔7个时刻重复一次;15个时刻重复一次。第4级与第3级模2相加产生的，序列，是从第2时刻开始每隔第5级与第4级模2相加结果如下：已知其为M序歹U。11111

4、9000100111110100011701010(30011111110001810101(4)000111201100191101000001(13)001102011101(6)10000141001121111100100015010012211111001001610100M序列：方波信号：逆重复M:1111101001011000101000100100010101111014 .画出广义最小二乘法的离线迭代算法的简单计算框图。答：广义最小二乘法的离线迭代算法的简单计算框图如下:5 .考虑如下数学模型y=asinx+bcosx,试用x(k),y(k),k=1,2,3,.,N输入输出

5、数据估计系统参数a,bo答：y1=asinx1bcosx1y2=asinx2bcosx2ayN=asinxNbcosxN二Y-二E二.I-(7-)-1-TY6 .利用最小二乘算法辨识如下模型参数z(k)-1.5z(k-1)+0.7z(k-2)=u(k-1)+0.5u(k-2)+v(k)其中，v(k)是零均值白噪声。当模型阶次为2时，可以获得准确的辨识结果，而模型阶次取3时，只能得到如下一组模型参数辨识结果(括号内为模型参数真值)：a1=-1.08884(-1.5)a2=0.08326(0.7)a3=0.28781(0.0)b1=1.00000(1.0)b2=0.91116(0.5)b3=0.2

6、0558(0.0)显然，辨识结果已经远远偏离了模型参数真值，试从理论上解释为什么会出现这种现象。答：又于n阶系统与n+1阶系统参数估计之间有如下的关系：对于n+1阶系统Az1y(k)=Bzuke(k)设其待估参数为Tn1：-|boaibi.an4%dbn.1-丁(2)则u(1)=u(n)-A/Y-中Mn)由题目知n=2时系统参数为准确值，则n=3时按照上式去计算，估算出的系数必远远偏离系统模型参数值。7 .请说明闭环系统不可辨识的原因。答：闭环系统不可辨识的原因：反馈使得一个闭环系统对不同的输入常产生差不多相同的输出，观测的输入输出数据所包含的信息比开环辨识少的多；输入信号与噪声因反馈而相关：

7、有偏估计，非一致性估计；在闭环条件下，用开环辨识方法系统的参数有时也是不可唯一辨识的。8 .设闭环系统前向通道模型为y(k)-1.4y(k-1)-0.45y(k-2)u(k-1)0.7u(k-2)；(k)反馈调节器为u(k)=y(k)0.2y(k-1)试画出其闭环系统框图，并判断系统是否可辨识？答：系统是可以辨识的，由于为非奇异，故在u(k)=y(k)+0.2y(k-1)条件下，参数是可以辨识的。闭环系统框图如下图所示：t(n,n1)二J(n)-J(n1)J(n1)9.对系统模型阶次进行辨识，得到1阶-4阶的参数估计，性能指标与系统模型阶次的关系如下表所示，利用F检验法判断系统模型的阶次。n=

8、1n=2n=3n=4n=5J51.814.6312.4612.4112.40解：由F检验法原理知若t(n,n三3.09则可以接受系统阶数。由计算得，t(1,2)=4.13,t(2,3)=0.49,t(3,4)=0.0034,t(4,5)=0所以系统的阶数为3。二.编程题1. (1)编程产生一组正态分布的白噪声信号，它的均值和方差以及长度可随意调整，将产生的白噪声信号存入数据文件data1.txt(2)编程产生一组M序列信号，它的幅值和长度可随意调整，将产生的M序列存入数据文件data2.txt(3)编程产生一组逆重复M序列信号，它的幅值和长度可随意调整，将产生的逆重复M序列存入数据文件data

9、3.txt解：(1)functiony=WNoise(N,E,VAR)%N为长度E为均值VAR为方差y=randn(1,N);y=y-mean(y);y=y/std(y);y=E+sqrt(VAR)*y;plot(y)title('严晓龙实验：产生一组正态分布的白噪声信号)savedata1.txty-ascii调用函数实验：WNoise(400,0,1),得到数据见data1.txt,如图所示:UFigure1".FileEditViewnsertToolsDesktopWindowHelpdlHIQlW装国口国IQ严晓龙实睑，产生一组正态分布的白噪声信号：1iij.t.J

10、I."0501C0150200度。打。350400(2) functionseq=mseq(a,L,N)%a为M序列幅值N为长度L为移位单位数register=randint(1,L)%寄存器初始化p=zeros(1,L);%特征向量p(L-1:L)=1;%默认最后两个寄存器相加temp=0;fori=1:Nseq(i)=a*register(L);temp=sum(register.*p);register(2:L)=register(1:L-1);%移位register(1)=mod(temp,2);endx=0:1/5:1.2;stairs(seq);grid;set(gca,

11、'ylim',-0.2,1.2);ylabel('M序列')title('严晓龙实验：移位寄存器产生的M序列')savedata2.txtseq-ascii调用函数实验：mseq(2,40,15),得到数据见data2.txt,和下图:RieEditViewInsertToolsDesktopWindowHelpdd会昵Xb3I0严晓龙实验：移位寄存器产生的M序列251oI0.武当IrB1.i'Ih-.1-J051015(3) functionseq=invM(a,L,N)%a为M序列幅值N为长度L为移位单位数register=randi

12、nt(1,L)%寄存器初始化p=zeros(1,L);%特征向量p(L-1:L)=1;%默认最后两个寄存器相加temp=0;fori=1:2AL-1seq(i)=register(L);temp=sum(register.*p);register(2:L)=register(1:L-1);%移位register(1)=mod(temp,2);endseq=seqseq;fori=1:2*(2AL-1)ifmod(i,2)=1invm(i)=1;elseinvm(i)=0;endseq(i)=a*xor(seq(i),invm(i);endfori=1:Nifmod(i,2*(2AL-1)=0m

13、seq(i)=seq(2*(2AL-1);elsemseq(i)=seq(mod(i,2*(2AL-1);endendseq=mseq;stairs(seq);grid;set(gca,'ylim',-0.2,1.2);title('严晓龙实验：产生一组逆重复M序列信号,)savedata3.txtseq-ascii调用函数实验：invM(1,10,40),得到数据见data3.txt,和下图2.12.mat中的数据是单输入单输出系统进行采样后100对输入输出数据，其中input表示系统的输入数据，output表示受到噪声污染后的系统的输出数据。在辨识过程中，可以认为噪

14、声具有正态分布，其均值为0。(1) 判断该系统的阶次(方法不限)(2) 利用递推最小二乘法进行参数估计。解:J是模型阶数的辨识，一般说来低阶模型描述粗糙，高阶模型精度高。残差平方总和NT9模型阶数的函数J(n)八.(y(k)-二三)2kX-在不同的模型阶数的假设下，参数估计得到的J(n)值亦不同。讨论如下当n=1时程序如下：启动matlab,打开12.mat;运行下面程序u=zeros(100,1);%构造输入矩阵z=zeros(100,1);%构造输出矩阵i=1:1:100;u(i,1)=input(i);z=zeros(100,1);%构造输出矩阵i=1:1:100;z(i,1)=outp

15、ut(i);r=100;forp=1:(r-2)%利用循环生成观测矩阵h(p,:)=-z(p+1)u(p+1);%endhl=h;forb=1:(r-2)%生成输出矩阵zl(b,:)=z(b+2);zl'endzl'%根据最小二乘法公式进行参数辩识c1=hl'*hl;c2=inv(c1);c3=hl'*zl;c=c2*c3;a1=1)a2=c(2)j=0;fork=4:100;hl=-z(k-1);u(k-1)'x=hl*c;y=z(k)-x;s=y*y;j=j+s;end仿真结果如下a1=-0.2576a2=0.6985j=0.8556al-0,257

16、6a2=0.69851 =0.855K(2)当n=2时程序如下(输入输出数据同上，只给出不同于一阶系统的程序不同之处)其中U、Z分别是作业要求给出得的输入输出，数据输入同上。启动matlab,打开12.mat;运行下面程序u=zeros(100,1);%构造输入矩阵z=zeros(100,1);%构造输出矩阵i=1:1:100;u(i,1)=input(i);z=zeros(100,1);%构造输出矩阵i=1:1:100;z(i,1)=output(i);r=100;%利用循环生成观测矩阵。forp=1:(r-2)h(p,:)=-z(p+1)-z(p)u(p+1)u(p);endhl=h;%生

17、成输出矩阵。forb=1:(r-2)zl(b,:)=z(b+2);zl'endzl'%根据最小二乘法公式进行参数辩识c1=hl'*hl;c2=inv(c1);c3=hl'*zl;c=c2*c3;%输出辩识参数a1=31)a2=c(2)b1=c(3)b2=c(4)j=0；%求(n)fork=4:100;%开始求Kh1=-z(k-1),-z(k-2),u(k-1),u(k-2)'x=h1'*c;y=z(k)-x;s=y*y;j=j+s;end仿真结果如下：a1=-0.4362a2=0.2407b1=1.8844b2=-1.1313j=0.5977当n

18、=3时程序如下启动matlab,打开12.mat;运行下面程序u=zeros(100,1);%构造输入矩阵z=zeros(100,1);%构造输出矩阵i=1:1:100;u(i,1)=input(i);z=zeros(100,1);%构造输出矩阵i=1:1:100;z(i,1)=output(i);r=100;forp=2:(r-1)h(p,:)=-z(p+1)-z(p)-z(p-1)u(p+1)u(p)u(p-1);endhl=h;forb=2:(r-1)zl(b,:)=z(b+1);zl'endzl'c1=hl'*hl;c2=inv(c1);c3=hl'*zl;c=c2*c3;a1=c(1)a2=32)a3=c(3)b1=c(4)b2=c(5)b3=c(6)j=0；