Lezione 12

In questa lezione affronteremo l'uso delle tecniche Monte Carlo per la simulazione di una misura che avete scolto nel laboratorio di fisica (spettrometro a prisma). Come vedremo la simulazione ci aiuterà a capire meglio il comportamento dell'apparato sperimentale permettendoci per esempio di stimare l'incertezza attesa sulle grandezze da misurare.

ESERCIZIO 12.0 - Simulazione dell'esperienza dello spettrometro a prisma (da consegnare) :

La dipendenza dell'indice di rifrazione dalla lunghezza d'onda della luce che incide su un materiale viene descritta dalla legge di Cauchy: $n(\lambda) = \sqrt{ A + \frac{B}{\lambda^{2}} }$ dove $A$ e $B$ sono due costanti che dipendono dalle caratteristiche del materiale che attraversa. L'esperienza dello spettrometro a prisma si propone di misurare l'indice di rifrazione del materiale di un prisma per le diverse lunghezze d'onda di una lampada al mercurio per determinare i parametri $A$ e $B$ che caratterizzano tale materiale. L'apparato sperimentale consiste in un goniometro sul quale viene posizionato il prisma. Una lampada a vapori di mercurio viene posizionata da un lato del canocchiale con due collimatori per produrre un fascio luminoso che incide sul prisma. Il fascio di luce riflesso o rifratto viene osservato tramite un altro canocchiale. Gli angoli corrispondeni all'orientamento del supporto sul goniometro ed alla posizione dei canocchiali sono leggibili su di un nonio.

Image title

La misura sperimentale consiste nella determinazione:

dell'angolo corrispondente al fascio non deflesso in assenza del prisma: $\theta_{0}$ .
dell'angolo corrispondente alla deviazione minima della riga del giallo: $\theta_{m}(\lambda_{1})$ .
dell'angolo corrispondente alla deviazione minima della riga del viola: $\theta_{m}(\lambda_{2})$ ..

L'analisi dati consiste nella seguente procedura:

determinazione degli angoli di deviazione minima: $\delta_{m}(\lambda_{1}) = \theta_{m}(\lambda_{1}) - \theta_{0}$ $\delta_{m}(\lambda_{2}) = \theta_{m}(\lambda_{2}) - \theta_{0}$
calcolo degli indici di rifrazione $n(\lambda )$ dalla relazione: $n(\lambda_{1}) = \frac{ sin \frac{ \delta_{m}(\lambda_{1}) + \alpha }{2} } { sin \frac{\alpha}{2} }$ $n(\lambda_{2}) = \frac{ sin \frac{ \delta_{m}(\lambda_{2}) + \alpha }{2} } { sin \frac{\alpha}{2} }$
calcolo dei parametri $A$ e $B$ dalle formule: $A = \frac{ \lambda_{2}^{2} n^{2}(\lambda_{2}) - \lambda_{1}^{2} n^{2}(\lambda_{1}) }{\lambda_{2}^{2} - \lambda_{1}^{2} }$ $B = \frac{ n^{2}(\lambda_{2}) - n^{2}(\lambda_{1}) }{ \frac{1}{\lambda_{2}^{2}} - \frac{1}{\lambda_{1}^{2}}}$

La simulazione della misura parte dall'assumere dei valori verosimili dei parametri della legge di Cauchy e determinare come questi si traducono in quantità osservabili. A partire da queste quantità osservabili simuliamo una misura aggiungendo una perturbazione gaussiana ai valori nominali e ne determiniamo l'impatto sulle grandezze da misurare ( $A$ e $B$ ).

Nell'esperimento l'unico tipo di grandezze misurate sono gli angoli, per cui possiamo assumere un'incertezza uguale per tutte le misure angolari e pari a $\sigma_{\theta}=0.3$ mrad. Nell'esperienza di laboratorio l'angolo di apertura del prisma vale $\alpha=60^{\circ}$ ed il materiale del prisma ha valori dei parametri di Cauchy $A$ =2.7, $B$ =60000 nm $^{2}$ . Consideriamo le due lunghezze d'onda estreme della lampada al mercurio, il giallo, $\lambda_{1}$ =579.1 nm, ed il viola, $\lambda_{2}$ =404.7 nm.

