Esercizio statistica - INFN Sezione di Roma

Prova in itinere di statistica da svolgere a casa
Consegna 25 Giugno 2014
Simulazione Monte Carlo dell’interazione di una particella con un gas
1. Simulare lo spettro di energia di un fascio di particelle: L’energia di un
certo fascio di particelle è distribuita secondo un’esponenziale:
pdf(E) = 1/E0*exp(-E/E0)
=>
dN ~ N0/E0 exp(-E/E0) dE
(dove E0 corrisponde all’energia caratteristica dello spettro di particelle e N0 al
numero di particelle nel fascio)
Partendo dalla funzione di generazione di numeri casuali equiprobabili in un
dato intervallo, ricavare la distribuzione esponenziale (in R: runif(min=0,
max=1, n=1)). Suggerimento: la funzione cumulativa F = integral(pdf(x’)dx’, 0,
x) è distribuita uniformemente
• Fare l’istrogramma di pdf(E) avendo generato 1000 particelle
• Fare un fit con una distribuzione esponenziale
2. Simulare il numero di urti di una particella che attraversa un gas:
Il numero di collisioni cui una particella incorre attraversando un gas è ben
rappresentato da una distribuzione poissoniana. Simulare questo processo
ricorrendo alla definizione della distribuzione di poisson:
“Si immagini di contare il numero di volte in cui si verifica un certo evento (di qualsiasi tipo
purché ben definito) in un intervallo di tempo finito Dt; e si immagini di suddividere tale
intervallo in intervallini di tempo dt “sufficientemente piccoli”. Il processo è poissoniano se
posso trovare una dimensione di intervallino dt per cui valgono le seguenti proprietà:
(a) la probabilità di avere un unico conteggio in un tempo dt è proporzionale a dt
(b) la probabilità di avere un conteggio in un tempo dt è molto minore della probabilità di non
conteggi nello stesso dt, e la probabilità di avere più di un conteggio è trascurabile;
(c) il numero di conteggi che osservo in dt è indipendente dal numero di conteggi che
osservo in un altro intervallo da questo disgiunto.
Cfr. Bini, Lezioni di Statistica per la Fisica Sperimentale, pag. 105 ”
Nel caso in questione invece di avere un conteggio in un tempo dt si ha una
collisione in uno spazio dx.
Assumere di avere un tubo di lunghezza pari a x = 1 cm riempito con del gas,
una probabilità che una particella abbia una collisione con le molecole del gas
in dx = 0.01 mm pari a p = c*dx, con c = 2 mm-1 e che la probabilità per una
collisione di avvenire in un dato intervallino dx sia indipendente da quella negli
altri intervallini dx. (Suggerimento: utilizzare il generatore di numeri random
“runif(min=0, max=1, n=1)” per valutare la probabilità di successo dell’evento
“collisione”. )
3. Scrivere una funzione che simuli il passaggio di N elettroni attraverso
il gas contenuto in un tubo:
Una particella che attraversa un gas è soggetta a collisioni con gli elettroni
delle molecule del gas. In queste collisioni gli elettroni vengono espulsi dagli
orbitali atomici, e gli atomi vengono ionizzati. Le scale tipiche delle energie in
gioco sono MeV per l’energia della particella incidente ed eV per l’energia
perduta in ogni collisione e tipicamente in 1 cm di materiale gassoso a
pressione atmosferica si hanno qualche decina di collisioni. Supponiamo di
voler studiare il passaggio di un elettrone nel tubo catodico di un televisore.
3.1) assumere che la perdita media di energia in ciascuna collisione sia
<Epersa >= 15 eV, che il tubo riempito di gas abbia una lunghezza totale di 5
cm, che ci sia un’interazione ogni 0.05 cm (per il momento assumere che
questo sia un numero fisso), e che l’energia della particella incidente sia 10
MeV. Calcolare l’energia persa in funzione della distanza percorsa all’interno
del tubo assumendo che l’energia persa in ogni collisione sia distribuita
uniformemente tra 12 e 18 eV, ie pdf(Epersa) = unif(12, 18). (L’energia persa
dopo aver percorso una certa distanza d sarà uguale alla somma dell’energia
persa in ciascuna collisione avvenute in d). Si chiede di studiare la la forma
della distribuzione dell’energia totale persa in funzione dello percorso d della
particella nel tubo d = 0.05, 0.1, 0.5 cm. Commentare il risultato. Come
cambia il risultato se l’energia persa fosse descritta da pdf(Epersa) =
gauss(mu=15, sigma=2). (Riflettere sul teorema del limite centrale)
3.2) Assumere che l’energia persa in ciascuna collisione sia costante e pari a
a dE= 15 eV, ricollegandosi al punto 2. In questo caso, assumere che il
numero di interazioni non sia più fisso ma segua la statistica poissoniana con
valore medio pari a Nint=1, 5, 10, 20, 50, 100 (ovvero distanze pari a d = 0.05,
0.1, 0.25 0.5, 1, 2.5, 5 cm). Si chiede di
a) studiare la distribuzione e la deviazione standard dell’energia totale persa
dalla particella attraversando il gas in funzione del numero medio di
interazioni. Fare l’istogramma dell’energia persa per alcuni valori del numero
medio di interazioni: 1, 5, 10, 50
b) fare un grafico della deviazione standard dell’energia persa in funzione del
numero di interazioni e confrontare con un andamento del tipo
sigma(Etotpersa) = Sqrt(A*Nint + B)
c) Commentare il risultato.
d) Assumere di far attraversare un foglietto di materiale al fascio di elettroni
prima di raggiungere il volume di gas. Questo foglietto assorbe o “cede”
energia secondo la legge:
Epersa= gauss(mu=dE, sigma=3 dE)
Rifare il grafico ed il fit del punto b) e commentare il risultato.
e) Facoltativo. Per Nint=10 fare la convoluzione dell’esponenziale dell’energia
iniziale derivata nel punto 1) con la perdita di energia calcolata nel punto 3.2)
e studiare la distribuzione dell’energia finale in uscita dal tubo.
Esempio di energia persa attraversando il gas (valori numerici arbitrari):
Esempio di energia persa attraversando il gas ed il foglietto (valori
numerici arbitrari):
Suggerimenti per lo svolgimento dell’esercizio con R
- per creare un vettore: int = c(1,3,5,7,10,15,20,30,50,100)
- per inizializzare un vettore: energiaFinale <- rep(0, particelle)
- per ridefinire un vettore che contiene 3 liste di 1000 valori come una matrice
1000 righe per 3 colonne:
MyEnergia <- c(energiaFinale, energiaFinale, energiaFinale)
dim(MyEnergia) <-c(1000,3)
la seconda componente e’ MyEnergia[ ,2]
- for(I in seq(0,100)) per fare un ciclo da 1 a 100 con step 1
- if (a == 1) b=32]
- Generatore di numeri casuali equiprobabili: runif(min=0, max=1, n=1)
- Generatore di numeri casuali secondo una Gaussiana: rnorm(mu=1, sd=1,
n=1)
- Generatore di numeri casuali secondo Poisson: rpois(lambda=1, n=1)
(n e’ sempre il numero di numeri casuali da generare)
esempio di fit per y ~ sqrt(x)
x= c(1,3,5,7,10,15,20,30,50,100)
y= sqrt(int)
plot(x, y)
tab <- data.frame(x=x, r=y)
fit <- nls(r ~ sqrt(a*x+b), start=c(a=1, b=1), data = tab, algorithm='port',
lower=c(a=0, b=0), upper=c(a=Inf,b=Inf))
summary(fit)
p <-summary(fit)$parameters[,'Estimate']
p
perr <-summary(fit)$parameters[,'Std. Error']
perr
plot(function(x) sqrt(p[1]*x+p[2]), col=2, add=T, xlim=range(tab$x))

Download Report