SDOF time history analysis

how I can write a code to do time history analysis for a specific ground motion

[Partha Pratim]

k=2

[Mathieu NOE]

hello

there is quite a lot of available ressources on this topic :

processing_groundmotion - File Exchange - MATLAB Central

PEER Ground Motion Editor - File Exchange - MATLAB Central

OpenSeismoMatlab - File Exchange - MATLAB Central

GitHub - bakerjw/GMMs: A collection of ground motion model functions

Eigenquake simulation - File Exchange - MATLAB Central

Answer 1

[Mathieu NOE]

if you need a code to perform some spectral analysis and integration (to get velocity and displacement)

you can try this :

filename = "bcj.xlsx"; 
data = readmatrix(filename);
acc = data(6:end,2); % select here which channel to process
dt = 0.01;
fs = 1/dt; % sampling rate 
t = (0:numel(acc)-1)*dt;

acc = acc/100; % convert from cm/s² to m/s²

figure(1)
plot(t,acc)
ylabel('Accel [m/s²]')
xlabel('Time [s]')
title(' Acceleration')

% fft of the Acceleration signal 
L = length(acc); % Length of signal
NFFT = 2^nextpow2(L); % Next power of 2 from length of y
accSpectrum = fft(acc, NFFT)/L;  
f = fs/2*linspace(0,1,NFFT/2+1);
accSpectrum = abs(accSpectrum(1:NFFT/2+1));
figure(2);
semilogx(f, accSpectrum);
title('FFT of acceleration')
xlabel('Frequency (Hz)')
ylabel('Amplitude')
hold on 

% add findpeaks to display dominant frequencies
NP = 20;
[PKS,LOCS] = findpeaks(accSpectrum,f,'SortStr','descend','NPeaks',NP);
plot(LOCS,PKS,'dr');
for ck = 1:NP
    text(LOCS(ck)*1.1,PKS(ck)*1.2,[num2str(LOCS(ck),3) 'Hz']);
end
hold off

●

% main code
fc = 0.25;  % Hz , this cut off frequency must be below first significant peak frequency in your signal
[vel,displ] = double_int2(acc,fc,fs); % double integration (cumtrapz) with high pass filtering

figure(3)
plot(t,vel)
ylabel('Vel. [m/s]')
xlabel('Time [s]')
title('Estimated Velocity')

figure(4)
plot(t,displ)
ylabel('Dis. [m]')
xlabel('Time [s]')
title('Estimated Displacement')

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function [vel,displ] = double_int2(x,fc,Fs)

    dt = 1/Fs;

    % DC blocking filter (discrete)
    fn = fc/(Fs/2);  % normalized cutt off frequency
    p = (sqrt(3) - 2*sin(pi*fn))/(sin(pi*fn) + sqrt(3)*cos(pi*fn));
    b = [1 -1];                                  % (numerator)
    a = [1 -p];                                  % (denominator)

    accel = filtfilt(b,a,x);                  % filter the test data with filtfilt

    % accel to velocity integration
    vel = dt*cumtrapz(accel);           % integration 
    vel = filtfilt(b,a,vel);     % detrend data with filtfilt

    % velocity to displ integration
    displ = dt*cumtrapz(vel);           % integration
    displ = filtfilt(b,a,displ);     % detrend data with filtfilt

end

[Mathieu NOE]

hello

problem solved ?