[Python] fitting theory with experimental data

I want to do fitting on experimental data and my theory I want to try to find the parameters that fit the same as the theory and experiment do I have to make an algorithm? what algorithm do I use to get a good fitting and mse for each parameter which is very small

I avoid using modules to do my fitting as much as possible trying to understand how the code works.

I use this mse

N = len(y)
error = (np.sum(y**2) - 2 * np.sum(y * y_t_fit) + np.sum(y_t_fit**2)) / N


this my code

import numpy as np
import matplotlib.pyplot as plt


#data
time = np.array([
0.000,
0.034,
0.068,
0.101,
0.135,
0.169,
0.203,
0.236,
0.270,
0.304,
0.338,
0.372,
0.405,
0.439,
0.473,
0.507,
0.541,
0.574,
0.608,
0.642,
0.676,
0.709,
0.743,
0.777,
0.811,
0.845,
0.878,
0.912,
0.946,
0.980,
1.013,
1.047,
1.081,
1.115,
1.149,
1.182,
1.216,
1.250,
1.284,
1.317,
1.351,
1.385,
1.419,
1.453,
1.486,
1.520,
1.554,
1.588,
1.622,
1.655,
1.689,
1.723,
1.757,
1.790,
1.824,
1.858,
1.892,
1.926,
1.959,
1.993,
2.027,
2.061,
2.094,
2.128,
2.162,
2.196,
2.230,
2.263,
2.297,
2.331,
2.365,
2.398,
2.432,
2.466,
2.500,
2.534,
2.567,
2.601,
2.635,
2.669,
2.703,
2.736,
2.770,
2.804,
2.838,
2.871,
2.905,
2.939,
2.973,
3.007,
3.040,
3.074,
3.108,
3.142,
3.145
])

y = np.array([
0.827,
0.828,
0.827,
0.817,
0.798,
0.769,
0.731,
0.692,
0.651,
0.607,
0.561,
0.510,
0.469,
0.435,
0.397,
0.368,
0.344,
0.325,
0.316,
0.314,
0.324,
0.333,
0.346,
0.372,
0.399,
0.429,
0.457,
0.485,
0.519,
0.545,
0.566,
0.585,
0.602,
0.618,
0.628,
0.650,
0.636,
0.619,
0.608,
0.592,
0.576,
0.559,
0.538,
0.513,
0.493,
0.478,
0.460,
0.444,
0.429,
0.422,
0.415,
0.413,
0.413,
0.417,
0.430,
0.441,
0.450,
0.465,
0.478,
0.493,
0.506,
0.518,
0.529,
0.536,
0.543,
0.550,
0.548,
0.546,
0.544,
0.539,
0.530,
0.520,
0.510,
0.504,
0.496,
0.486,
0.476,
0.471,
0.469,
0.463,
0.457,
0.456,
0.458,
0.459,
0.464,
0.471,
0.479,
0.484,
0.489,
0.497,
0.505,
0.507,
0.511,
0.517,
0.516
])

y -= 0.5

# find the max
def max_p (y):
result = []
for i in range(1,len(y)-1):
if y[i-1] < y > y[i+1]:
result.append(i)
return result
max = max_p(y)

#T1
y0 = y[max[0]]
yT = y[max[1]]
T = time[max[1]] - time[max[0]]

#T2
y1 = y[max[1]]
yT1 = y[max[2]]
T = time[max[2]] - time[max[1]]

max = max_p(y)

peak_values = [y for i in max]



#Parameter
A = 0.32
b = 0.5
k = 4
m = 1

w = np.sqrt(k/m - (b/2*m)**2)
T = 2*np.pi/w

y_t = np.exp((-b/(2*m))*t) * (A*np.cos(w*t))


plt.plot(time,y,label='experiment',linestyle='-')
plt.scatter(time[max], peak_values, color='red', marker='o', s=100, label="Max")
plt.plot(t,y_t,"r--")
plt.legend()
plt.grid()
plt.show()

Continue reading...