子空间系统辨识matlab,MIMO系统子空间系统辨识算法求助 - 信息科学 - 小木虫 - 学术 科研 互动社区...

我写的一段程序，可以参考下

/step 5 利用最小二乘计算B、D///

/*参考[Peter Van Overschee,96]Subspace Identification for Linear Systems:

Theory Implementation Applications,pp.56*/

Matrix<Type> Mi,Li,MM(I*(p*I-n),m),LL(I*(p*I-n),p*I),IGamd(p*I,p+n),BD;

Mi=trT(U2)*R31*pinv(R11);

Li=trT(U2);

for (int i=1;i<=I;i++)

{

for(int j=1;j<=p*I-n;j++)

{

for(int k=1;k<=m;k++)

MM((i-1)*(p*I-n)+j,k)=Mi(j,(i-1)*m+k);

}

}

for(int i=1;i<=I*(p*I-n);i++)

for(int j=1;j<=p*I;j++)

LL(i,j)=double(0);

for (int i=1;i<=I;i++)

{

for(int j=i;j<=I;j++)

{

for(int k=1;k<=p*I-n;k++)

{

for(int t=1;t<=p;t++)

LL((i-1)*(p*I-n)+k,(j-i)*p+t)=Li(k,(j-1)*p+t);

}

}

}

//设置矩阵IGamd=[Ip,0;0,Gamd]

for (int i=1;i<=p;i++)

{

for(int j=1;j<=p+n;j++)

{

if (j==i)

IGamd(i,j)=double(1);

else

IGamd(i,j)=double(0);

}

}

for (int i=1;i<=p*(I-1);i++)

{

for(int j=1;j<=p;j++)

IGamd(p+i,j)=double(0);

for(int j=1;j<=n;j++)

IGamd(p+i,p+j)=Gamd(i,j);

}

//QRD<Type> qrBD;qrBD.dec(LL*IGamd);BD=qrBD.solve(MM);

BD=pinv(LL*IGamd)*MM; //最小二乘计算BD=[D;B]

for(int i=0;i<n+p;i++)

{

if(i<p)

D.setRow(BD.getRow(i),i); //提取D矩阵

else

B.setRow(BD.getRow(i),i-p); //提取B矩阵

}

//计算算法耗时，