import numpy as np
import matplotlib.pyplot as plt
import scipy.optimize as opt

D0 = 0.0065

eps = np.array([65,65,65,65,65,65, 69, 71, 72]) * 1e-4
t = np.array([0,0,0,0,0,0, 1, 60, 24*60]) / 60

print(eps)

def prony(t, *vars):
    D = vars[:2]
    tau = vars[2:]

    print(D, tau)
    D_val = D0
    for i in range(len(D)):
        D_val += D[i] * (1 - np.exp(- t/tau[i]))

    return D_val*25

coeff, cov = opt.curve_fit(prony, t, eps, p0=[10 for i in range(4)], bounds=(0, np.inf))
print(coeff)


t_list = np.logspace(-1, 5, 100)/3600
plt.scatter(0.1, D0, color='blue')
plt.scatter(t, eps, color='blue')
plt.plot(t_list, prony(t_list, *coeff))
plt.xscale('log')
#plt.xlim(0.1, 24*60*60)
#plt.ylim(min(eps), max(eps))
plt.show()
