luyunjun023 发表于 2009-1-3 23:17

ssi程序有点问题,请高手检查!!!

%hankel矩阵
clear
clc
close all hidden
d=100   %采样频率
n1=20   %行(2*n1)
n2=40   %列
load Y.txt %响应数据数列
y=Y';
yy=zeros(2*n1,n2);    %空hankel矩阵
for i=1:n2
    yy(:,i)=y(i:i+2*n1-1) %hankel矩阵第i列
end
hankel_a=yy/((n2)^0.5)   %最后的hankel矩阵
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%求奇异值
yp=hankel_a(1:n1,:);
yf=hankel_a(n1+1:2*n1,:);
pref=yf*yp'*pinv(yp*yp')*yp;
=size(pref);
pi=pref(1:r-1,:);   
=svd(pi); %奇异值分解得d0
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%定阶
%假如阶数为i
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%参数识别
pr=p0(:,1:i);    %U1
dd=d0(1:i,1:i);%S1
dr=dd^0.5;       %S1^0.5
qr=q0(:,1:i);    %V1’

oi=pr*dr;      %U1*S1^0.5
xi=pinv(oi)*pi;

%compuet xi1
clear yf p0 d0 q0 pr dd dr qr

pi=pref(2:r,:);   
=svd(pi);
xi1=pinv(oi)*pi;

a=xi1*pinv(xi);%A

=eig(a);
z1=diag(z);
z2=log(z1);
fr=abs(z2)*d;
damp=real(z2)./fr*d;

luyunjun023 发表于 2009-1-3 23:17

算出来结果和实际值差别很大

luyunjun023 发表于 2009-1-6 09:14

怎么没有人来了????都回家过年了吗

luyunjun023 发表于 2009-1-8 15:29

高手,你在哪里???

huazi071783 发表于 2010-4-26 16:09

楼主不知道你这个信号Y是几个通道的数据,看你的语句Y就是单通道的信号了,单个通道数据怎么能识别出阵型呢?
页: [1]
查看完整版本: ssi程序有点问题,请高手检查!!!