8.3. Examples: Lecture 8#

8.3.1. Generation of several toy experiments#

  • toys_01.py

#!/usr/bin/python
'''
generazione di N_toy toy experiment contenenti ciascuno N_evt campioni,
ciascuno dei quali segue una distribuzione uniforme fra 0 ed 1
python3 toys_01.py 100 1000
'''

import sys
from myrand import generate_uniform

def main () :
    '''
    Funzione che implementa il programma principale
    '''

    if len(sys.argv) < 3 :
        print ('usage: ', sys.argv[0], 'N_toy N_evt')
        exit ()

    N_toy = int (argv[1])
    N_evt = int (argv[2])
 
    # loop over toys
    for i in range (N_toy):
        randlist = generate_uniform (N_evt)


# ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- 


if __name__ == "__main__":
    main ()

8.3.2. Generation of several toy experiments and visualization of the means#

  • toys_02.py

#!/usr/bin/python
'''
generazione di N_toy toy experiment contenenti ciascuno N_evt campioni,
ciascuno dei quali segue una distribuzione uniforme fra 0 ed 1
e disegno della distribuzione delle medie dei campioni al variare dei toy
python3 toys_02.py 100 1000
'''

import sys
import numpy as np
from math import floor
import matplotlib.pyplot as plt

from myrand import generate_uniform
from stats import stats

def main () :
    '''
    Funzione che implementa il programma principale
    '''

    if len(sys.argv) < 3 :
        print ('usage: ', sys.argv[0], 'N_toy N_evt')
        exit ()

    N_toy = int (sys.argv[1])
    N_evt = int (sys.argv[2])
 
    means = []

    # loop over toys
    for i in range (N_toy):
        randlist = generate_uniform (N_evt)
        toy_stats = stats (randlist)
        means.append (toy_stats.mean ())

    means_stats = stats (means)
    xMin = means_stats.mean () - 5. * means_stats.sigma ()
    xMax = means_stats.mean () + 5. * means_stats.sigma ()
    nBins = floor (len (means) / 10.) + 1     # number of bins of the histogram
    bin_edges = np.linspace (xMin, xMax, nBins + 1)  # edges o the histogram bins
    fig, ax = plt.subplots ()
    ax.set_title ('Histogram of the mean over ' + str (N_toy) + ' toys', size=14)
    ax.set_xlabel ('mean value')
    ax.set_ylabel ('toys in bin')
    ax.hist (means,      # list of numbers
             bins = bin_edges,
             color = 'orange',
             # normed = True,
            )
    plt.show ()


# ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- 


if __name__ == "__main__":
    main ()

8.3.3. Generation of several toy experiments and visualization of the error on the means#

  • toys_03.py

#!/usr/bin/python
'''
generazione di N_toy toy experiment contenenti ciascuno N_evt campioni,
ciascuno dei quali segue una distribuzione uniforme fra 0 ed 1
e disegno della distribuzione delle medie dei campioni al variare dei toy
python3 toys_03.py 100 1000
'''

import sys
import numpy as np
from math import floor
import matplotlib.pyplot as plt

from myrand import generate_uniform
from stats import stats

def main () :
    '''
    Funzione che implementa il programma principale
    '''

    if len(sys.argv) < 3 :
        print ('usage: ', sys.argv[0], 'N_toy N_evt')
        exit ()

    N_toy = int (sys.argv[1])
    N_evt = int (sys.argv[2])
 
    means = []
    sigma_means = []

    # loop over toys
    for i in range (N_toy):
        randlist = generate_uniform (N_evt)
        toy_stats = stats (randlist)
        means.append (toy_stats.mean ())
        sigma_means.append (toy_stats.sigma_mean ())

    # compare the distribution of the sigma on the mean
    # calculated for each toy to the sigma of the mean distribution

    sigma_means_stats = stats (sigma_means)
    means_stats = stats (means)

    print ('sigma of the means disitribution:             ', means_stats.sigma ())
    print ('mean of the sigma of the means disitribution: ', sigma_means_stats.mean ())

    # plot the distribution of the sigma on the mean
    # calculated for each toy
    xMin = sigma_means_stats.mean () - 5. * sigma_means_stats.sigma ()
    xMax = sigma_means_stats.mean () + 5. * sigma_means_stats.sigma ()
    nBins = floor (len (sigma_means) / 10.) + 1     # number of bins of the histogram
    bin_edges = np.linspace (xMin, xMax, nBins + 1)  # edges o the histogram bins
    fig, ax = plt.subplots ()
    ax.set_title ('Histogram of the sigma on the mean over ' + str (N_toy) + ' toys', size=14)
    ax.set_xlabel ('mean value')
    ax.set_ylabel ('toys in bin')
    ax.hist (sigma_means,      # list of numbers
             bins = bin_edges,
             color = 'orange',
             # normed = True,
            )

    # get the sigma of the means distribution
    ax.plot ([means_stats.sigma (), means_stats.sigma ()], plt.ylim ())
    plt.show ()


# ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- 


if __name__ == "__main__":
    main ()

8.3.4. Uniform mapping of the (x,y) plane#

  • toys_04.py

#!/usr/bin/python
'''
Mappatura uniforme del piano con numeri pseudo-casuali
python3 toys_04.py 2. 4. 0. 4. 1000
'''

import sys
import numpy as np
from math import floor
import matplotlib.pyplot as plt

from myrand import generate_range
from stats import stats

def main () :
    '''
    Funzione che implementa il programma principale
    '''

    if len(sys.argv) < 6 :
        print ('usage: ', sys.argv[0], 'xMin xMax yMin yMax N_evt')
        exit ()

    xMin = float (sys.argv[1])
    xMax = float (sys.argv[2])
    yMin = float (sys.argv[3])
    yMax = float (sys.argv[4])
    N_evt = int (sys.argv[5])
 
    x_coord = generate_range (xMin, xMax, N_evt)
    y_coord = generate_range (yMin, yMax, N_evt)

    fig, ax = plt.subplots ()
    ax.set_title ('uniform mapping', size=14)
    ax.set_xlabel ('x axis')
    ax.set_ylabel ('y axis')
    ax.scatter (x_coord, y_coord)
    plt.show ()


# ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- 


if __name__ == "__main__":
    main ()

8.3.5. Hit-or-miss integration method#

  • toys_05.py

#!/usr/bin/python
'''
Metodo di integrazione hit-or-miss
python3 toys_05.py 0. 6.28 2. 1000
'''

import sys
import numpy as np
from math import floor
import matplotlib.pyplot as plt

from myrand import generate_range
from stats import stats


def func (x) : 
    return 1. + np.sin (x) ; 


# ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- 


def main () :
    '''
    Funzione che implementa il programma principale
    '''

    if len(sys.argv) < 5 :
        print ('usage: ', sys.argv[0], 'xMin xMax yMax N_evt')
        exit ()

    xMin = float (sys.argv[1])
    xMax = float (sys.argv[2])
    yMax = float (sys.argv[3])
    N_evt = int (sys.argv[4])
 
    x_coord = generate_range (xMin, xMax, N_evt)
    y_coord = generate_range (0., yMax, N_evt)

    x_under = []
    y_under = []
    for x, y in zip (x_coord, y_coord):
        if (func (x) > y) : 
            x_under.append (x)
            y_under.append (y)

    points_under = len (x_under)
    A_rett = (xMax - xMin) * yMax
    frac = float (points_under) / float (len (x_coord))
    integral = A_rett * frac
    integral_unc = A_rett**2 * frac * (1 - frac) / len (x_coord)

    print ('integral: ', integral)
    print ('integral uncertainty: ', integral_unc)

    fig, ax = plt.subplots ()
    ax.set_title ('try and catch', size=14)
    ax.set_xlabel ('x coord')
    ax.set_ylabel ('y coord')
    ax.scatter (x_coord, y_coord, color='lightskyblue')
    ax.scatter (x_under, y_under, color='steelblue')

    x = np.linspace(xMin, xMax, 1000)
    y = func (x)
    ax.plot (x, y, color = 'red')

    plt.show ()


# ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- ---- 


if __name__ == "__main__":
    main ()