Per costruire la nostra simulazione procediamo per gradi:

Costruire una classe EsperimentoPrisma con le seguenti caratteristiche:
- come data membri deve avere sia i valori 'veri' che i valori misurati di tutte le quantità ed in più un generatore di numeri casuali RandomGen (vedi lezione 10);
- nel costruttore deve definire tutti i valori 'veri' delle quantità misurabili a partire dai parametri $A$ , $B$ ed $\alpha$ e dalle lunghezze d'onda; N.B.: il valore di $\theta_{0}$ è arbitrario ma, una volta definito, i $\theta_{m}$ sono fissati.
Aggiungere un metodo Esegui() che effettua la misura sperimentale e determina dei valori misurati di $\theta_{0}$ , $\theta_{m}(\lambda_{1})$ e $\theta_{m}(\lambda_{2})$ ; N.B.: il valore misurato di un angolo si ottiene estraendo un numero distribuito in gaussianamente intorno al suo valore \'vero\' e deviazione standard $\sigma_{\theta}$ .
Aggiungere alla classe un metodo Analizza() che a partire dalle pseudomisure svolga l'analisi dei dati sino alla determinazione di A e B (misurati)
i metodi necessari per accedere ai valori dei data membri, sia quelli 'veri' sia quelli misurati.

Scrivere un programma che esegua 10000 volte l'esperimento, faccia un istogramma dei valori misurati, e calcoli media e deviazione standard di tali valori. Per il calcolo di medie e deviazioni standard potete decidere di immagazzinare i dati in un contenitore std::vector e utilizzare le funzioni sviluppate nelle prime lezioni oppure accedervi direttamente dagli istogrammi di ROOT come mostrato qui (ricordandosi di aggiungere histo.StatOverflows( kTRUE ); in modo da forzare l'utilizzo di eventuali underflow e overflow per calcoli statistici).

In particolare è utile studiare le seguenti distribuzioni :

distribuzione delle differenze tra i valori misurati e quelli attesi per $\theta_{0}$ , $\theta_{1}$ e $\theta_{2}$ . Calcolare media e deviazione standard delle distribuzioni.
distribuzione della differenza tra i valori misurati e quelli attesi di $\delta_{m}(\lambda_{1})$ e $\delta_{m}(\lambda_{2})$ , calcolare media e deviazione standard delle distribuzioni. Produrre un istogramma bidimensionale delle differenze per le due lunghezze d'onda e calcolare il coefficiente di correlazione.
distribuzione della differenza tra i valori misurati e quelli attesi di $n(\lambda_{1})$ e $n(\lambda_{2})$ , calcolare media e deviazione standard delle distribuzioni. Produrre un istogramma bidimensionale delle differenze per le due lunghezze d'onda, e calcolare il coefficiente di correlazione.
distribuzione delle differenza tra i valori misurati e quelli attesi di $A$ e $B$ , calcolare media e deviazione standard delle distribuzioni. Produrre un istogramma bidimensionale delle differenze e calcolare il coefficiente di correlazione.

In caso provate a confrontare l'errore stimato su $A$ e $B$ dalla simulazione con quello ottenuto dalla propagazione degli errori.

La classe EsperimentoPrisma :

La classe EsperimentoPrisma deve contenere al suo interno un generatore di numeri casuali per la simulazione del processo di misura, tutti i parametri che definiscono l'esperimento e, per le quantità misurate, sia il valore assunto sia il valore ottenuto dal processo di misura ed analisi dati dell'esperimento. Il suo header file, potrà pertanto essere del tipo:

#ifndef _Esperimento_Prisma_h_
#define _Esperimento_Prisma_h_

#include "RandomGen.h"

class EsperimentoPrisma {

 public :

  EsperimentoPrisma( unsigned int seed ) ;
  ~EsperimentoPrisma() {;} ;

  void Esegui() ;
  void Analizza() ;

  double getAmis() { return m_A_misurato ; } ;

  // .... aggiungere i metodi necessari

 private:

  // generatore di numeri casuali

  RandomGen m_rgen ;

  // parametri dell'apparato sperimentale

