# Damped Oscillation Equation Fitting

18 visualizzazioni (ultimi 30 giorni)
Haardhik il 29 Ott 2023
Commentato: Haardhik il 9 Nov 2023
I am attempting to fit Data Tables in Excel to a Sin curve in order to determine its Omega for a physics experiment, where I need the equation to be in the form
x' = Ae^-bt sin(ωt' + φ)
Unfortunately I have no idea about Matlab but when I plot them as Column vectors, it is the same as I would recieve in excel and nowhere close to a sin curve. Any and all help is appreciated, as I have been through the previous questions that are similar and they pick up from a sin curve, unable to find instructions to that step, I am posting this here. Attached is the different trials I wish to do this for.
##### 2 CommentiMostra NessunoNascondi Nessuno
Star Strider il 6 Nov 2023
Your data may be undersampled and could show aliasing. If you can, repeat the experiment with a much faster sampling frequency (100 times the current sampling frequency would work), then do the analysis.
Haardhik il 7 Nov 2023
This was the fastest possible sampling frequency on the equipment, it's just a high school physics lab with a basic Vernier Caliper Motion Sensor.

Accedi per commentare.

### Risposte (2)

Sam Chak il 29 Ott 2023
You can probably estimate the value of omega (angular frequency) by counting the number of zero-crossing events. If you want to fit the damped sinusoidal equation to the data, you can try this:
Warning: Column headers from the file were modified to make them valid MATLAB identifiers before creating variable names for the table. The original column headers are saved in the VariableDescriptions property.
Set 'VariableNamingRule' to 'preserve' to use the original column headers as table variable names.
Tab = 112×6 table
Time_s_ Position_m_ Time_s__1 Position_m__1 Time_s__2 Position_m__2 _______ ___________ _________ _____________ _________ _____________ 0.04 -0.0193 0.06 -0.01769 0.06 -0.013135 0.1 0.022134 0.1 0.020726 0.1 0.016225 0.14 -0.022867 0.14 -0.01769 0.14 -0.014782 0.18 0.022683 0.2 0.01661 0.18 0.01403 0.24 -0.017379 0.24 -0.01769 0.24 -0.013135 0.28 0.02186 0.28 0.019902 0.28 0.015677 0.32 -0.021495 0.32 -0.015495 0.32 -0.013135 0.36 0.020488 0.38 0.017158 0.38 0.012658 0.42 -0.017654 0.42 -0.017416 0.42 -0.01341 0.46 0.021311 0.46 0.01853 0.46 0.014853 0.5 -0.020398 0.5 -0.013849 0.5 -0.011763 0.54 0.018293 0.56 0.017433 0.56 0.013207 0.6 -0.017928 0.6 -0.016593 0.6 -0.013135 0.64 0.021037 0.64 0.017158 0.64 0.01403 0.68 -0.017928 0.7 -0.012751 0.68 -0.010117 0.74 0.016372 0.74 0.017433 0.74 0.013481
%% Trial 1
t = Tab.Time_s_(1:end-2); % Time
y = Tab.Position_m_(1:end-2); % Position
[~, count] = zerocrossrate(y, Method="comparison"); % count zero-crossing
T = t(end)/(count/2); % Period of a sine wave
omega = 2*pi/T % Estimate the angular frequency
omega = 68.4867
% Curve-Fitting
fo = fitoptions('Method','NonlinearLeastSquares',...
'Lower',[ 0, 0, omega-5],...
'Upper',[0.03, 1, omega+5],...
'StartPoint',[max(y) 0.5, omega]);
ft = fittype('a*exp(-b*x)*cos(c*x)', 'options',fo);
[curve, gof] = fit(t, y, ft) % check if omega ≈ c
curve =
General model: curve(x) = a*exp(-b*x)*cos(c*x) Coefficients (with 95% confidence bounds): a = 0.02388 (0.02315, 0.0246) b = 0.2507 (0.2357, 0.2657) c = 69.05 (69.03, 69.07)
gof = struct with fields:
sse: 1.7646e-04 rsquare: 0.9907 dfe: 107 adjrsquare: 0.9905 rmse: 0.0013
% Making plots
plot(t, y, '.'), hold on
plot(curve), grid on
xlabel('Time'), ylabel('Position'), title('Trial 1')
legend('Data', 'Fitted Sine')
##### 7 CommentiMostra 5 commenti meno recentiNascondi 5 commenti meno recenti
Sam Chak il 6 Nov 2023
The reason I asked for the initial velocity is to verify the values of a, b, c, and d found by Alex. Additionally, your data is undersampled because there are fewer than 2 data points per cycle (period). Please follow @Star Strider's advice to obtain more data.
% Using Alex Sha's finding
a = 0.023540378407752; % amplitude
b = 0.247177601901692; % decay rate constant
c = 69.1369475428933; % omega
d = 1.32644349058437; % phi
% period
T = 2*pi/c
T = 0.0909
% number of points per second
numPts = 110/5 % 110 data points over 5 sec
numPts = 22
% How many point in a period?
T*numPts % less than two points
ans = 1.9994
% Proposed fitting model
f = @(t) a*exp(-b*t).*sin(c*t + d);
t = linspace(0, 5, 19811); % minimum 19810 points over 5 sec to have 360 points per period
plot(t, f(t)), xlabel('t'), title('f(t) = a*exp(-b*t).*sin(c*t + d)')
% Velocity by formula (df/dt)
df = a*c*exp(-b*t).*cos(c*t + d) - a*b*exp(-b*t).*sin(c*t + d);
plot(t, df), xlabel('t'), title('df/dt')
% initial velocity based on df formula
df0 = df(1)
df0 = 0.3881
% initial velocity computed directly from {a, b, c, d}
df0 = a*c*cos(d) - a*b*sin(d)
df0 = 0.3881
% Test if initial velocity is -0.939057777778
True_df0 = -0.939057777778;
ToF = logical(abs(df0 - True_df0) < 1e-4) % True (1) / False (0)
ToF = logical
0
Haardhik il 9 Nov 2023
I realized that this data is sifted only for the peaks, If I were to apply the same to the actual data that sampled every 0.02 seconds, would it be possible?
I also noticed that you used a Cos wave in ur Fittype and hence will attach the original data herewith.

