FittingDemo.C File Reference

Detailed Description

Example for fitting signal/background.

This example can be executed with:

root > .x FittingDemo.C  (using the cling interpreter)
root > .x FittingDemo.C+ (using the native complier via ACLIC)

 
****************************************
Minimizer is Minuit2 / Migrad
Chi2                      =      58.9284
NDf                       =           54
Edm                       =  9.73335e-07
NCalls                    =          606
p0                        =    -0.864746   +/-   0.887832    
p1                        =      45.8428   +/-   2.65453     
p2                        =     -13.3213   +/-   0.980305    
p3                        =      13.8087   +/-   2.24457     
p4                        =     0.172313   +/-   0.0374057   
p5                        =     0.987278   +/-   0.0112928   
Minuit2Minimizer: Minimize with max-calls 1780 convergence for edm < 0.01 strategy 1
Minuit2Minimizer : Valid minimum - status = 0
FVAL  = 58.9283860619491762
Edm   = 4.48436860430199946e-13
Nfcn  = 163
p0   = -0.864713   +/-  0.891795
p1   = 45.8434  +/-  2.64223
p2   = -13.3214    +/-  0.976969
p3   = 13.8074  +/-  2.1776
p4   = 0.172309    +/-  0.0358282
p5   = 0.987281    +/-  0.0112683
 
Covariance Matrix:
 
                         p0          p1          p2          p3          p4          p5
p0                   0.7953     -1.2054     0.34842    -0.15946  -0.0037284  0.00042265
p1                  -1.2054      6.9814     -2.5255     -3.0272   -0.037043  -0.0018129
p2                  0.34842     -2.5255     0.95447      1.1587    0.014337  0.00058758
p3                 -0.15946     -3.0272      1.1587      4.7419     0.05537   0.0017371
p4               -0.0037284   -0.037043    0.014337     0.05537   0.0012837  2.8395e-05
p5               0.00042265  -0.0018129  0.00058758   0.0017371  2.8395e-05  0.00012697
 
Correlation Matrix:
 
                         p0          p1          p2          p3          p4          p5
p0                        1    -0.51155     0.39991    -0.08211    -0.11669    0.042058
p1                 -0.51155           1    -0.97834    -0.52614     -0.3913   -0.060891
p2                  0.39991    -0.97834           1     0.54464     0.40958    0.053373
p3                 -0.08211    -0.52614     0.54464           1      0.7097    0.070794
p4                 -0.11669     -0.3913     0.40958      0.7097           1    0.070334
p5                 0.042058   -0.060891    0.053373    0.070794    0.070334           1
****************************************
Minimizer is Minuit2 / Migrad
Chi2                      =      58.9284
NDf                       =           54
Edm                       =  4.48437e-13
NCalls                    =          163
p0                        =    -0.864713   +/-   0.891795    
p1                        =      45.8434   +/-   2.64223     
p2                        =     -13.3214   +/-   0.976969    
p3                        =      13.8074   +/-   2.1776      
p4                        =     0.172309   +/-   0.0358282   
p5                        =     0.987281   +/-   0.0112683   

 
#include "TH1.h"
#include "TMath.h"
#include "TF1.h"
#include "TLegend.h"
#include "TCanvas.h"
 
// Quadratic background function
double background(double *x, double *par) {
   return par[0] + par[1]*x[0] + par[2]*x[0]*x[0];
}
 
 
// Lorenzian Peak function
double lorentzianPeak(double *x, double *par) {
  return (0.5*par[0]*par[1]/TMath::Pi()) /
    TMath::Max( 1.e-10,(x[0]-par[2])*(x[0]-par[2])
   + .25*par[1]*par[1]);
}
 
// Sum of background and peak function
double fitFunction(double *x, double *par) {
  return background(x,par) + lorentzianPeak(x,&par[3]);
}
 
void FittingDemo() {
 //Bevington Exercise by Peter Malzacher, modified by Rene Brun
 
   const int nBins = 60;
 
   double data[nBins] = { 6, 1,10,12, 6,13,23,22,15,21,
                           23,26,36,25,27,35,40,44,66,81,
                           75,57,48,45,46,41,35,36,53,32,
                           40,37,38,31,36,44,42,37,32,32,
                           43,44,35,33,33,39,29,41,32,44,
                           26,39,29,35,32,21,21,15,25,15};
   TCanvas *c1 = new TCanvas("c1","Fitting Demo",10,10,700,500);
   c1->SetFillColor(33);
   c1->SetFrameFillColor(41);
   c1->SetGrid();
 
   TH1F *histo = new TH1F("histo",
      "Lorentzian Peak on Quadratic Background",60,0,3);
   histo->SetMarkerStyle(21);
   histo->SetMarkerSize(0.8);
   histo->SetStats(false);
 
   for(int i=0; i < nBins;  i++) histo->SetBinContent(i+1,data[i]);
 
   // create a TF1 with the range from 0 to 3 and 6 parameters
   TF1 *fitFcn = new TF1("fitFcn",fitFunction,0,3,6);
   fitFcn->SetNpx(500);
   fitFcn->SetLineWidth(4);
   fitFcn->SetLineColor(kMagenta);
 
   // first try without starting values for the parameters
   // This defaults to 1 for each param.
   // this results in an ok fit for the polynomial function
   // however the non-linear part (lorenzian) does not
   // respond well.
   fitFcn->SetParameters(1,1,1,1,1,1);
   histo->Fit("fitFcn","0");
 
   // second try: set start values for some parameters
   fitFcn->SetParameter(4,0.2); // width
   fitFcn->SetParameter(5,1);   // peak
 
   histo->Fit("fitFcn","V+","ep");
 
   // improve the picture:
   TF1 *backFcn = new TF1("backFcn",background,0,3,3);
   backFcn->SetLineColor(kRed);
   TF1 *signalFcn = new TF1("signalFcn",lorentzianPeak,0,3,3);
   signalFcn->SetLineColor(kBlue);
   signalFcn->SetNpx(500);
   double par[6];
 
   // writes the fit results into the par array
   fitFcn->GetParameters(par);
 
   backFcn->SetParameters(par);
   backFcn->Draw("same");
 
   signalFcn->SetParameters(&par[3]);
   signalFcn->Draw("same");
 
   // draw the legend
   TLegend *legend=new TLegend(0.6,0.65,0.88,0.85);
   legend->SetTextFont(72);
   legend->SetTextSize(0.04);
   legend->AddEntry(histo,"Data","lpe");
   legend->AddEntry(backFcn,"Background fit","l");
   legend->AddEntry(signalFcn,"Signal fit","l");
   legend->AddEntry(fitFcn,"Global Fit","l");
   legend->Draw();
 
}

Author

Rene Brun

Definition in file FittingDemo.C.