I have the following code to fit the curve and show me the time constant tau. However, when I try to fit the experimental data it doesn't work. Anybody can help?
from scipy.signal import lti
import pylab as plt
from scipy.optimize import curve_fit
import math as m
def model1(x, gain1, tau1):
y = lti(gain1, [tau1, 1]).step(T=x)[1]
return y
def get_time_constants(time_vector, power_vector):
time_constant = []
for i in range(len(power_vector)):
# fit to output and estimate parameters - gain and tau
par1 = curve_fit(model, time_vector[i], power_vector[i])
y_fit = model(time_vector[i], par1[0][0], par1[0][1])
time_constant.append(par1[0][1])
time_constant_mean = sum(time_constant) / len(time_constant)
plt.plot(time_vector[i], power_vector[i])
plt.plot(time_vector[i], y_fit, label='Time Constant: %.2f s' %par1[0][1])
plt.show()
power = [[0.0, 68.13, 108.445, 133.43, 158.56, 164.575, 168.26, 172.035, 172.94, 173.795, 173.96, 174.145, 174.165, 174.195, 174.215, 174.29, 174.305, 174.325]]
time = [[0.0, 1.0, 2.0, 3.0999999999999943, 5.099999999999994, 6.099999999999994, 7.099999999999994, 9.099999999999994, 10.199999999999989, 12.199999999999989, 13.199999999999989, 14.199999999999989, 16.19999999999999, 17.299999999999983, 19.299999999999983, 20.30000000000001, 21.30000000000001, 23.30000000000001]]
get_time_constants(time, power)
As so many times ... just a question of starting values
def get_time_constants(time_vector, power_vector):
time_constant = []
for t, p in zip( time_vector, power_vector ):
guess_amp = max( p )
guess_tau = guess_amp / ( ( p[1] - p[0]) / ( t[1] - t[0] ) )
par1, _ = curve_fit(model, t, p, p0=[ guess_amp, guess_tau ])
y_fit = model(t, *par1)
time_constant.append(par1[1])
time_constant_mean = sum(time_constant) / len(time_constant)
plt.plot( t, p, ls='', marker='+' )
plt.plot( t, y_fit, label='Time Constant: %.2f s' %par1[1])
plt.legend(loc=0)
plt.show()
works just fine.