 |
 |
 |
 |
 |
 |
||
 |
 |
 |
 |
 |
| 5 | 1/1 | 返回列表 |
| 查看: 2570 | 回復(fù): 9 | ||
| 當(dāng)前只顯示滿足指定條件的回帖，點擊這里查看本話題的所有回帖 | ||
書尋玉鐵桿木蟲 (正式寫手)
|
[求助]
MIMO系統(tǒng)子空間系統(tǒng)辨識算法求助 已有1人參與
|
|
| 我想問下，大家有沒有做過MIMO系統(tǒng)的子空間系統(tǒng)辨識啊，最好是雙輸入雙輸出系統(tǒng)的辨識，有沒有相關(guān)的資料共享下，最好是MATLAB的例程。我按書上的方法編制了MATLAB程序，算出來不對，其中B、D怎么用最小二乘法求呢，最好有詳細(xì)的資料提供下啊，不勝感激啊。 |
回復(fù)此樓» 猜你喜歡
384求調(diào)劑
已經(jīng)有6人回復(fù)
-
生物學(xué)學(xué)碩求調(diào)劑
已經(jīng)有4人回復(fù)
-
081700 調(diào)劑 267分
已經(jīng)有6人回復(fù)
-
一志愿河北工業(yè)大學(xué)0817化工278分求調(diào)劑
已經(jīng)有8人回復(fù)
-
269求調(diào)劑
已經(jīng)有4人回復(fù)
-
一志愿211 初試270分 求調(diào)劑
已經(jīng)有5人回復(fù)
-
265求調(diào)劑
已經(jīng)有10人回復(fù)
-
327求調(diào)劑
已經(jīng)有3人回復(fù)
-
材料專業(yè)求調(diào)劑
已經(jīng)有9人回復(fù)
-
291求調(diào)劑
已經(jīng)有7人回復(fù)
southisland
新蟲 (初入文壇)
- 應(yīng)助: 0 (幼兒園)
- 金幣: 88
- 帖子: 10
- 在線: 3.7小時
- 蟲號: 2987725
- 注冊: 2014-02-23
- 專業(yè): 控制理論與方法
9樓2014-03-21 09:21:33
southisland
新蟲 (初入文壇)
- 應(yīng)助: 0 (幼兒園)
- 金幣: 88
- 帖子: 10
- 在線: 3.7小時
- 蟲號: 2987725
- 注冊: 2014-02-23
- 專業(yè): 控制理論與方法
【答案】應(yīng)助回帖
|
我寫的一段程序，可以參考下 /////////////////////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); } } } //設(shè)置矩陣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矩陣 } //計算算法耗時 |
2樓2014-02-23 16:50:39
southisland
新蟲 (初入文壇)
- 應(yīng)助: 0 (幼兒園)
- 金幣: 88
- 帖子: 10
- 在線: 3.7小時
- 蟲號: 2987725
- 注冊: 2014-02-23
- 專業(yè): 控制理論與方法
【答案】應(yīng)助回帖
|
找到了當(dāng)時寫的matlab函數(shù) function [A,B,C,D]=mysubid(Y,U,I,obqr) %基本的線性開環(huán)子空間辨識算法 %Y：辨識輸出信號， U：辨識輸入信號， I：Hankel矩陣行數(shù) %obqr:如果使用QR分解的斜向投影，則置于非0的任意數(shù) %obqr==0,使用斜向投影 %計算時可選直接投影或者基于QR分解的斜向投影 %extract matrices A & C from matric Gam,then Solve least square: %M=L*X*[B;D] % 2012-4-11 rewrite 2012-5-11 %2012-7-11 增加可測噪聲模型的辨識 %%%START%%% if (nargin < 4);obqr = 0;end [m,N]=size(U);p=size(Y,1); J=N-2*I+1; % obqr = 1; %% STEP 1: 構(gòu)造輸入輸出數(shù)據(jù)hankel矩陣塊% HklU=blkhank(U,I,J ); Up=[];Yp=[];Uf=[];Yf=[]; for i=1:I Up=[Up;U(:,i:J+i-1)]; Uf=[Uf;U(:,I+i:J+I+i-1)]; Yp=[Yp;Y(:,i:J+i-1)]; Yf=[Yf;Y(:,I+i:J+I+i-1)]; end %% STEP 2：投影并加權(quán)（直接投影和QR分解斜向投影） %直接正交投影 Wp=[Yp;Up]; if obqr==0 PI_Uf=eye(J)-Uf'*pinv(Uf*Uf')*Uf; O=Yf*PI_Uf*Wp'*pinv(Wp*PI_Uf*Wp')*Wp; else %采用QR分解計算斜向投影 YuW=[Uf;Wp;Yf]; [yu_q,yu_r]=qr(YuW'); R= yu_r';% R = R(1:2*I*(m+p),1:2*I*(m+p)); % Truncate Q=yu_q'; O=R(I*p+2*I*m+1:2*I*(m+p),p*I+1:I*(2*p+m))*inv(R(m*I+1  2*m+p)*I,m*I+12*m+p)*I))...*R(m*I+1 2*m+p)*I,12*m+p)*I)*Q(12*m+p)*I, ;end %% STEP 3：SVD分解，選擇系統(tǒng)狀態(tài)階數(shù) [Uu,S,V] = svd(O); ss = diag(S)'; lnss=log(ss); h = bar(1:I*p,lnss,'b'); grid on;title('Singular Value Histogram');xlabel('The Number of Singular Value ');ylabel('Singular Value'); n = input(' System order ? ');%選擇系統(tǒng)階數(shù) %% STEP 4：由Gam利用其移不變性計算C和A U1=Uu(:,1:n);V1=V(1:n, ;r=sqrt(S(1:n,1:n)); Gam=U1*r; Xf=r*V1; C_h=Gam(1:p,1:n); A_h=pinv(Gam(1:p*(I-1),1:n))*Gam(p+1:p*I,1:n); %% STEP 5：計算B和D: 對于M=L*X*[B;D]采用LS計算B,D tao=Gam';%注意：轉(zhuǎn)置矩陣 PI_Gam=eye(size(tao,2))-tao'*inv(tao*tao')*tao;%Gam轉(zhuǎn)置的正交投影算子，PI_Gam*Gam=0 Ms=PI_Gam*Yf*pinv(Uf); for k=1:I M(I*p*(k-1)+1:I*p*k, =Ms(:,m*(k-1)+1:m*k);%構(gòu)造M的列塊矩陣L(I*p*(k-1)+1:I*p*k, =[PI_Gam(:,p*(k-1)+1:p*I),zeros(p*I,p*(k-1))];%構(gòu)造L的塊方陣end IG = [eye(p),zeros(p,n);zeros(p*(I-1),p),Gam(1 I-1)*p,];%構(gòu)造【I,0;0,gamma】矩陣% Solve least squares sol_bd = (L*IG)\M; D_h = sol_bd(1:p, ;B_h = sol_bd(p+1:p+n,; |
3樓2014-02-23 16:55:40
書尋玉
鐵桿木蟲 (正式寫手)
- 應(yīng)助: 2 (幼兒園)
- 金幣: 5543.2
- 散金: 3
- 紅花: 22
- 帖子: 478
- 在線: 241.1小時
- 蟲號: 2209924
- 注冊: 2012-12-27
- 性別: GG
- 專業(yè): 制造系統(tǒng)與自動化
4樓2014-02-23 18:48:04
| 最具人氣熱帖推薦 [查看全部] | 作者 | 回/看 | 最后發(fā)表 | |
|---|---|---|---|---|
|
[考研] 291求調(diào)劑 +7 | hhhhxn.. 2026-03-23 | 7/350 |
|
|---|---|---|---|---|
|
[考研] 材料專碩英一數(shù)二306 +8 | z1z2z3879 2026-03-18 | 8/400 |
|
|
|
[考研] 336求調(diào)劑 +4 | 收到VS 2026-03-20 | 4/200 |
|
|
|
[考研] 328求調(diào)劑，英語六級551，有科研經(jīng)歷 +7 | 生物工程調(diào)劑 2026-03-17 | 12/600 |
|
|
|
[考研] 316求調(diào)劑 +7 | 梁茜雯 2026-03-19 | 7/350 |
|
|
|
[考研] 招08考數(shù)學(xué) +6 | laoshidan 2026-03-20 | 14/700 |
|
|
|
[考研] 310求調(diào)劑 +4 | baibai1314 2026-03-16 | 4/200 |
|
|
|
[考研] 307求調(diào)劑 +11 | 冷笙123 2026-03-17 | 11/550 |
|
|
|
[考研] 一志愿華中農(nóng)業(yè)071010，總分320求調(diào)劑 +5 | 困困困困坤坤 2026-03-20 | 6/300 |
|
|
|
[考研] 298求調(diào)劑一志愿211 +3 | 上岸6666@ 2026-03-20 | 3/150 |
|
|
|
[考研] 269專碩求調(diào)劑 +6 | 金恩貝 2026-03-21 | 6/300 |
|
|
|
[考研] 求調(diào)劑 +5 | Zhangbod 2026-03-21 | 7/350 |
|
|
|
[考研] 考研調(diào)劑 +4 | 來好運來來來 2026-03-21 | 4/200 |
|
|
|
[考研] 材料工程專碩 348分求調(diào)劑 +3 | 冬辭. 2026-03-17 | 5/250 |
|
|
|
[考研] 297求調(diào)劑 +3 | 喜歡還是不甘心 2026-03-20 | 3/150 |
|
|
|
[考研] 一志愿南昌大學(xué)，327分，材料與化工085600 +9 | Ncdx123456 2026-03-19 | 9/450 |
|
|
|
[考研] 308求調(diào)劑 +3 | 阿姐阿姐家啊 2026-03-18 | 3/150 |
|
|
|
[考研] 一志愿西南交通 專碩 材料355 本科雙非 求調(diào)劑 +5 | 西南交通專材355 2026-03-19 | 5/250 |
|
|
|
[考研] 求調(diào)劑 +3 | eation27 2026-03-20 | 3/150 |
|
|
|
[考研] 【同濟軟件】軟件（085405）考研求調(diào)劑 +3 | 2026eternal 2026-03-18 | 3/150 |
|