系统辨识之经典辨识法Word文档格式.docx
《系统辨识之经典辨识法Word文档格式.docx》由会员分享,可在线阅读,更多相关《系统辨识之经典辨识法Word文档格式.docx(22页珍藏版)》请在冰豆网上搜索。
−1(𝑡
)
𝑎
𝑑
𝑡
(1-2)面积法原则上可以求出n为任意阶的个系数。
以n为3为例。
有:
3𝑦
)𝑑
2𝑦
{
|𝑡
→∞=
→∞=0(1-3)
)|𝑡
→∞=1
将式
(1)中的y(t)移至右边,在[0,t]上积分,得
3
(1-4)定义:
𝐹
1(𝑡
)=∫0𝑡
[1−𝑦
)]𝑑
(1-5)由式(1-3)条件可知,当t→∞时,
(1-6)
同理,定义
𝐹
2
(1-7)
由式(1-,3)条件可知,当t→∞时,
(1-8)
因此,可得
[𝐹
−1(𝑡
)−𝑎
−1𝑦
)]dt(1-9)
𝑎
=𝐹
(∞)(1-10)
②当系统的传递函数存在零点时,传递函数如下:
G(s)=k
bsmmn+basmn-1-1smn-1-1++LL++asbs11+1+1,(nm)(1-11)
asn+
其中,Kh=()/U0定义
1
G(s)=K
P(s)
其中,
P(s)=
bsasnmn++basmn-1-1smn-1-1++LL++asbs11+1+1=+1i=1Csii(1-12)m
根据[1−h*(t)]的Laplace变换,求出一阶面积A1,确定L[h(*1t]),并定义二
阶面积A2,以此类推,得到i阶面积Ai。
进一步利用e−st拉氏变换,得到
L[1−h*(t])=
Msii,进而得到Ai的值:
i=0
Ai=01−h*(t)(i1)!
−−t)i−1dt+tj−=20Ai−−j101−h*(t)
(−j!
t)jdt(1-13)(
根据ACi=i,可得:
+𝑎
−1+⋯+𝑎
+1
=(𝑏
𝑚
+𝑏
−1+⋯+𝑏
+1)(1+∑∞𝑖
=1𝐴
𝑖
)。
比较上式两边s的各次幂,便可得到a,b,A之间的关系,如下:
b1AnAn−1LAnm−+1−1An+1
b2An+1AnLAnm−+2An+2=−
MLLLLM
bmAnm+−1Anm+−2LAnAnm+
b1
a1110LL0000bM2+AA12(1-14)
a2=A1
MLLLLLM
bmAn
anAn−1An−2LA110
由此可知,根据式(1-12)、(1-13)、(1-14)便可得到辨识传递函数的参数a,b。
1.2实验过程1.2.1无零点模型系统假设系统的传递函数模型为G(s)=
21,为无零点的模型,利用
10𝑠
+6.5𝑠
+1
Matlab编程,分别在没有噪声和有噪声两种情况下进行辨识,比较辨识结果。
1.没有噪声时,程序如下:
clear;
%==================获得原传递函数方程=======================%num=[1];
den=[106.51];
%=====================产生阶跃采样序列======================%T=0.2;
%采样周期t=0:
T:
30;
%采样时间L=length(t);
%数据长度h=step(num,den,t);
%原传递函数的阶跃响应
K=h(L)%系统增益
%======================面积法求解参数======================%s1=0;
fori=1:
Ls1=s1+(1-h(i))*T;
F(i)=s1;
enda1=s1;
s2=0;
fori=1:
Ls2=s2+(F(i)-a1*h(i))*T;
enda2=s2;
num1=[1];
den1=[a2a11];
disp('
原传递函数为:
'
)
G1=tf(num,den)
通过辨识得到的传递函数为:
G2=tf(num1,den1)
%=============原传递函数和辨识函数的阶跃响应对比图=============%step(G1,'
b-'
G2,'
r-.'
title('
原系统与辨识后所得到系统阶跃响应对比'
)legend('
原响应曲线'
'
辨识响应曲线'
(1)当采样周期T=0.2秒,采样时间t=30s时,行程序后得到原传递函数G1和辨识得到的传递函数G2如图1.1:
图1.1
原系统和辨识后系统的阶跃响应对比图如下:
图1.2
(2)当采样周期T=0.2秒,采样时间t=50s时,行程序后得到原传递函数G1和辨识得到的传递函数G2如下:
图1.3
(3)当采样周期T=0.02秒,采样时间t=50s时,行程序后得到原传递函数G1和辨识得到的传递函数G2如下:
图1.5
2.有噪声的情况下,系统程序如下:
主程序还是用面积法,在程序中添加以下代码:
%产生期望为0,方差为0.01的噪声
figure
(1)w=randn(1,L);
%建立服从正态分布的随机矩阵。
w=w/std(w);
w=w-mean(w);
qw=0;
fc=0.01;
w=qw+sqrt(fc)*w;
%=====================阶跃采样序列中加入白噪声==================%h=h+w;
plot(t,w);
(1)加入的噪声如下图所示:
图1.7
(2)当采样周期T=0.02s,采样时间t=50s时,辨识结果如下:
图1.8
原系统与辨识系统阶跃响应如图所示:
结合上述无测量噪声和有测量噪声两种情况下的辨识结果,列出如下所示的表格:
表1-1
噪声情况
条件
增益
a1
a2
参数
采样时间
数据长度
1.0
10
6.5
真值
无测量噪声
0.2
30
0.9985
11.43
6.594
估计值
50
1.0000
11.31
6.6
0.02
10.13
6.51
有测量噪
(方差为0.01)
1
7.784
6.709
分析:
通过对比不同的采样周期和不同的采样时间在无测量噪声情况的辨识结果可知,在相同的采样周期下,适当的增加采样时间,可以提高辨识精度,尤其是对增益的提高有很大影响;
而在相同的采样时间下,适当的减小采样时间,对于系统参数的辨识精度有很大的提高。
因此,可以发现合理采样时间和数据长度,可以提高辨识的精度,令辨识后的传递函数系数与原传递函数系数更接近,差距小,从而得到满意的辨识结果。
通过对比无测量噪声和有测量噪声两种情况下的辨识结果,我们可以发现在白噪声的情况下,曲线拟合较无噪声情况下要差,说明白噪声对于面积法辨识系统存在较大的干扰,会对辨识结果产生一定的影响。
1.2.2有零点模型系统
17.5𝑆
2+7.5𝑆
+1假设系统的传递函数为G(s)=4
𝑆
3+5𝑆
2+8𝑆
+1,为有零点的模型,其中n=3,m=2,用面积法需要求解𝐴
1~𝐴
5,利用Matlab编程,分别在没有噪声和有噪声两种情况下进行辨识,比较辨识结果。
(1)没有噪声时,程序如下:
%==================获得原传递函数方程=======================%num=[17.57.51];
den=[4581];
%=====================产生阶跃采样序列======================%T=0.02;
100;
%数据长度y=step(num,den,t);
k=y(L)%系统增益
%======================面积法求解参数======================%sum1=0;
L-1;
sum1=sum1+(1-(y(i)+y(i+1))/2)*T;
A(i)=sum1;
endA1=sum1sum2=0;
sum2=sum2+(A(i)-A1*(y(i)+y(i+1))/2)*T;
B(i)=sum2;
endA2=sum2sum3=0;
sum3=sum3+(B(i)-A2*(y(i)+y(i+1))/2)*T;
C(i)=sum3;
endA3=sum3sum4=0;
sum4=sum4+(C(i)-A3*(y(i)+y(i+1))/2)*T;
D(i)=sum4;
endA4=sum4sum5=0;
sum5=sum5+(D(i)-A4*(y(i)+y(i+1))/2)*T;
end
A5=sum5
%==============根据所得A(i),利用公式求取a、b=================%
M=(-1)*(inv([A3,A2;
A4,A3]))*[A4;
A5];
b1=M(1,1);
b2=M(2,1);
N=[100;
A110;
A2A11]*[b1;
b2;
0]+[A1;
A2;
A3];
a1=N(1,1);
a2=N(2,1);
a3=N(3,1);
%================根据所求a、b,得到辨识后传递函数==============%num1=[b2b11];
den1=[a3a2a11];
%=============原传递函数和辨识函数的阶跃响应对比图=============%
step(G1,'
当采样时间取0.02秒,数据长度为100时,辨识结果如下:
图1.10
原系统与辨识后的系统阶跃响应对比图:
当采样时间为0.02,数据长度为400时,系统辨识结果如下:
图1-12
图1-13
当采样时间为0.2秒,数据长度为400时,系统辨识结果如下:
图1-14原系统与辨识后的系统阶跃响应对比图:
综上所述,结果如表
表1-2
a3
a2
a1
b2
b1
4
5
8
17.5
7.5
100
3.559
5.224
8.051
17.71
7.515
400
4.016
4.914
7.979
17.42
7.479
0.2
4.129
4.168
7.792
16.61
7.278
通过对比不同的采样周期和不同的采样时间在无测量噪声情况的辨识结果可知,对于存在有零点的系统来说,通过面积法辨识系统必须合理的选择分子分母的阶次,否则不能得出正确的辨识结果。
从表格中也可以发现,在相同的采样时间下,由于系统收敛过程较长从,增加数据长度对系统的辨识精度有了很大的提高;
同样在相同的数据长度下,合理的减小采样时间,也可以提高系统的辨识精度。
2.脉冲响应法
脉冲响应法(impulseresponsemethod)是指线性系统在零初始条件输入脉冲信号,信号后输出的瞬态响应,即输出响应叫脉冲响应。
建立系统的非参数模来观测系统脉冲响应,以求得系统数学模型的待定参数,进而实现系统辨识的目的。
2.1Hankel矩阵法
一个n阶过程的脉冲传递函数为
G(z)=−1
1−−1+bzaz22−2−2++LL++bznazn−n−n(2-1)bz1
1+az1+
将传递函数转化成为状态方程后,进一步推导,可知
G(z)=g
(1)z−1−1+g
(2)z−2+g(3)z−3+L
根据上述公式,可得
bz1−1+bz2−2+L+bzn−n
=g
(1)z−1+g
(2)z−2+g(3)z−3+L1+az1−1+az2−2+L+azn−n
=g
(1)z−1+L+g
(2)+ag1
(1)z−2+L+g(n)+n−1g(i)ani−z−n(2-2)
i=1
n2n−1
+g(n++1)g(i)an+−iz−+(n1)+L+g(2n)+g(i)a2ni−z−2n
i=1i=1脉冲相应法参数计算公式如下:
g
(1)g
(2)Lgn()angn(+1)
g
(2)g(3)Lgn(+1)an−1=−gn(+2)(2-3)
MMOMMM
gn()gn(+!
)Lgn(2−1)a1gn
(2)
b110L00g
(1)
b2=a11L00g
(2)
MMMOMMM
bnan−1an−2La11gn()(2-4)
Hankel矩阵的定义:
𝑔
(𝑘
)𝑔
+1)𝑔
+𝑙
−1)
⋯
+2)𝑔
H(l,k)=()(2-5)
⋮⋱⋮
−1)𝑔
)⋯𝑔
+2𝑙
−2)
公式(2-3)左边的矩阵是一种特定的Hankel矩阵。
记作H(n,l),并且它是可逆的。
因此,只需要将获得的脉冲响应值g(k),k=1,2,3,…,2n填入式(3-3)、
(3-4)便可以求出脉冲传递函数的估计值,即传递函数系数:
1,𝑎
2…,𝑎
𝑏
1,𝑏
2…,𝑏
。
2.2.实验过程
假设系统模型的传递函数为G(s)=𝑠
3+155𝑠
2𝑠
++61𝑠
+100,利用Matlab编程,分别在没有噪声和有噪声两种情况下进行辨识,比较辨识结果。
(1)无噪声的情况下,程序如下:
clear
%==================获得原传递函数方程=======================%
num=[51];
den=[11550100];
%=====================产生脉冲采样序列======================%T=0.02;
%采样时间gk=impulse(num,den,t);
%原传递函数的脉冲响应
%=====================Hamkel矩阵求解参数=====================%H=[gk
(1)gk
(2)gk(3);
gk
(2)gk(3)gk(4);
gk(3)gk(4)gk(5)];
%构造Hankel
矩阵求a、b系数矩阵
H1=H^-1;
%H求逆G1=-[gk(4);
gk(5);
gk(6)];
a=H1*G1;
%这里a=[a3;
a2;
a1]a=flipud(a).'
;
%a中的值倒转,a=[a1,a2,a3]
A=[100;
a
(1)10;
a
(2)a
(1)1];
%构造下三角阵G2=[gk
(1);
gk
(2);
gk(3)];
b=A*G2;
a=[1a];
原传递函数及其Z函数'
dsys1=c2d(G1,T)
辨识得到的Z函数和传递函数'
dsys2=filt(b.'
a,0.02)%求得辨识得到的Z函数
G2=T*d2c(dsys2,'
tustin'
0.02)%双线性变换
%=============原传递函数和辨识函数的脉冲响应对比图=============%
impulse(G1,'
)%绘出原系统与辨识得到的系统的脉冲
响应曲线
原模型与辨识后所得到系统脉冲响应对比'
原始响应曲线'
辨识曲线'
)程序运行结果如下:
原系统的传递函数和脉冲阶跃函数如下:
图2.1
辨识之后系统的传递函数和脉冲函数如下:
图2.2
在没有噪声的情况下,原系统与辨识系统脉冲响应对比图:
图2.3
(2)在噪声期望为0,方差为0.01的情况下,对系统进行参数辨识,在Mtlab中编程实现,程序如下:
clear
%==================获得原传递函数方程=======================%num=[51];
L=length(t);
%产生均值为0,方差为0.01的白噪声
w=randn(1,L);
w=w/std(w);
%=====================阶跃采样序列中加入白噪声==================%gk=gk+w;
%=====================Hamkel矩阵求解参数=====================%
H=[gk
(1)gk
(2)gk(3);
%构造Hankel
a=H1*G1;
%a中的值倒转,a=[a1;
a2;
a3]
G1=tf(num,den)dsys1=c2d(G1,T)
%=============原传递函数和辨识函数的脉冲响应对比图=============%impulse(G1,'
)%绘出原系统与辨识得到的系统的脉冲响应曲线title('
图2.4
原系统与辨识之后的系统脉冲响应图如下:
2.3结果分析
通过对比在没有噪声情况下原传递函数和辨识后系统的脉冲响应,我们发现两者的脉冲响应较为接近,但传递函数有一定差别。
而在加入白噪声的情况下,辨识后的系统脉冲响应与原传递函数的脉冲响应有较大的差别,存在很大的不确