共享ERA/DC算法和SSI算法程序

幽幽深谷

function [fr,damp,v,d0]=eradc(n1,n2,h,d);
%  fr:frequecy matrix
%  damp:damp ratio matrix
%  n2:the order of hankel matrix

% generate R-matrix from h
[r,c]=size(h);
r(1:c-n1)=0;
for q=1:c-n1;
    r(q)=0;
    for i=1:n1;
        r(q)=h(q+i)*h(i)+r(q);
    end;
end;

% generate Hankel-matrix from R-matrix
h0=hankel(r(1:n2),r(n2:2*n2-1));
h1=hankel(r(2:n2+1),r(n2+1:2*n2));

% ERA start
[p0 d0 q0]=svd(h0);

%judge value's validlity
i=4;
cr1=3;
cr2=3;
cr3=3;
while (cr1>1.5 | cr2>1.5 | cr3>1.5) & i<20
    cr1=d0(i,i)/d0(i+2,i+2);
    cr2=d0(i+2,i+2)/d0(i+4,i+4);
    cr3=d0(i+4,i+4)/d0(i+6,i+6);
    i=i+2;
end;
i=i-4;

pr=p0(:,1:i);
dd=d0(1:i,1:i);
qr=q0(:,1:i);
dr=dd^(-0.5);

a=dr*pr'*h1*qr*dr;
[v,z]=eig(a);

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



n1是自相关函数的平均次数
21是hankel矩阵的阶数
h是脉冲的时程响应
d是采样频率

幽幽深谷

function [fr,damp,v,d0]=ssi2(n1,n2,h,d);

% generate R-matrix from h
[r,c]=size(h);
r(1:c-n1)=0;
for q=1:c-n1;
    r(q)=0;
    for i=1:n1;
        r(q)=h(q+i)*h(i)+r(q);
    end;
end;

% generate Hankel-matrix from R-matrix
h0=hankel(r(1:n2),r(n2:2*n2-1));

%compute xi
yp=h0(1:n2/2,:);
yf=h0(n2/2+1:n2,:);
pref=yf*yp'*pinv(yp*yp')*yp;
[r,c]=size(pref);
pi=pref(1:r-1,:);
[p0 d0 q0]=svd(pi);

%judge qiyi value
i=4;
cr1=3;
cr2=3;
cr3=3;
while (cr1>1.5 | cr2>1.5 | cr3>1.5) & i<20
    cr1=d0(i,i)/d0(i+2,i+2);
    cr2=d0(i+2,i+2)/d0(i+4,i+4);
    cr3=d0(i+4,i+4)/d0(i+6,i+6);
    i=i+2;
end;
i=i-4;

pr=p0(:,1:i);
dd=d0(1:i,1:i);
dr=dd^0.5;
qr=q0(:,1:i);

oi=pr*dr;
xi=pinv(oi)*pi;

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

pi=pref(2:r,:);
[p0 d0 q0]=svd(pi);
xi1=pinv(oi)*pi;

a=xi1*pinv(xi);

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

幽幽深谷

这两个程序是我自己编写用来做论文的。
若有不当之处让各位大虾指正

Galaxy163

能否介绍一下是做什么用的，谢谢！

幽幽深谷

可以识别频率和阻尼，如果在识别出振型，还需要加入几句语句。
因为我只需要识别频率和阻尼，所以没写关于振型的程序。

Galaxy163

谢谢共享！

wxywan368

把识别振型语句加进去才算完整的SSI，仁兄不妨把那几个语句也加进去，我们学习下！^_^

幽幽深谷

没时间加了
做论文的程序，本不应该贴出来的
如果哪位兄弟想用它发论文，请与我联系一下
体凉

ap0108220

:@)
楼主可否发一份完整的程序给我啊，正在赶论文～急需

luyunjun023

楼主能否提供详细的n1,n2,h,d
想学习！！！

cammer534

楼主极大的促进了国内振动学术研究的进展速度，一点都不夸张！！顶起。SSI 程序 我找了那么久  竟然在这里有。汗啊······

jgbjgb2010

谢谢！刚刚接触，先试试，非常感谢

secondye

楼主，如果有输入数据h就好了，现在不知道输入数据格式是怎么样的

huazi071783

楼主你好，稳定图是怎么画的啊？画稳定图的程序有么？

secondye

楼上的，你也是在随机子空间吗？用到稳定图法
我最近也在搞稳定图这方面，方便的话加我：79857660，我们交流一下，
我是华南理工大学的