




版权说明:本文档由用户提供并上传,收益归属内容提供方,若内容存在侵权,请进行举报或认领
文档简介
1、精选优质文档-倾情为你奉上系统辨识与参数估计大作业 第一题 递推最小二乘估计参数考虑如上图所示的仿真对象,选择模型结构为:,其中是服从正态分布的不相关随机噪声;输入信号采用4阶逆M序列,特征多项式取,幅度为1,循环周期为;控制值,使数据的噪信比分别为10,73,100三种情况。加权因子;数据长度L=500;初始条件取,(1) 利用递推最小二乘算法在线估计参数,(2) 利用模型阶次辨识方法(AIC准则),确定模型的阶次。(3) 估计噪声的方差和模型静态增益(4) 作出参数估计值随时间的变化图答:设过程的输入输出关系可以描述成是输出量,是可观测的数据向量,是均值为0的随机噪声选取的模型为结构是加权
2、最小二乘参数估计递推算法RWLS的公式如下,为了把p(k)的对称性,可以把p(k)写成如果把设成1的时候,加权最小二乘法就退化成最小二乘法。用AIC准则定阶法来定阶,所用公式其中模型参数和 噪声方差的极大似然估计值为 ,AIC的定阶公式写成取分别计算,找到使最小的那个作为模型的阶次。一般说来,这样得到的模型阶次都能比较接近实际过程的真实阶次。信噪比为10%时:参数a1a2b1b2噪声方差静态增益模型阶次真值-1.50.710.51估计值-1.5190.722591.03140.509231.09517.56612信噪比为73%时:参数a1a2b1b2噪声方差静态增益模型阶次真值-1.50.71
3、0.51估计值-1.5190.722591.03140.509231.09517.56612信噪比为100%时:参数a1a2b1b2噪声方差静态增益模型阶次真值-1.50.710.51估计值-1.5190.722591.03140.50923估计值7.56612源程序:%function a1 a2 b1 b2 na nb fangcha Kk=rwls(L,syn,Np) %na,nb为模型阶次,fangcha为噪声方差,Kk为静态增益a1=0;a2=0;b1=0;b2=0;na=0;nb=0;fangcha=0;Kk=0;L=500;Np=62;syn=1;x(1:4)=1 0 1 0;f
4、or i=1:Np temp=xor(x(1),x(4); M(i)=x(4); for j=4:-1:2 x(j)=x(j-1); end x(1)=temp;endS=ones(1,Np);%先产生一个全是1的序列Sif mod(Np,2)=0%判断Np是奇数还是偶数 p=Np/2;else p=(Np-1)/2;endfor j=1:p S(2*j)=0;%将S序列的偶数位值均置为0,从而使S序列是0或1的方波序列end IM=xor(M,S); %使用M序列与方波序列S复合生成逆M序列IMu=IM*2-1;for i=(Np+1):L u(i)=u(i-Np);endrandn(
5、9;seed',2);v=randn(1,L);syms c;y(1)=0;y(2)=u(1);e(1)=c*v(1);e(2)=1.5*e(1)+c*v(2);for i=3:L y(i)=1.5*y(i-1)-0.7*y(i-2)+u(i-1)+0.5*u(i-2); e(i)=1.5*e(i-1)-0.7*e(i-2)+c*v(i);endm=sum(e.2);n=sum(y.2);%c=solve('m/n-syn*syn','c');c=solve('m/n-syn*syn');c1=abs(double(c(1);z(1)=c
6、1*v(1);z(2)=u(1)+1.5*z(1)+c1*v(2);for i=3:L z(i)=1.5*z(i-1)-0.7*z(i-2)+u(i-1)+0.5*u(i-2)+c1*v(i);endcta=zeros(4,L);cta(:,1)=0.001 0.001 0.001 0.001'P=diag(106 106 106 106);%k=2时h=-z(1) 0 u(1) 0'K= P*h*inv(h'*P*h+1);P=P-K*K'*(h'*P*h+1);cta(:,2)=cta(:,1)+K*(z(2)-h'*cta(:,1);for
7、 k=3:L h=-z(k-1) -z(k-2) u(k-1) u(k-2)' K=P*h*inv(h'*P*h+1); P=P-K*K'*(h'*P*h+1); cta(:,k)=cta(:,(k-1)+K*(z(k)-h'*cta(:,(k-1);end %以上为参数估计值 z=z' for na=1:4 for nb=1:4 A=zeros(L,na); B=zeros(L,nb);for i=1:L for j=1:na if i>j A(i,j)=z(i-j); end endendfor i=1:L for j=1:nb if
8、i>j B(i,j)=u(i-j); end endendH=A B;cta1=inv(H'*H)*H'*z;%cta1为模型参数极大似然估计值cgm(na,nb)=(z-H*cta1)'*(z-H*cta1)/L;%cgm为噪声方差极大似然估计值AIC(na,nb)=L*log(cgm(na,nb)+2*(na+nb); endendna nb=find(AIC=min(min(AIC);fangcha=cgm(na,nb)/(c12);a1=cta(1,500);a2=cta(2,500);b1=cta(3,500);b2=cta(4,500);Kk=(cta
9、(3,L)+cta(4,L)/(1+cta(1,L)+cta(2,L);m=1:L;plot(m,cta(1,:),'b-',m,cta(2,:),'k-',m,cta(3,:),'y-',m,cta(4,:),'r-')第二大题 卡尔曼滤波一个系统模型为同时有下列条件:(1) 初始条件已知且有。(2) 是一个标量零均值白高斯序列,且自相关函数已知为,另外,我们有下列观测模型,即,且有下列条件:(3) 和是独立的零均值白高斯序列,且有 ,(4) 对于所有的和,与观测噪声过程和是不相关的,即,我们希望得到由观测矢量估计状态矢量的卡尔
10、曼滤波的公式表示,并求解以下问题:(a) 求出卡尔曼增益矩阵,并得出最优估计和观测矢量之间的递推关系。(b) 用模拟数据确定状态矢量的估计值,并画出当时,的图(c) 通常,状态矢量的真实值是得不到的,但是为了用作图来说明问题,表1和表2给出了状态矢量元素的真实值。对于,在同一幅图中画出真实值和在(b)中确定的的估计值。对重复这一过程。当从1变到10时,对每个元素,计算并画出各自的误差图,即。(d) 当从1变到10时,通过用由卡尔曼滤波器决定的状态误差协方差矩阵画出和,而(e) 讨论一下(C)中计算的误差与(d)中方差之间的关系。 表1观测值时间下标观测值观测值13223037349251126
11、15672258283930510384表2由模拟得到的实际状态值时间下标实际状态值实际状态值00.0.11.1.23.1.35.2.49.3.512.3.616.4.721.4.825.4.931.5.1036.5.答:(a)卡尔曼增益矩阵:估计值与观测值之间的递归关系为:(b)状态矢量估计值的计算框图:+(c)和的图:(d)真实值与估计值的比较图:各自的误差图:(e)通过用卡尔曼滤波器的状态误差协方差矩阵画出的和:(f)分析:(e)中的方差是(d)中的误差平方后取均值,是均方误差。误差直接由真实值减去估计值,有正有负,而均方误差没有这个缺陷,更能综合的表示滤波的效果。源程序:%P K X
12、Z均为2*2矩阵,用三维矩阵表示他们,记录每一个k值下的情况,如P(:,:,k)%k从1开始,k=1时取初值clear all;P(:,:,1)=10 0;0 10;%P(0|0)初值,P表示P(k|k)X(:,:,1)=0;0 ; %X(0)初值fy= 1 1;0 1; %tao=1;1;%H= 1 0;0 1;Q=1;R=1 0;0 2;Z(:,:,2)=3.;2.;Z(:,:,3)=3.;0.;Z(:,:,4)=7.;3.;Z(:,:,5)=9.;2.;Z(:,:,6)=11.;2.;Z(:,:,7)=15.;6.;Z(:,:,8)=22.;5.;Z(:,:,9)=28.;3.;Z(:,
13、:,10)=30.;5.;Z(:,:,11)=38.;4.;for k=1:10 P1(:,:,k)=fy*P(:,:,k)*fy'+tao*Q*tao' %P1(:,:,k)表示P(k+1|k) K(:,:,k+1)=P1(:,:,k)*H'*inv(H*P1(:,:,k)*H'+R); X(:,:,k+1)=fy*X(:,:,k)+K(:,:,k+1)*(Z(:,:,k+1)-H*fy*X(:,:,k); P(:,:,k+1)=(eye(2,2)-K(:,:,k+1)*H)*P1(:,:,k);endgx1(1:11)=X(1,1,1:11); %gx1,g
14、x2分别表示X的估计值gx2(1:11)=X(2,1,1:11);plot(0:10,gx1,0:10,gx2);legend('x1估计值','x2估计值','location','northwest');title('状态矢量X的估计值x1(k)和x2(k)');xlabel('k');ylabel('X');x1=0,1.,3.,5.,9.,12.,. 16.,21.,25.,31.,36.;x2=0,1.,1.,2.,3.,3.,. 4.,4.,4.,5.,5.;figure
15、;plot(0:10,gx1,0:10,x1,0:10,gx2,0:10,x2);legend('x1估计值','x1真实值','x2估计值','x2真实值','location','northwest');title('状态矢量X的估计值与真实值之间的关系');xlabel('k');ylabel('X');delta1=x1-gx1;%delta1,delta2分别表示真值与估计值的差值delta2=x2-gx2;figure;plot(1:10,delta1(2:11),1:10,delta2(2:11);legend('x1的误差','x2的误差','location','northwest');title('状态矢量X
温馨提示
- 1. 本站所有资源如无特殊说明,都需要本地电脑安装OFFICE2007和PDF阅读器。图纸软件为CAD,CAXA,PROE,UG,SolidWorks等.压缩文件请下载最新的WinRAR软件解压。
- 2. 本站的文档不包含任何第三方提供的附件图纸等,如果需要附件,请联系上传者。文件的所有权益归上传用户所有。
- 3. 本站RAR压缩包中若带图纸,网页内容里面会有图纸预览,若没有图纸预览就没有图纸。
- 4. 未经权益所有人同意不得将文件中的内容挪作商业或盈利用途。
- 5. 人人文库网仅提供信息存储空间,仅对用户上传内容的表现方式做保护处理,对用户上传分享的文档内容本身不做任何修改或编辑,并不能对任何下载内容负责。
- 6. 下载文件中如有侵权或不适当内容,请与我们联系,我们立即纠正。
- 7. 本站不保证下载资源的准确性、安全性和完整性, 同时也不承担用户因使用这些下载资源对自己和他人造成任何形式的伤害或损失。
最新文档
- 企业绩效评估的指标体系试题及答案
- 2024人力资源管理师考试知识提示试题及答案
- 2024年福建三支一扶考试相关试题及答案
- 2024年高中地理全球性问题探讨试题及答案
- 常见计算机应用案例试题及答案
- 2025年-海南省建筑安全员B证考试题库及答案
- Photoshop数字影像处理案例教程 习题及答案 Chapter 7 产品颜色处理
- 2024年陪诊师考试实务案例试题及答案
- 如何备战2024年陪诊师考试的试题及答案
- 2024年马工学管理的实践与挑战试题与答案
- MOOC 写作与表达-常熟理工学院 中国大学慕课答案
- 农贸市场应急预案
- 肥胖患者麻醉管理专家共识2023年版中国麻醉学指南与专家共识
- 中药饮片处方点评表-副本(文档良心出品)
- DL-T 5605-2021太阳能热发电厂蒸汽发生系统设计规范-PDF解密
- 学校双重预防体系建设指导书
- 螺蛳粉出口贸易的现状及策略分析
- 2024年江苏省苏州市中考数学一模练习卷
- 家政聘用合同模板
- 防汛应急预案培训课件
- 桩基施工安全培训课件
评论
0/150
提交评论