#J.COURTIN 01/2023         MONTE CARLO  pour BESSEL

from random import gauss




#formule de Bessel
def focale(D,e):
    return (D**2-e**2) / (4*D)





#Génération des mesuree aléatoire fictive
def mesureAleatoire(D,e,DD,De,N=100):
    alea_D=[]
    alea_e=[]
    maFocale=[] #initialisation liste vides

    for i in range(N):
        alea_D+=[gauss(D,DD)] #création des variables aléatoires
        alea_e+=[gauss(e,De)]
        maFocale+=[focale(alea_D[-1], alea_e[-1])] #focale simulée

    return (alea_D, alea_e, maFocale)




###### MAIN ##########

D = 80.0   #paramètres de distance
e = 50.7   #mesurés D et e

DD = 0.1   #incertitude mesurée
De = 0.7

N=1000     # nombre de mesures simulées


#Simulation :
mes_D, mes_e, ma_focale = mesureAleatoire(D,e,DD,De,N)

print("D = ", mean(mes_D), " +/- ", Sigma(mes_D) )
print("e = ", mean(mes_e), " +/- ", Sigma(mes_e) )


print("focale : f' = ", mean(ma_focale), " +/- ", Sigma(ma_focale) )