  double m_lambda1, m_lambda2, m_alpha, m_sigmat;

  // valori delle quantita' misurabili :
  // input    : valori assunti come ipotesi nella simulazione
  // misurato : valore dopo la simulazione di misura 

  double m_A_input, m_A_misurato;
  double m_B_input, m_B_misurato;
  double m_n1_input, m_n1_misurato;
  double m_n2_input, m_n2_misurato;
  double m_dm1_input, m_dm1_misurato;
  double m_dm2_input, m_dm2_misurato;
  double m_th0_input, m_th0_misurato;
  double m_th1_input, m_th1_misurato;
  double m_th2_input, m_th2_misurato;

};

#endif

N.B.: in questo header mancano i metodi tipo Get per accedere ai data membri.

La configurazione dell'esperimento, con il calcolo di tutti i valori 'veri' per le quantità misurabili, può venire fatta nel costruttore della classe, che in questo caso risulta più complicato del solito:

#include "EsperimentoPrisma.h"

EsperimentoPrisma::EsperimentoPrisma( unsigned int seed ) :
  m_rgen( seed ),
  m_lambda1(579.1E-9),
  m_lambda2(404.7E-9),
  m_alpha(60.*M_PI/180.),
  m_sigmat(0.3E-3),
  m_A_input(2.7),
  m_B_input(60000E-18)
{

  // calcolo degli indici di rifrazione attesi

  m_n1_input = sqrt( m_A_input + m_B_input / (m_lambda1*m_lambda1) ) ;
  m_n2_input = sqrt( m_A_input + m_B_input / (m_lambda2*m_lambda2) ) ;

  // theta0 e' arbitrario, scelgo M_PI/2.

  m_th0_input = M_PI/2. ;

  // determino theta1 e theta2 

  m_dm1_input = 2.*asin( m_n1_input * sin (0.5 * m_alpha) ) - m_alpha ;
  m_th1_input = m_th0_input + m_dm1_input ;
  m_dm2_input = 2.*asin( m_n2_input * sin (0.5 * m_alpha) ) - m_alpha ;
  m_th2_input = m_th0_input + m_dm2_input ;

}

Notate l'uso della lista di inizializzazione nel costruttore che abbiamo già visto nella scorsa lezione: questa permette di inizializzare direttamente i valori dei data membri anzichè procedere all'assegnazione dei valori nel costruttore dopo che i data membri siano stati costruiti. Nel nostro caso particolare diventa utile per inizializzare l'oggetto m_rgen invocando il costruttore opportuno ( altrimenti avrebbe usato il costruttore senza argomenti che nel nostro caso non esiste ). Potete trovare una buona spiegazione dell'utilizzo delle liste di inizializzazione qui.

Istrogrammi bidimensionali :

Per verificare le correlazioni tra due variabili si possono utilizzare istogrammi bidimensionali, in cui i bins sono definiti da range di valori sia di una variabile che dell'altra Per costruire questi istogrammi, utilizzeremo la classe TH2F di ROOT. Questa ha un costruttore che premette di dividere in bins sia la coordinata x che la y dell'istogramma:

TH2F::TH2F(char* nome, char* titolo, int n_bins_x, double xmin, double xmax, int n_bins_y, double ymin, double ymax);

La classe TH2F ha gli stessi metodi Fill() e Draw() della classe TH1F, l'unica differenza è che ora Fill ha bisogno della coppia completa di valori per decidere il canale in cui viene collocato l'evento:

TH2F::Fill(double x, double y);

Calcolo del coefficiente di correlazione :

Il coefficiente di correlazione $\rho$ tra due variabili $x$ e $y$ , può venire espresso come: $(<xy> - <x> <y>)/\sigma_{x}\sigma_{y}$ Si noti che dalla relazione: $\sigma_{x}^{2} = <x^{2}> - <x>^{2}$ risulta che può essere più conveniente nel ciclo di esperimenti calcolare le sommatorie di $x$ , $y$ , $x^{2}$ , $y^{2}$ e xy, ed al termine del ciclo effettuare le medie e trovare incertezze e correlazioni. Il tutto non richiede strettamente di immagazinare tutti i dati dei 10000 esperimenti.

Risultati attesi :

Image title