|
马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
%最小二乘法模态参数识别(复模态—频响函数实虚部)
clear
clc
close all hidden
format long
global mn
%fni=input('最小二乘法模态参数识别-输入数据文件名','s');
%fid=fopen(fni,'r')
%mn=fscanf(fid,'%f',1); %模态阶数
%df=fscanf(fid,'%d',1); %频率间隔
%f0=fscanf(fid,'%d',mn); %模态频率初值数组
%d0=fscanf(fid,'%d',mn); %模态阻尼比初值数组
%fno=fscanf(fid,'%d',1); %读入输出数据文件名
%b=fscanf(fid,'%f',[2 inf]); %实测频响函数实部虚部数据
%status=fclose(fid);
mn=3;df=0.1953125;f0=[5 10 20];d0=[0.05 0.05 0.05];fno='out8_2a.mat';
%tt=0:0.01:20;
%y=(1-0.2*tt).*sin(2*pi*5*tt);
%yy=fft(y);y=yy(1:1000);
%H0=[real(y);imag(y)];
load y
b=[y(1:10);y(11:20)];
f=0:df:(length(b(1,:))-1)*df;
w=2*pi*f;
H=b(1,:)+b(2,:)*i;
w0=2*pi*f0;
for j=1:mn
l=4*(j-1);
x0(l+1:l+4)=[-w0(j)*d0(j),w0(j)*sqrt(1-d0(j)^2),1,1];
end
x=lsqcurvefit('fun82',x0,w,H);
for j=1:mn
l=4*(j-1);
c=x(l+1)+i*x(l+2);
d=x(l+3)+i*x(l+4);
F(j)=abs(c)/(2*pi);
D(j)=-real(c)/abs(c);
S(j)=d;
end
H1=fun82(x,w);
subplot(2,1,1);
plot(f,real(H),':',f,real(H1));
xlabel('频率(Hz)');
ylabel('实部');
legend('实测','拟合');
grid on;
subplot(2,1,2);
plot(f,imag(H),':',f,imag(H1),'r');
xlabel('频率(Hz)');
ylabel('虚部');
legend('实测','拟合');
grid on;
fid=fopen(fno,'w');
fprintf(fid,' 频率(Hz) 阻尼比(%%) 振型系数\n');
for k=1:mn
fprintf(fid,'%10.4f %10.4f %10.6f\n',F(k),D(k)*100,imag(S(k)));
end
status=fclose(fid);
其中b=[y(1:10);y(11:20)];一旦改变y响亮的长度,比如改成b=[y(1:20);y(21:40)];就会出现错误。请教一下是怎么回事?
|
|