Damped harmonic motion curve fit

45 views (last 30 days)
Stashu Kozlowski on 9 Oct 2019
Commented: Star Strider on 21 Oct 2020
Hey,
I have a data set in matlab, when plotted it looks like this: My goal is to determen a damped sinusoidal equation that would fit this data set, I honestly dont even know how to start. I have included my code, but it isn't much.Any help is much appreciated. Thank you!

Star Strider on 10 Oct 2019
Try this:
D = load('Stashu Kozlowski DHM.mat'); % File Attached
x = D.x;
y = D.y;
y = detrend(y); % Remove Linear Trend
yu = max(y);
yl = min(y);
yr = (yu-yl); % Range of ‘y’
yz = y-yu+(yr/2);
zci = @(v) find(v(:).*circshift(v(:), [-1 0]) <= 0); % Returns Approximate Zero-Crossing Indices Of Argument Vector
zt = x(zci(y));
per = 2*mean(diff(zt)); % Estimate period
ym = mean(y); % Estimate offset
fit = @(b,x) b(1) .* exp(b(2).*x) .* (sin(2*pi*x./b(3) + 2*pi/b(4))) + b(5); % Objective Function to fit
fcn = @(b) norm(fit(b,x) - y); % Least-Squares cost function
[s,nmrs] = fminsearch(fcn, [yr; -10; per; -1; ym]) % Minimise Least-Squares
xp = linspace(min(x),max(x), 500);
figure
plot(x,y,'b', 'LineWidth',1.5)
hold on
plot(xp,fit(s,xp), '--r')
hold off
grid
xlabel('Time')
ylabel('Amplitude')
legend('Original Data', 'Fitted Curve')
text(0.3*max(xlim),0.7*min(ylim), sprintf('$y = %.3f\\cdot e^{%.0f\\cdot x}\\cdot sin(2\\pi\\cdot x\\cdot %.0f%.3f)$', [s(1:2); 1./s(3:4)]), 'Interpreter','latex')
The estimated parameters are:
s =
-1.398211481931498e+00
-6.142349926864338e+02
2.591368008158479e-04
-5.442228857001487e+00
-3.075267405594925e-15
and the fit is nearly perfect: Star Strider on 21 Oct 2020
My pleasure!
It shouldn’t be difficult to generate a time-domain model of your circuit. Use that to define ‘fit’ in my code, instead of the funciton I posted. It may be necessary to change the number or order of initial parameter estimates (the ‘[yr; -10; per; -1; ym]’ in the fminsearch call) to matrch those in your model.
An updated (and slightly more robust) version of ‘zci’ is:
zci = @(v) find(diff(sign(v)));
The rest of the code is unchanged, except as necessary to work with your expression of ‘fit’.

Alex Sha on 10 Oct 2019
The fellow results are little better.
Parameter Best Estimate
---------- -------------
b1 -1.36099782974822
b2 -599.110824641553
b3 0.000259106153388041
b4 1.22915310606227
b5 0.0138196119722517