Fitting Example


In [1]:
import ROOT


Welcome to JupyROOT 6.10/06

Inject into the interpreter the functions.


In [2]:
%%cpp -d
//Define functions for fitting 
// 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]);
}

Construct the histogram containing the input data


In [3]:
nbins = 60
data = [ 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 ]
xlow = 0
xup = 3

histo = ROOT.TH1F('histo', 'Lorentzian Peak on Quadratic Background', nbins, xlow, xup)

for i in range(nbins):
   histo.SetBinContent(i+1, data[i])

Create the function and try to fit it without setting any parameter


In [4]:
%jsroot on
c = ROOT.TCanvas()
nparams = 6
fitFcn = ROOT.TF1('fitFcn', ROOT.fitFunction, xlow, xup, nparams)

histo.Fit(fitFcn)
c.Draw()


 FCN=140.05 FROM HESSE     STATUS=FAILED         17 CALLS         354 TOTAL
                     EDM=4.67928e-18    STRATEGY= 1  ERROR MATRIX UNCERTAINTY 100.0 per cent
  EXT PARAMETER                APPROXIMATE        STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  p0          -8.60512e-01   5.51405e-01   3.19762e-03   2.31987e-09
   2  p1           5.42112e+01   4.03368e-01   2.33915e-03  -3.02547e-09
   3  p2          -1.64868e+01   1.64971e-01   9.56676e-04  -1.04872e-08
   4  p3          -3.35916e+05   4.24264e-01  -0.00000e+00   0.00000e+00
   5  p4           0.00000e+00   1.41421e+00  -0.00000e+00   0.00000e+00
   6  p5           0.00000e+00   1.41421e+00  -0.00000e+00   0.00000e+00

Less than optimal. Set parameters and fit again, draw histogram with error bars


In [5]:
fitFcn.SetParameter(4, 0.2); # width
fitFcn.SetParameter(5, 1);   # peak
histo.Fit(fitFcn)
histo.Draw('E')
c.Draw()


 FCN=58.9284 FROM MIGRAD    STATUS=CONVERGED     330 CALLS         331 TOTAL
                     EDM=3.88837e-07    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                   STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  p0          -8.64911e-01   8.91778e-01   2.08428e-03  -4.57833e-04
   2  p1           4.58438e+01   2.64186e+00   1.52471e-03  -3.74735e-04
   3  p2          -1.33216e+01   9.76821e-01   6.23584e-04  -1.32301e-03
   4  p3           1.38076e+01   2.17658e+00   5.05237e-03  -5.14136e-05
   5  p4           1.72316e-01   3.58134e-02   9.42174e-05   1.38074e-02
   6  p5           9.87284e-01   1.12684e-02   4.23955e-05   5.14294e-02

Much better. Now time to beautify the plot. Construct a TF1 for the background and Lorentzian functions and draw them in the same canvas. We save the fit results and set the parameters of the functions accordingly


In [6]:
pars = fitFcn.GetParameters()
backFcn = ROOT.TF1('backFcn', ROOT.background, xlow, xup, 3)
backFcn.SetLineColor(ROOT.kGreen)
backFcn.Draw('Same')
backFcn.SetParameters(pars[0], pars[1], pars[2])
c.Draw()



In [7]:
signalFcn = ROOT.TF1('signalFcn', ROOT.lorentzianPeak, xlow, xup, 3)
signalFcn.SetLineColor(ROOT.kBlue)
signalFcn.SetParameters(pars[3], pars[4], pars[5])
signalFcn.Draw('Same')
c.Draw()


We can now add a legend


In [11]:
legend = ROOT.TLegend(0.45, 0.65, 0.73, 0.85)
legend.SetTextFont(72)
legend.SetTextSize(0.04)
legend.AddEntry(histo, 'Data', 'LE')
legend.AddEntry(backFcn, 'Background fit', 'L')
legend.AddEntry(signalFcn, 'Signal fit', 'L')
legend.AddEntry(fitFcn, 'Global Fit', 'L')
legend.Draw('Same')
c.Draw()



In [ ]: