disp(' ')
disp(' srs.m ver 2.0 July 3, 2006')
disp(' by Tom Irvine Email: tomirvine@aol.com') disp(' ')
disp(' This program calculates the shock response spectrum') disp(' of an acceleration time history, which is pre-loaded into Matlab.') disp(' The time history must have two columns: time(sec) & acceleration') disp(' ') %
clear t; clear y; clear yy; clear n; clear fn; clear a1; clear a2 clear b1; clear b2; clear jnum; clear THM; clear resp; clear x_pos; clear x_neg; %
iunit=input(' Enter acceleration unit: 1= G 2= m/sec^2 '); %
disp(' ')
disp(' Select file input method '); disp(' 1=external ASCII file ');
disp(' 2=file preloaded into Matlab '); file_choice = input(''); %
if(file_choice==1)
[filename, pathname] = uigetfile('*.*'); filename = fullfile(pathname, filename); %
fid = fopen(filename,'r');
THM = fscanf(fid,'%g %g',[2 inf]); THM=THM'; else
THM = input(' Enter the matrix name: '); end %
t=double(THM(:,1));
y=double(THM(:,2)); %
tmx=max(t); tmi=min(t); n = length(y); %
out1 = sprintf('\\n %d samples \\n',n); disp(out1) %
dt=(tmx-tmi)/(n-1); sr=1./dt; %
out1 = sprintf(' SR = %g samples/sec dt = %g sec \\n',sr,dt); disp(out1) %
fn(1)=input(' Enter the starting frequency (Hz) '); if fn(1)>sr/30. fn(1)=sr/30.; end %
idamp=input(' Enter damping format: 1= damping ratio 2= Q '); %
disp(' ') if(idamp==1)
damp=input(' Enter damping ratio (typically 0.05) '); else
Q=input(' Enter the amplification factor (typically Q=10) '); damp=1./(2.*Q); end %
disp(' ')
disp(' Select algorithm: ')
disp(' 1=Kelly-Richman 2=Smallwood '); ialgorithm=input(' '); %
tmax=(tmx-tmi) + 1./fn(1); limit = round( tmax/dt ); n=limit;
yy=zeros(1,limit); for i=1:length(y) yy(i)=y(i); end %
disp(' ')
disp(' Calculating response..... ') %
% SRS engine %
for j=1:1000 %
omega=2.*pi*fn(j);
omegad=omega*sqrt(1.-(damp^2)); cosd=cos(omegad*dt); sind=sin(omegad*dt); domegadt=damp*omega*dt; %
if(ialgorithm==1)
a1(j)=2.*exp(-domegadt)*cosd; a2(j)=-exp(-2.*domegadt); b1(j)=2.*domegadt;
b2(j)=omega*dt*exp(-domegadt);
b2(j)=b2(j)*( (omega/omegad)*(1.-2.*(damp^2))*sind -2.*damp*cosd ); b3(j)=0; %
else E=exp(-damp*omega*dt); K=omegad*dt; C=E*cos(K); S=E*sin(K); Sp=S/K; % a1(j)=2*C; a2(j)=-E^2; b1(j)=1.-Sp; b2(j)=2.*(Sp-C); b3(j)=E^2-Sp; end
forward=[ b1(j), b2(j), b3(j) ]; back =[ 1, -a1(j), -a2(j) ]; %
resp=filter(forward,back,yy); %
x_pos(j)= max(resp); x_neg(j)= min(resp); %
jnum=j;
if fn(j) > sr/8.
break end
fn(j+1)=fn(1)*(2. ^ (j*(1./12.))); end %
% Output options %
disp(' ')
disp(' Select output option ');
choice=input(' 1=plot only 2=plot & output text file ' ); disp(' ') %
if choice == 2 %%
[writefname, writepname] = uiputfile('*','Save SRS data as'); writepfname = fullfile(writepname, writefname); writedata = [fn' x_pos' (abs(x_neg))' ]; fid = fopen(writepfname,'w'); fprintf(fid,' %g %g %g\\n',writedata'); fclose(fid); %%
% disp(' Enter output filename '); % SRS_filename = input(' ','s'); %
% fid = fopen(SRS_filename,'w'); % for j=1:jnum
% fprintf(fid,'%7.2f .3f .3f \\n',fn(j),x_pos(j),abs(x_neg(j))); % end
% fclose(fid); end %
% Plot SRS %
disp(' ')
disp(' Plotting output..... ') %
% Find limits for plot %
srs_max = max(x_pos);
if max( abs(x_neg) ) > srs_max srs_max = max( abs(x_neg )); end
srs_min = min(x_pos);
if min( abs(x_neg) ) < srs_min srs_min = min( abs(x_neg )); end %
figure(1);
plot(fn,x_pos,fn,abs(x_neg),'-.'); %
if iunit==1
ylabel('Peak Accel (G)'); else
ylabel('Peak Accel (m/sec^2)'); end
xlabel('Natural Frequency (Hz)'); Q=1./(2.*damp);
out5 = sprintf(' Acceleration Shock Response Spectrum Q=%g ',Q); title(out5); grid;
set(gca,'MinorGridLineStyle','none','GridLineStyle',':','XScale','log','YScale','log');
legend ('positive','negative',2); %
ymax= 10^(round(log10(srs_max)+0.8)); ymin= 10^(round(log10(srs_min)-0.6)); %
fmax=max(fn); fmin=fmax/10.; %
fmax= 10^(round(log10(fmax)+0.5)); %
if fn(1) >= 0.1 fmin=0.1; end
if fn(1) >= 1 fmin=1; end
if fn(1) >= 10 fmin=10; end
if fn(1) >= 100 fmin=100; end
axis([fmin,fmax,ymin,ymax]); %
disp(' ')
百度搜索“77cn”或“免费范文网”即可找到本站免费阅读全部范文。收藏本站方便下次阅读，免费范文网，提供经典小说教育文库冲击响应谱计算的matlab程序在线全文阅读。
相关推荐：