系统辨识作业和答案.doc

系统辨识作业和答案一 问答题1. 介绍系统辨识的步骤。答：（1）先验知识和建模目的的依据；（2）实验设计；（3）结构辨识；（4）参数估计；（5）模型适用性检验。2. 考虑单输入单输出随机系统，状态空间模型转换成ARMA模型。答：ARMA模型的特点是u(k)=0,3. 设有一个五级移位寄存器，反馈取自第2级和第3级输出的模2加法和。试说明：（1） 其输出序列是什么？（2） 是否是M序列？（3） 它与反馈取自第4级与第3级输出模2加法和所得的序列有何不同？（4） 其逆M序列是什么？答：（1）设设输入序列1 1 1 1 1 其输出序列为：1 1 1 1 1 0 0 1 0 1 不是M序列第4级与第3级模2相加结果不同点：第2级和第3级模二相加产生的序列，是从第4时刻开始，每隔7个时刻重复一次；第4级与第3级模2相加产生的,序列，是从第2时刻开始每隔15个时刻重复一次。第5级与第4级模2相加结果如下：已知其为M序列。M序列: 1 1 1 1 1 0 0 0 0 1 0 0 0 1 1 0 0 1 0 1 0 1方波信号: 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0逆重复M: 0 1 0 1 0 0 1 0 1 1 1 0 1 1 0 0 1 1 1 1 1 14. 画出广义最小二乘法的离线迭代算法的简单计算框图。答：广义最小二乘法的离线迭代算法的简单计算框图如下：输入u（k），z（k）设di=0，并计算产生，并计算产生，并计算收敛准则满足否输出5. 考虑如下数学模型，试用输入输出数据估计系统参数。答：6. 利用最小二乘算法辨识如下模型参数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.20558(0.0)显然，辨识结果已经远远偏离了模型参数真值，试从理论上解释为什么会出现这种现象。答：对于n阶系统与n+1阶系统参数估计之间有如下的关系：对于n+1阶系统 设其待估参数为则由题目知n=2时系统参数为准确值，则n=3时按照上式去计算，估算出的系数必远远偏离系统模型参数值。7. 请说明闭环系统不可辨识的原因。答：闭环系统不可辨识的原因：反馈使得一个闭环系统对不同的输入常产生差不多相同的输出，观测的输入输出数据所包含的信息比开环辨识少的多；输入信号与噪声因反馈而相关：有偏估计，非一致性估计；在闭环条件下，用开环辨识方法系统的参数有时也是不可唯一辨识的。8. 设闭环系统前向通道模型为反馈调节器为试画出其闭环系统框图，并判断系统是否可辨识？答：系统是可以辨识的，由于为非奇异，故在条件下，参数是可以辨识的。闭环系统框图如下图所示：9. 对系统模型阶次进行辨识，得到1阶-4阶的参数估计，性能指标与系统模型阶次的关系如下表所示，利用F检验法判断系统模型的阶次。n=1n=2n=3n=4n=5J51.814.6312.4612.4112.40解：由F检验法原理知 若则可以接受系统阶数。由计算得，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序列存入数据文件data3.txt解：（1）function y=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(严晓龙实验：产生一组正态分布的白噪声信号)save data1.txt y -ascii调用函数实验：WNoise(400,0,1)，得到数据见data1.txt，如图所示：（2）function seq=mseq(a,L,N)% a为M序列幅值 N为长度 L为移位单位数register=randint(1,L) %寄存器初始化p=zeros(1,L); %特征向量p(L-1:L)=1; %默认最后两个寄存器相加temp=0;for i=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,ylim,-0.2,1.2);ylabel(M序列)title(严晓龙实验 ：移位寄存器产生的M序列)save data2.txt seq -ascii调用函数实验：mseq(2,40,15)，得到数据见data2.txt，和下图：（3）function seq=invM(a,L,N)% a为M序列幅值 N为长度 L为移位单位数register=randint(1,L) %寄存器初始化p=zeros(1,L); %特征向量p(L-1:L)=1; %默认最后两个寄存器相加temp=0;for i=1:2L-1 seq(i)=register(L); temp=sum(register.*p); register(2:L)=register(1:L-1); %移位 register(1)=mod(temp,2);endseq=seq seq;for i=1:2*(2L-1) if mod(i,2)=1 invm(i)=1; else invm(i)=0; end seq(i)=a*xor(seq(i),invm(i);endfor i=1:N if mod(i,2*(2L-1)=0 mseq(i)=seq(2*(2L-1); else mseq(i)=seq(mod(i,2*(2L-1); endendseq=mseq;stairs(seq);grid;set(gca,ylim,-0.2,1.2);title(严晓龙实验：产生一组逆重复M序列信号)save data3.txt seq -ascii调用函数实验：invM (1,10,40)，得到数据见data3.txt，和下图2. 12.mat中的数据是单输入单输出系统进行采样后100对输入输出数据，其中input表示系统的输入数据，output表示受到噪声污染后的系统的输出数据。在辨识过程中，可以认为噪声具有正态分布，其均值为0。（1） 判断该系统的阶次（方法不限）（2） 利用递推最小二乘法进行参数估计。解：模型阶数的辨识,一般说来低阶模型描述粗糙，高阶模型精度高。残差平方总和J(n)是模型阶数的函数 在不同的模型阶数的假设下，参数估计得到的J(n)值亦不同。讨论如下(1)当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)=output(i);r=100;for p=1:(r-2) %利用循环生成观测矩阵 h(p,:)=-z(p+1) u(p+1); %endhl=h;for b=1:(r-2) %生成输出矩阵 zl(b,:)=z(b+2); zlendzl%根据最小二乘法公式进行参数辩识c1=hl*hl;c2=inv(c1);c3=hl*zl;c=c2*c3;a1=c(1)a2=c(2)j=0;for k=4:100; hl=-z(k-1);u(k-1); x=hl*c; y=z(k)-x; s=y*y; j=j+s;endj仿真结果如下a1 = -0.2576 a2 = 0.6985 j = 0.8556(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;%利用循环生成观测矩阵。for p=1:(r-2) h(p,:)=-z(p+1) -z(p) u(p+1) u(p);endhl=h;%生成输出矩阵。for b=1:(r-2) zl(b,:)=z(b+2); zlendzl%根据最小二乘法公式进行参数辩识c1=hl*hl;c2=inv(c1);c3=hl*zl;c=c2*c3;%输出辩识参数a1=c(1)a2=c(2)b1=c(3)b2=c(4)j=0;%求J（n）for k=4:100; %开始求K h1=-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;endj仿真结果如下:a1 = -0.4362 a2 = 0.2407 b1 = 1.8844 b2 = -1.1313j = 0.5977(3)当n=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;for p=2:(r-1) h(p,:)=-z(p+1) -z(p) -z(p-1) u(p+1) u(p) u(p-1);endhl=h;for b=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=c(2)a3=c(3)b1=c(4)b2=c(5)b3=c(6)j=0;for k=4:100; hl=-z(k-1);-z(k-2);-z(k-3);u(k-1);u(k-2);u(k-3); x=hl*c; y=z(k)-x; s=y*y; j=j+s;endj仿真结果如下a1 = -1.0000 a2 = 0 a3 = 2.1316e-14 b1 = -2.2737e-12b2 = 1.3642e-12 b3 = 1.8190e-12 j = 2.6845数据分析如下利用LS法先对系统参数进行初步辩识并根据其确定残差平方总和J(n) 结果如下：阶数辩识naibiJ (n)1-0.25760.69850.85562-0.43620.24071.8844-1.13130.59773-1.00000