Soluzioni: Lezione 11

Contents

Soluzioni: Lezione 11¶

../_images/linea.png

11.1 generazione di coppie di numeri pseudo-casuali¶

esercizio_01.cpp

/*
c++ -o esercizio_01 esercizio_01.cpp funzioni.cc `root-config --glibs --cflags`
*/

#include <iostream>
#include <vector>
#include <cmath>

#include "TF1.h"
#include "TGraph.h"
#include "TCanvas.h"

#include "funzioni.h"


using namespace std ;

int main (int argc, char ** argv)
{

  int nPoints = 5 ;
  double min = 0.5 ;
  double max = 10.5 ;
  double epsilon = 0.15 ;

  TF1 f_true ("f_true", retta, 0., 2. * nPoints, 2) ;
  f_true.SetParameter (0, 0.4) ;
  f_true.SetParameter (1, 1.5) ;

  double x = 0. ;
  for (int i = 0 ; i < nPoints ; ++i)
    {
      x+= rand_range (0., 1.) ;
      double y = f_true.Eval (x) + rand_TCL (0, epsilon) ;
    }
  
  return 0 ;
}

funzioni.h

#ifndef funzioni_h
#define funzioni_h

// retta secondo il prototipo utilizzato da root
double retta (double * x, double * par) ;

// generazione di numeri casuali su un intervallo 
// con distribuzione uniforme
double rand_range (double min, double max) ;

// generazione di numeri casuali
// distribuiti secondo una Gaussiana
// con il metodo del teorema centrale del limite
float rand_TCL (double mean, double sigma, int N = 10) ;


#endif

funzioni.cc

#include "funzioni.h"
#include <cmath>


double retta (double * x, double * par)
  {
    return par[0] + par[1] * x[0] ;
  }


double rand_range (double min, double max)
  {
    return min + (max - min) * rand () / static_cast<double> (RAND_MAX) ;
  } 


float rand_TCL (double mean, double sigma, int N)
  {
    double y = 0. ; 
    double xMin = mean - sqrt (3 * N) * sigma ;
    double xMax = mean + sqrt (3 * N) * sigma ;
    for (int i = 0 ; i < N ; ++i)
      y += rand_range (xMin, xMax) ;
    y /= N ;
    return y ;
  }

../_images/linea.png

11.2 disegno della generazione dei punti¶

esercizio_02.cpp

/*
c++ -o esercizio_02 esercizio_02.cpp funzioni.cc `root-config --glibs --cflags`
*/

#include <iostream>
#include <vector>
#include <cmath>

#include "TF1.h"
#include "TGraphErrors.h"
#include "TCanvas.h"

#include "funzioni.h"


using namespace std ;

int main (int argc, char ** argv)
{

  int nPoints = 5 ;
  double min = 0.5 ;
  double max = 10.5 ;
  double epsilon = 0.15 ;

  TF1 f_true ("f_true", retta, 0., 2. * nPoints, 2) ;
  f_true.SetParameter (0, 0.4) ;
  f_true.SetParameter (1, 1.5) ;

  vector<double> v_x ;
  vector<double> v_y ;

  double x = 0. ;
  for (int i = 0 ; i < nPoints ; ++i)
    {
      x+= rand_range (0., 1.) ;
      v_x.push_back (x) ;
      double y = f_true.Eval (x) + rand_TCL (0, epsilon) ;
      v_y.push_back (y) ;
    }
  
  // errore sulle coordinate x dei punti
  vector<double> v_ex (v_x.size (), 0.) ;
  // errore sulle coordinate y dei punti
  vector<double> v_ey (v_y.size (), epsilon);

  TGraphErrors g_retta (v_x.size (), &v_x[0], &v_y[0], &v_ex[0], &v_ey[0]) ;
  g_retta.SetMarkerStyle (4) ;
  g_retta.SetMarkerColor (kRed) ;

  TCanvas c1 ;
  g_retta.Draw ("AP") ;
  c1.Print ("retta.png", "png") ; 

  return 0 ;
}

../_images/linea.png

11.3 fit dei punti generati¶

esercizio_03.cpp

/*
c++ -o esercizio_03 esercizio_03.cpp funzioni.cc `root-config --glibs --cflags`
*/

#include <iostream>
#include <vector>
#include <cmath>

