Logo ROOT  
Reference Guide
ErrorIntegral.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_fit
3/// \notebook -js
4/// Estimate the error in the integral of a fitted function
5/// taking into account the errors in the parameters resulting from the fit.
6/// The error is estimated also using the correlations values obtained from
7/// the fit
8///
9/// run the macro doing:
10///
11/// ~~~{.cpp}
12/// .x ErrorIntegral.C
13/// ~~~
14///
15/// \macro_image
16/// \macro_output
17/// \macro_code
18///
19/// \author Lorenzo Moneta
20
21#include "TF1.h"
22#include "TH1D.h"
23#include "TVirtualFitter.h"
24#include "TMath.h"
25#include <assert.h>
26#include <iostream>
27#include <cmath>
28
29TF1 * fitFunc; // fit function pointer
30
31const int NPAR = 2; // number of function parameters;
32
33//____________________________________________________________________
34double f(double * x, double * p) {
35 // function used to fit the data
36 return p[1]*TMath::Sin( p[0] * x[0] );
37}
38
39//____________________________________________________________________
40void ErrorIntegral() {
41 fitFunc = new TF1("f",f,0,1,NPAR);
42 TH1D * h1 = new TH1D("h1","h1",50,0,1);
43
44 double par[NPAR] = { 3.14, 1.};
45 fitFunc->SetParameters(par);
46
47 h1->FillRandom("f",1000); // fill histogram sampling fitFunc
48 fitFunc->SetParameter(0,3.); // vary a little the parameters
49 h1->Fit(fitFunc); // fit the histogram
50
51 h1->Draw();
52
53 /* calculate the integral*/
54 double integral = fitFunc->Integral(0,1);
55
57 assert(fitter != 0);
58 double * covMatrix = fitter->GetCovarianceMatrix();
59
60 /* using new function in TF1 (from 12/6/2007)*/
61 double sigma_integral = fitFunc->IntegralError(0,1);
62
63 std::cout << "Integral = " << integral << " +/- " << sigma_integral
64 << std::endl;
65
66 // estimated integral and error analytically
67
68 double * p = fitFunc->GetParameters();
69 double ic = p[1]* (1-std::cos(p[0]) )/p[0];
70 double c0c = p[1] * (std::cos(p[0]) + p[0]*std::sin(p[0]) -1.)/p[0]/p[0];
71 double c1c = (1-std::cos(p[0]) )/p[0];
72
73 // estimated error with correlations
74 double sic = std::sqrt( c0c*c0c * covMatrix[0] + c1c*c1c * covMatrix[3]
75 + 2.* c0c*c1c * covMatrix[1]);
76
77 if ( std::fabs(sigma_integral-sic) > 1.E-6*sic )
78 std::cout << " ERROR: test failed : different analytical integral : "
79 << ic << " +/- " << sic << std::endl;
80}
#define f(i)
Definition: RSha256.hxx:104
double cos(double)
double sqrt(double)
double sin(double)
1-Dim function class
Definition: TF1.h:210
virtual Double_t Integral(Double_t a, Double_t b, Double_t epsrel=1.e-12)
IntegralOneDim or analytical integral.
Definition: TF1.cxx:2505
virtual Double_t * GetParameters() const
Definition: TF1.h:514
virtual void SetParameters(const Double_t *params)
Definition: TF1.h:638
virtual void SetParameter(Int_t param, Double_t value)
Definition: TF1.h:628
virtual Double_t IntegralError(Double_t a, Double_t b, const Double_t *params=0, const Double_t *covmat=0, Double_t epsilon=1.E-2)
Return Error on Integral of a parametric function between a and b due to the parameter uncertainties.
Definition: TF1.cxx:2692
1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:614
virtual void FillRandom(const char *fname, Int_t ntimes=5000)
Fill histogram following distribution in function fname.
Definition: TH1.cxx:3445
virtual TFitResultPtr Fit(const char *formula, Option_t *option="", Option_t *goption="", Double_t xmin=0, Double_t xmax=0)
Fit histogram with function fname.
Definition: TH1.cxx:3808
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2998
Abstract Base Class for Fitting.
virtual Double_t * GetCovarianceMatrix() const =0
static TVirtualFitter * GetFitter()
static: return the current Fitter
Double_t x[n]
Definition: legend1.C:17
TH1F * h1
Definition: legend1.C:5
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
Double_t Sin(Double_t)
Definition: TMath.h:627