Accedi per commentare.

Alex Sha il 29 Ott 2023
if taking fitting function as: x=a*exp(-b*t)*sin(w*t+phi)
for trial-1
Sum Squared Error (SSE): 9.1948847955911E-5
Root of Mean Square Error (RMSE): 0.000914274913678052
Correlation Coef. (R): 0.999583948620075
R-Square: 0.999168070338902
Parameter Best Estimate
--------- -------------
a 0.023540378407752
b 0.247177601901692
w 69.1369475428933
phi 1.32644349058437
for trial-2
Sum Squared Error (SSE): 0.000115663498814737
Root of Mean Square Error (RMSE): 0.0010162233075687
Correlation Coef. (R): 0.999720103291284
R-Square: 0.999440284924735
Parameter Best Estimate
--------- -------------
a 0.0200937306323941
b 0.236486196043985
w 69.1712756320905
phi 0.896068977752904
for trial-3
Sum Squared Error (SSE): 9.64366701017117E-5
Root of Mean Square Error (RMSE): 0.000936320992461801
Correlation Coef. (R): 0.999614594177523
R-Square: 0.999229336892693
Parameter Best Estimate
--------- -------------
a 0.0158395749665591
b 0.215537219553407
w 69.1871168416771
phi 7.30505431226722
Note: phi is not unique
##### 4 CommentiMostra 2 commenti meno recentiNascondi 2 commenti meno recenti
Haardhik il 6 Nov 2023
I assume you used the CFTool function. I too attempted the same but in the custom equation part, I always got it wrong, could you please tell me how you reached here. I attempted the fitting fucntion you gave above and still failed.
Haardhik il 9 Nov 2023
If you could provide me the base code you used for this, it would be very helpful

Accedi per commentare.

### Categorie

Scopri di più su Linear Predictive Coding in Help Center e File Exchange

R2023b

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by