#-----------------------------------------------------------------------
# Calcul d'incertitudes, loi normale, méthode Monte Carlo
#-----------------------------------------------------------------------

# Bibliothèques utilisées
import numpy as np
import scipy.optimize
from scipy.optimize import curve_fit
import matplotlib
import matplotlib.pyplot as plt

# Paramètres modifiables
N_exp = 1000 # nombre de tirages aléatoires
n_pnt = 10 # nombre de points d'expérience

# Valeurs mesurées (attention, la taille des listes doit correspondre à n_pnt !)
valU = np.array([]) # a remplacer par vos valeurs
valI = np.array([]) # a remplacer par vos valeurs
sU = np.array([]) # a remplacer par vos valeurs
sI = np.array([]) # a remplacer par vos valeurs

# Définitions des fonctions utiles
def make_data(n) :
    tab1 = np.zeros(n)
    tab2 = np.zeros(n)

    for k in range(n) :
        tab1[k] = float(np.random.normal(valU[k],sU[k],1))
        tab2[k] = float(np.random.normal(valI[k],sI[k],1))

    return tab1,tab2

def func(x, a): # regression linéaire de la forme y = a*x
    return a*x

# Tirages aléatoires et regressions linéaires
valR = np.zeros(N_exp)

for k in range(N_exp) :
    U,I = make_data(n_pnt)
    popt, pcov = curve_fit(func, I, U)
    valR[k] = float(popt)
    plt.plot(I,func(I,*popt),'x-')

# Calcul des indicateurs statistiques
Rm = np.mean(valR)
sigma = np.std(valR,ddof=1)
sR = sigma/np.sqrt(N_exp)
print("Valeur moyenne : %2.4f"  %Rm)
print("Incertitude type : %2.4f"  %sR)

# Tracés pour vérification (à commenter si execution trop longue)
plt.plot(valI,valU,'o',color='black',markersize=10)

plt.show()