divide by zero encountered in true_divide error without having zeros in my data

Question

this is my code and this is my data, and this is the output of the code. I've tried adding one the values on the x axes, thinking maybe values so little can be interpreted as zeros. I've no idea what true_divide could be, and I cannot explain this divide by zero error since there is not a single zero in my data, checked all of my 2500 data points. Hoping that some of you could provide some clarification. Thanks in advance.

import pandas as pd
import matplotlib.pyplot as plt
from iminuit import cost, Minuit
import numpy as np

frame = pd.read_excel('/Users/lorenzotecchia/Desktop/Analisi Laboratorio/Analisi     dati/Quinta Esperienza/500Hz/F0000CH2.xlsx', 'F0000CH2')
data =  pd.read_excel('/Users/lorenzotecchia/Desktop/Analisi Laboratorio/Analisi     dati/Quinta Esperienza/500Hz/F0000CH1.xlsx', 'F0000CH1')
# tempi_500Hz = pd.DataFrame(frame,columns=['x'])
# Vout_500Hz = pd.DataFrame(frame,columns=['y'])
tempi_500Hz = pd.DataFrame(frame,columns=['x1'])
Vout_500Hz = pd.DataFrame(frame,columns=['y1']) 
# Vin_500Hz = pd.DataFrame(data,columns=['y'])

def fit_esponenziale(x, α, β):
    return α * (1 - np.exp(-x / β))

plt.xlabel('ω(Hz)')
plt.ylabel('Attenuazioni')
plt.title('Fit Parabolico')
plt.scatter(tempi_500Hz, Vout_500Hz)
least_squares = cost.LeastSquares(tempi_500Hz, Vout_500Hz, np.sqrt(Vout_500Hz), fit_esponenziale)
m = Minuit(least_squares, α=0, β=0)
m.migrad()
m.hesse()
plt.errorbar(tempi_500Hz, Vout_500Hz, fmt="o", label="data")
plt.plot(tempi_500Hz, fit_esponenziale(tempi_500Hz, *m.values), label="fit")
fit_info = [
f"$\\chi^2$ / $n_\\mathrm{{dof}}$ = {m.fval:.1f} / {len(tempi_500Hz) - m.nfit}",]

for p, v, e in zip(m.parameters, m.values, m.errors):
    fit_info.append(f"{p} = ${v:.3f} \\pm {e:.3f}$")

plt.legend()
plt.show()

input

output and example of data

Answer 1

Here is a working Minuit vs curve_fit example. I scaled the function such that the decay in the exponential is in the order of 1 (generally a good idea for non linear fits ). Eventually, both methods give very similar results.

Note:I leave it open whether the error makes sense like this or not. The starting values equal to zero was definitively a bad idea.

import pandas as pd
import matplotlib.pyplot as plt
from iminuit import cost, Minuit
from scipy.optimize import curve_fit
import numpy as np

def fit_exp(x, a, b, c):
    return a * (1 - np.exp(- 1000 * b * x) ) + c

nn = 170
xl = np.linspace( 0, 0.001, nn )
yl = fit_exp( xl, 15, 5.3, -8.1 ) + np.random.normal( size=nn, scale=0.05 )

#######################
### Minuit
#######################
least_squares = cost.LeastSquares(xl, yl, np.sqrt( np.abs( yl ) ), fit_exp )
print(least_squares)
m = Minuit(least_squares, a=1, b=5, c=-7)
print( "grad: ")
print( m.migrad() ) ### needs to be called to get fit values
print( m.values )### gives slightly different output
print("Hesse:")
print( m.hesse() )

#######################
### curve_fit
#######################

opt, cov = curve_fit(
    fit_exp, xl, yl, sigma=np.sqrt( np.abs( yl ) ),
    absolute_sigma=True
)
print( " curve_fit: ")
print( opt )
print( " covariance ")
print( cov )

fig = plt.figure()
ax = fig.add_subplot(1, 1, 1 )
ax.plot(xl, yl, marker ='+', ls='' )

ax.plot(xl, fit_exp(xl, *m.values), ls="--")
ax.plot(xl, fit_exp(xl, *opt), ls=":")


plt.show()

Answer 2

There is actually a way to fit this completely linear. See e.g.here

import matplotlib.pyplot as plt
import numpy as np

from scipy.integrate import cumtrapz



def fit_exp(x, a, b, c):
    return a * (1 - np.exp( -b * x) ) + c

nn = 170
xl = np.linspace( 0, 0.001, nn )
yl = fit_exp( xl, 15, 5300, -8.1 ) + np.random.normal( size=nn, scale=0.05 )

"""
with y = a( 1- exp(-bx) ) + c
we have Y = int y = -1/b y + d x + h ....try it out or see below
so we get a linear equation for b (actually 1/b ) to optimize
this goes as:
"""

Yl = cumtrapz( yl, xl, initial=0 )
ST = [xl, yl, np.ones( nn ) ]
S = np.transpose( ST )

eta = np.dot( ST, Yl )
A = np.dot( ST, S )

sol = np.linalg.solve( A, eta )
bFit = -1/sol[1]
print( bFit )
"""
now we can do a normal linear fit
"""

ST = [ fit_exp(xl, 1, bFit, 0), np.ones( nn ) ]
S = np.transpose( ST )
A = np.dot( ST, S )
eta = np.dot( ST, yl )
aFit, cFit = np.linalg.solve( A, eta )

print( aFit, cFit )
print(aFit + cFit, sol[0] ) ### consistency check

fig = plt.figure()
ax = fig.add_subplot(1, 1, 1 )
ax.plot(xl, yl, marker ='+', ls='' )
## at best a sufficient fit, at worst a good start for a non-linear fit
ax.plot(xl, fit_exp( xl, aFit, bFit, cFit ) ) 

plt.show()

"""
a ( 1 - exp(-b x)) + c = a + c - a exp(-b x) = d - a exp( -b x )
int y = d x + a/b exp( -b x ) + g
    = d x +a/b exp( -b x ) + a/b - a/b + c/b - c/b + g
    = d x - 1/b ( a - a exp( -b x ) + c ) + c/b + a/b + g
    = d x + k y + h
with k = -1/b and h = g + c/b + a/b.
d and h are fitted but not used, but as a+c = d we can check 
for consistency
"""