#include "TF1.h"
#include "TGraphErrors.h"
#include "TCanvas.h"
#include "TFitResult.h"

#include "funzioni.h"


using namespace std ;

int main (int argc, char ** argv)
{

  int nPoints = 10 ;
  double min = 0.5 ;
  double max = 10.5 ;
  double epsilon = 0.15 ;

  TF1 f_true ("f_true", retta, 0., 2. * nPoints, 2) ;
  f_true.SetParameter (0, 0.4) ;
  f_true.SetParameter (1, 1.5) ;

  vector<double> v_x ;
  vector<double> v_y ;

  double x = 0. ;
  for (int i = 0 ; i < nPoints ; ++i)
    {
      x+= rand_range (0., 1.) ;
      v_x.push_back (x) ;
      double y = f_true.Eval (x) + rand_TCL (0, epsilon) ;
      v_y.push_back (y) ;
    }
  
  // errore sulle coordinate x dei punti
  vector<double> v_ex (v_x.size (), 0.) ;
  // errore sulle coordinate y dei punti
  vector<double> v_ey (v_y.size (), epsilon);

  TGraphErrors g_retta (v_x.size (), &v_x[0], &v_y[0], &v_ex[0], &v_ey[0]) ;
  g_retta.SetMarkerStyle (4) ;
  g_retta.SetMarkerColor (kRed) ;

  TF1 f_fit ("f_fit", retta, 0., 2. * nPoints, 2) ;
  TFitResultPtr fit_result = g_retta.Fit (&f_fit, "S") ;

  cout << endl ;
  cout.precision (3) ;
  cout << "risultato del fit: " << fit_result->IsValid () << endl ;
  cout << "termine noto : " << f_fit.GetParameter (0) << "\t+- " << f_fit.GetParError (0) << endl ;
  cout << "pendenza     : " << f_fit.GetParameter (1) << "\t+- " << f_fit.GetParError (1) << endl ;

  TCanvas c1 ("c1", "", 800, 800) ;
  g_retta.Draw ("AP") ;
  c1.Print ("retta.png", "png") ; 

  return 0 ;
}

../_images/linea.png

11.4 il calcolo del Q<sup>2<sub> minimo¶

esercizio_04.cpp

/*
c++ -o esercizio_04 esercizio_04.cpp funzioni.cc `root-config --glibs --cflags`
*/

#include <iostream>
#include <vector>
#include <cmath>

#include "TF1.h"
#include "TGraphErrors.h"
#include "TCanvas.h"
#include "TFitResult.h"

#include "funzioni.h"


using namespace std ;

int main (int argc, char ** argv)
{

  int nPoints = 10 ;
  double min = 0.5 ;
  double max = 10.5 ;
  double epsilon = 0.15 ;

  TF1 f_true ("f_true", retta, 0., 2. * nPoints, 2) ;
  f_true.SetParameter (0, 0.4) ;
  f_true.SetParameter (1, 1.5) ;

  vector<double> v_x ;
  vector<double> v_y ;

  double x = 0. ;
  for (int i = 0 ; i < nPoints ; ++i)
    {
      x+= rand_range (0., 1.) ;
      v_x.push_back (x) ;
      double y = f_true.Eval (x) + rand_TCL (0, epsilon) ;
      v_y.push_back (y) ;
    }
  
  // errore sulle coordinate x dei punti
  vector<double> v_ex (v_x.size (), 0.) ;
  // errore sulle coordinate y dei punti
  vector<double> v_ey (v_y.size (), epsilon);

  TGraphErrors g_retta (v_x.size (), &v_x[0], &v_y[0], &v_ex[0], &v_ey[0]) ;
  g_retta.SetMarkerStyle (4) ;
  g_retta.SetMarkerColor (kRed) ;

  TF1 f_fit ("f_fit", retta, 0., 2. * nPoints, 2) ;
  TFitResultPtr fit_result = g_retta.Fit (&f_fit, "S") ;

  cout << endl ;
  cout.precision (3) ;
  cout << "risultato del fit: " << fit_result->IsValid () << endl ;
  cout << "termine noto : " << f_fit.GetParameter (0) << "\t+- " << f_fit.GetParError (0) << endl ;
  cout << "pendenza     : " << f_fit.GetParameter (1) << "\t+- " << f_fit.GetParError (1) << endl ;

  cout << "numero di gradi di libertà: " << fit_result->Ndf () << endl ;
  cout << "valore di Q quadro: " << fit_result->Chi2 () << endl ;

  double Q2 = 0 ;
  for (int i = 0 ; i < v_x.size () ; ++i)
    {
      Q2 += (v_y.at (i) - f_fit.Eval (v_x.at (i))) * (v_y.at (i) - f_fit.Eval (v_x.at (i))) / 
            (v_ey.at (i) * v_ey.at (i)) ; 
    }

  cout << "valore di Q quadro ricalcolato: " << Q2 << endl ;


  TCanvas c1 ;
  g_retta.Draw ("AP") ;
  c1.Print ("retta.png", "png") ; 

  return 0 ;
}

