#############################
######### Z - score #########
############ J.B ############
########### 2021  ###########
#############################
#############################


import matplotlib.pyplot as plt
import numpy as np
from statistics import stdev, mean
plt.rcParams['font.size'] = 12
plt.close('all')

## mesures

L = 30.0 # cm distance laser-reseau
l = 13.5 # cm distance laser-ordre interference
p = 1
longueur_onde = 650 #nm / vert : 532 nm
# distances
u_L = .1  #cm incertitude
u_l = .1  #cm incertitude

#pour affichage
nom_resultat = "$a$" # pas du réseau
unite = "nm" #


## Méthode de Monté-Carlo
N =  10000 # nombre de tirages
L_alea = L + np.random.uniform(-u_L,u_L,N)
l_alea = l + np.random.uniform(-u_l,u_l,N)

tan_theta = l_alea/L_alea
sin_theta = np.sin(np.arctan(tan_theta))
a_alea =  p * longueur_onde / sin_theta


## représentation graphique
plt.figure(figsize=(10,6))

_,bins,_ = plt.hist(a_alea , 50,
                            density = True,
                            color ='green',
                            alpha = 0.2)
# Répartition gaussienne
sigma, mu = stdev(a_alea ), mean(a_alea )
y = ((1 / (np.sqrt(2 * np.pi) * sigma)) *
     np.exp(-0.5 * (1 / sigma * (bins - mu))**2))
plt.plot(bins, y, '--', color ='black')
# Affichage moyenne
plt.plot([mu,mu],[0,max(y)],'--r')
texte = "$\mu$ = {:.3e} ".format(mu) +  unite
plt.text(mu+sigma/10,max(y)/20,texte,color='r')
# Affichage écart type
plt.plot([mu,mu+sigma],[max(y)*.6,max(y)*.6],'k')
texte = "$\sigma =  {:.0e} $".format(sigma) + unite
plt.text(mu+sigma/10,max(y)*.56,texte)
# noms des axes et titre
plt.xlabel(nom_resultat + ' ('+unite + ')')
plt.ylabel('Fréquence')
plt.title('Méthode de Monté-Carlo')
plt.show()