../_images/linea.png

11.5 la distribuzione del Q<sup>2<sub> minimo con toy-experiment¶

esercizio_05.cpp

/*
c++ -o esercizio_05 esercizio_05.cpp funzioni.cc `root-config --glibs --cflags`
*/

#include <iostream>
#include <vector>
#include <cmath>

#include "TF1.h"
#include "TGraphErrors.h"
#include "TCanvas.h"
#include "TFitResult.h"
#include "TH1F.h"

#include "funzioni.h"


using namespace std ;

int main (int argc, char ** argv)
{

  int nPoints = 10 ;
  int nToys = 10000 ;
  double min = 0.5 ;
  double max = 10.5 ;
  double epsilon = 0.15 ;

  TF1 f_true ("f_true", retta, 0., 2. * nPoints, 2) ;
  f_true.SetParameter (0, 0.4) ;
  f_true.SetParameter (1, 1.5) ;

  vector<double> v_x ;
  vector<double> v_y ;
  // errore sulle coordinate x dei punti
  vector<double> v_ex (nPoints, 0.) ;
  // errore sulle coordinate y dei punti
  vector<double> v_ey (nPoints, epsilon);

  TH1F h_Q2 ("h_Q2", "histogram of Q2", 20 * nPoints, 0., 4 * nPoints) ;

  // loop di toy experiment
  for (int iToy = 0 ; iToy < nToys ; ++iToy)
    {
      v_x.clear () ;
      v_y.clear () ;
      double x = 0. ;
      for (int i = 0 ; i < nPoints ; ++i)
        {
          x+= rand_range (0., 1.) ;
          v_x.push_back (x) ;
          double y = f_true.Eval (x) + rand_TCL (0, epsilon) ;
          v_y.push_back (y) ;
        }
      
      TGraphErrors g_retta (v_x.size (), &v_x[0], &v_y[0], &v_ex[0], &v_ey[0]) ;
      TF1 f_fit ("f_fit", retta, 0., 2. * nPoints, 2) ;
      TFitResultPtr fit_result = g_retta.Fit (&f_fit, "SQ") ;
      h_Q2.Fill (fit_result->Chi2 ()) ;
    
    } // loop di toy experiment

  TCanvas c1 ;
  h_Q2.SetFillColor (kOrange) ;
  h_Q2.Draw ("hist") ;
  c1.Print ("Q2.png", "png") ; 

  return 0 ;
}

../_images/linea.png

11.6 la determinazione dell’incetezza sulle y<sub>i</sub>¶

esercizio_06.cpp

/*
c++ -o esercizio_06 esercizio_06.cpp funzioni.cc `root-config --glibs --cflags`
*/

#include <iostream>
#include <vector>
#include <cmath>

#include "TF1.h"
#include "TGraphErrors.h"
#include "TCanvas.h"
#include "TFitResult.h"
#include "TH1F.h"

#include "funzioni.h"


using namespace std ;

int main (int argc, char ** argv)
{

  int nPoints = 10 ;
  int nToys = 10000 ;
  double min = 0.5 ;
  double max = 10.5 ;
  double epsilon = 0.15 ;
  double modificatore = 0.33 ;

  TF1 f_true ("f_true", retta, 0., 2. * nPoints, 2) ;
  f_true.SetParameter (0, 0.4) ;
  f_true.SetParameter (1, 1.5) ;

  vector<double> v_x ;
  vector<double> v_y ;
  // errore sulle coordinate x dei punti
  vector<double> v_ex (nPoints, 0.) ;
  // errore sulle coordinate y dei punti
  vector<double> v_ey (nPoints, epsilon);

  TH1F h_Q2 ("h_Q2", "histogram of Q2", 20 * nPoints, 0., 4 * nPoints) ;
  int ndf = 0 ;

  // loop di toy experiment
  for (int iToy = 0 ; iToy < nToys ; ++iToy)
    {
      v_x.clear () ;
      v_y.clear () ;
      double x = 0. ;
      for (int i = 0 ; i < nPoints ; ++i)
        {
          x+= rand_range (0., 1.) ;
          v_x.push_back (x) ;
          double y = f_true.Eval (x) + rand_TCL (0, modificatore * epsilon) ;
          v_y.push_back (y) ;
        }
      
      TGraphErrors g_retta (v_x.size (), &v_x[0], &v_y[0], &v_ex[0], &v_ey[0]) ;
      TF1 f_fit ("f_fit", retta, 0., 2. * nPoints, 2) ;
      TFitResultPtr fit_result = g_retta.Fit (&f_fit, "SQ") ;
      h_Q2.Fill (fit_result->Chi2 ()) ;
      ndf = fit_result->Ndf () ;
    
    } // loop di toy experiment

  TCanvas c1 ;
  h_Q2.SetFillColor (kOrange) ;
  h_Q2.Draw ("hist") ;
  c1.Print ("Q2.png", "png") ; 

  double scale = h_Q2.GetMean () / ndf ;

  cout << "sigma realistica: " << sqrt (scale) * epsilon << endl ;

  return 0 ;
}

../_images/linea.png

11.7 la matrice di covarianza del fit¶

esercizio_07.cpp

/*
c++ -o esercizio_07 esercizio_07.cpp funzioni.cc `root-config --glibs --cflags`
*/

#include <iostream>
#include <vector>
#include <cmath>

#include "TF1.h"
#include "TGraphErrors.h"
#include "TCanvas.h"
#include "TFitResult.h"

#include "funzioni.h"


using namespace std ;

int main (int argc, char ** argv)
{

  int nPoints = 10 ;
  double min = 0.5 ;
  double max = 10.5 ;
  double epsilon = 0.15 ;

  TF1 f_true ("f_true", retta, 0., 2. * nPoints, 2) ;
  f_true.SetParameter (0, 0.4) ;
  f_true.SetParameter (1, 1.5) ;

  vector<double> v_x ;
  vector<double> v_y ;

  double x = 0. ;
  for (int i = 0 ; i < nPoints ; ++i)
    {
      x+= rand_range (0., 1.) ;
      v_x.push_back (x) ;
      double y = f_true.Eval (x) + rand_TCL (0, epsilon) ;
      v_y.push_back (y) ;
    }
  
  // errore sulle coordinate x dei punti
  vector<double> v_ex (v_x.size (), 0.) ;
  // errore sulle coordinate y dei punti
  vector<double> v_ey (v_y.size (), epsilon);

  TGraphErrors g_retta (v_x.size (), &v_x[0], &v_y[0], &v_ex[0], &v_ey[0]) ;
  g_retta.SetMarkerStyle (4) ;
  g_retta.SetMarkerColor (kRed) ;

  TF1 f_fit ("f_fit", retta, 0., 2. * nPoints, 2) ;
  TFitResultPtr fit_result = g_retta.Fit (&f_fit, "S") ;

  cout << endl ;
  cout.precision (3) ;
  cout << "risultato del fit: " << fit_result->IsValid () << endl ;
  cout << "termine noto : " << f_fit.GetParameter (0) << "\t+- " << f_fit.GetParError (0) << endl ;
  cout << "pendenza     : " << f_fit.GetParameter (1) << "\t+- " << f_fit.GetParError (1) << endl ;

  cout << "\nmatrice di covarianza del risultato\n\n:" ;
  cout << "\t " << fit_result->CovMatrix (0,0) ;
  cout << "\t " << fit_result->CovMatrix (1,0) << "\n" ;
  cout << "\t " << fit_result->CovMatrix (0,1) ;
  cout << "\t " << fit_result->CovMatrix (1,1) << "\n" ; 

  TCanvas c1 ;
  g_retta.Draw ("AP") ;
  c1.Print ("retta.png", "png") ; 

  return 0 ;
}

../_images/linea.png

previous

Esercizi per la lezione 11

next

Lezione 12: fit di distribuzioni binnate con ROOT