rf201_composite¶
Addition and convolution: composite pdf with signal and background component
pdf = f_bkg * bkg(x,a0,a1) + (1-fbkg) * (f_sig1 * sig1(x,m,s1 + (1-f_sig1) * sig2(x,m,s2)))
Author: Clemens Lange, Wouter Verkerke (C++ version)
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Sunday, February 05, 2023 at 11:14 AM.
In [1]:
import ROOT
Welcome to JupyROOT 6.29/01
Setup component pdfs¶
Declare observable x
In [2]:
x = ROOT.RooRealVar("x", "x", 0, 10)
Create two Gaussian PDFs g1(x,mean1,sigma) anf g2(x,mean2,sigma) and their parameters
In [3]:
mean = ROOT.RooRealVar("mean", "mean of gaussians", 5)
sigma1 = ROOT.RooRealVar("sigma1", "width of gaussians", 0.5)
sigma2 = ROOT.RooRealVar("sigma2", "width of gaussians", 1)
sig1 = ROOT.RooGaussian("sig1", "Signal component 1", x, mean, sigma1)
sig2 = ROOT.RooGaussian("sig2", "Signal component 2", x, mean, sigma2)
[#0] WARNING:InputArguments -- The parameter 'sigma1' with range [-1e+30, 1e+30] of the RooGaussian 'sig1' exceeds the safe range of (0, inf). Advise to limit its range. [#0] WARNING:InputArguments -- The parameter 'sigma2' with range [-1e+30, 1e+30] of the RooGaussian 'sig2' exceeds the safe range of (0, inf). Advise to limit its range.
Build Chebychev polynomial pdf
In [4]:
a0 = ROOT.RooRealVar("a0", "a0", 0.5, 0.0, 1.0)
a1 = ROOT.RooRealVar("a1", "a1", -0.2, 0.0, 1.0)
bkg = ROOT.RooChebychev("bkg", "Background", x, [a0, a1])
Method 1 - Two RooAddPdfs¶
Add signal components
Sum the signal components into a composite signal pdf
In [5]:
sig1frac = ROOT.RooRealVar("sig1frac", "fraction of component 1 in signal", 0.8, 0.0, 1.0)
sig = ROOT.RooAddPdf("sig", "Signal", [sig1, sig2], [sig1frac])
Add signal and background¶
Sum the composite signal and background
In [6]:
bkgfrac = ROOT.RooRealVar("bkgfrac", "fraction of background", 0.5, 0.0, 1.0)
model = ROOT.RooAddPdf("model", "g1+g2+a", [bkg, sig], [bkgfrac])
Sample, fit and plot model¶
Generate a data sample of 1000 events in x from model
In [7]:
data = model.generate({x}, 1000)
Fit model to data
In [8]:
model.fitTo(data)
Out[8]:
<cppyy.gbl.RooFitResult object at 0x(nil)>
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization
[#1] INFO:Minimization -- The following expressions have been identified as constant and will be precalculated and cached: (sig1,sig2)
[#1] INFO:Minimization -- The following expressions will be evaluated in cache-and-track mode: (bkg)
**********
** 1 **SET PRINT 1
**********
**********
** 2 **SET NOGRAD
**********
PARAMETER DEFINITIONS:
NO. NAME VALUE STEP SIZE LIMITS
1 a0 5.00000e-01 1.00000e-01 0.00000e+00 1.00000e+00
2 a1 0.00000e+00 1.00000e-01 0.00000e+00 1.00000e+00
MINUIT WARNING IN PARAM DEF
============== STARTING VALUE IS AT LIMIT.
MINUIT WARNING IN PARAMETR
============== VARIABLE2 IS AT ITS LOWER ALLOWED LIMIT.
MINUIT WARNING IN PARAMETR
============== VARIABLE2 BROUGHT BACK INSIDE LIMITS.
3 bkgfrac 5.00000e-01 1.00000e-01 0.00000e+00 1.00000e+00
4 sig1frac 8.00000e-01 1.00000e-01 0.00000e+00 1.00000e+00
**********
** 3 **SET ERR 0.5
**********
**********
** 4 **SET PRINT 1
**********
**********
** 5 **SET STR 1
**********
NOW USING STRATEGY 1: TRY TO BALANCE SPEED AGAINST RELIABILITY
**********
** 6 **MIGRAD 2000 1
**********
FIRST CALL TO USER FUNCTION AT NEW START POINT, WITH IFLAG=4.
MINUIT WARNING IN MIGrad
============== VARIABLE2 IS AT ITS LOWER ALLOWED LIMIT.
START MIGRAD MINIMIZATION. STRATEGY 1. CONVERGENCE WHEN EDM .LT. 1.00e-03
FCN=1926.2 FROM MIGRAD STATUS=INITIATE 31 CALLS 32 TOTAL
EDM= unknown STRATEGY= 1 NO ERROR MATRIX
EXT PARAMETER CURRENT GUESS STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 a0 5.00000e-01 1.00000e-01 0.00000e+00 2.31850e+00
2 a1 7.50936e-02 1.00000e-01 5.53213e-01 -4.13475e-03
3 bkgfrac 5.00000e-01 1.00000e-01 0.00000e+00 2.67741e+01
4 sig1frac 8.00000e-01 1.00000e-01 0.00000e+00 3.58316e-01
ERR DEF= 0.5
MIGRAD MINIMIZATION HAS CONVERGED.
MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=1924.17 FROM MIGRAD STATUS=CONVERGED 119 CALLS 120 TOTAL
EDM=2.77609e-07 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 a0 5.07642e-01 7.96101e-02 4.71991e-03 -9.40696e-04
2 a1 2.66515e-01 1.34621e-01 5.72412e-03 -6.75205e-05
3 bkgfrac 4.27793e-01 3.57760e-02 1.24946e-03 -1.00877e-02
4 sig1frac 6.41234e-01 9.72166e-02 4.23189e-03 1.73342e-03
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 4 ERR DEF=0.5
6.392e-03 2.434e-03 -5.449e-04 -1.352e-03
2.434e-03 1.871e-02 -3.804e-03 -8.768e-03
-5.449e-04 -3.804e-03 1.282e-03 2.494e-03
-1.352e-03 -8.768e-03 2.494e-03 9.583e-03
PARAMETER CORRELATION COEFFICIENTS
NO. GLOBAL 1 2 3 4
1 0.22587 1.000 0.223 -0.190 -0.173
2 0.79303 0.223 1.000 -0.777 -0.655
3 0.82183 -0.190 -0.777 1.000 0.712
4 0.73007 -0.173 -0.655 0.712 1.000
**********
** 7 **SET ERR 0.5
**********
**********
** 8 **SET PRINT 1
**********
**********
** 9 **HESSE 2000
**********
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=1924.17 FROM HESSE STATUS=OK 23 CALLS 143 TOTAL
EDM=2.75649e-07 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER INTERNAL INTERNAL
NO. NAME VALUE ERROR STEP SIZE VALUE
1 a0 5.07642e-01 7.95738e-02 1.88796e-04 1.52838e-02
2 a1 2.66515e-01 1.33780e-01 2.28965e-04 -4.85861e-01
3 bkgfrac 4.27793e-01 3.55511e-02 2.49893e-04 -1.44921e-01
4 sig1frac 6.41234e-01 9.67457e-02 8.46378e-04 2.86365e-01
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 4 ERR DEF=0.5
6.386e-03 2.397e-03 -5.349e-04 -1.328e-03
2.397e-03 1.847e-02 -3.740e-03 -8.613e-03
-5.349e-04 -3.740e-03 1.266e-03 2.455e-03
-1.328e-03 -8.613e-03 2.455e-03 9.489e-03
PARAMETER CORRELATION COEFFICIENTS
NO. GLOBAL 1 2 3 4
1 0.22394 1.000 0.221 -0.188 -0.171
2 0.78997 0.221 1.000 -0.773 -0.651
3 0.81931 -0.188 -0.773 1.000 0.708
4 0.72690 -0.171 -0.651 0.708 1.000
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
Plot data and PDF overlaid
In [9]:
xframe = x.frame(Title="Example of composite pdf=(sig1+sig2)+bkg")
data.plotOn(xframe)
model.plotOn(xframe)
Out[9]:
<cppyy.gbl.RooPlot object at 0xef5d7c0>
Overlay the background component of model with a dashed line
In [10]:
model.plotOn(xframe, Components={bkg}, LineStyle="--")
Out[10]:
<cppyy.gbl.RooPlot object at 0xef5d7c0>
[#1] INFO:Plotting -- RooAbsPdf::plotOn(model) directly selected PDF components: (bkg) [#1] INFO:Plotting -- RooAbsPdf::plotOn(model) indirectly selected PDF components: ()
Overlay the background+sig2 components of model with a dotted line
In [11]:
model.plotOn(xframe, Components={bkg, sig2}, LineStyle=":")
Out[11]:
<cppyy.gbl.RooPlot object at 0xef5d7c0>
[#1] INFO:Plotting -- RooAbsPdf::plotOn(model) directly selected PDF components: (bkg,sig2) [#1] INFO:Plotting -- RooAbsPdf::plotOn(model) indirectly selected PDF components: (sig)
Print structure of composite pdf
In [12]:
model.Print("t")
0xb648d50 RooAddPdf::model = 0.885987/1 [Auto,Clean]
0xa9c3f70/V- RooChebychev::bkg = 0.733485 [Auto,Dirty]
0xd2790c0/V- RooRealVar::x = 5
0xdd8e4a0/V- RooRealVar::a0 = 0.507642 +/- 0.0795738
0xddc1250/V- RooRealVar::a1 = 0.266515 +/- 0.13378
0x96328c0/V- RooRealVar::bkgfrac = 0.427793 +/- 0.0355511
0x962df70/V- RooAddPdf::sig = 1/1 [Auto,Clean]
0xdd934a0/V- RooGaussian::sig1 = 1 [Auto,Dirty]
0xd2790c0/V- RooRealVar::x = 5
0xda376e0/V- RooRealVar::mean = 5
0x1019740/V- RooRealVar::sigma1 = 0.5
0xe100e20/V- RooRealVar::sig1frac = 0.641234 +/- 0.0967457
0xdddef00/V- RooGaussian::sig2 = 1 [Auto,Dirty]
0xd2790c0/V- RooRealVar::x = 5
0xda376e0/V- RooRealVar::mean = 5
0x10112f0/V- RooRealVar::sigma2 = 1
Method 2 - One RooAddPdf with recursive fractions¶
Construct sum of models on one go using recursive fraction interpretations
model2 = bkg + (sig1 + sig2)
In [13]:
model2 = ROOT.RooAddPdf("model", "g1+g2+a", [bkg, sig1, sig2], [bkgfrac, sig1frac], True)
NB: Each coefficient is interpreted as the fraction of the left-hand component of the i-th recursive sum, i.e.
sum4 = A + ( B + ( C + D) with fraction fA, and fC expands to
sum4 = fAA + (1-fA)(fBB + (1-fB)(fCC + (1-fC)D))
Plot recursive addition model¶
In [14]:
model2.plotOn(xframe, LineColor="r", LineStyle="--")
model2.plotOn(xframe, Components={bkg, sig2}, LineColor="r", LineStyle="--")
model2.Print("t")
[#1] INFO:Plotting -- RooAbsPdf::plotOn(model) directly selected PDF components: (bkg,sig2)
[#1] INFO:Plotting -- RooAbsPdf::plotOn(model) indirectly selected PDF components: ()
0xefd7d30 RooAddPdf::model = 0.885987/1 [Auto,Clean]
0xa9c3f70/V- RooChebychev::bkg = 0.733485 [Auto,Dirty]
0xd2790c0/V- RooRealVar::x = 5
0xdd8e4a0/V- RooRealVar::a0 = 0.507642 +/- 0.0795738
0xddc1250/V- RooRealVar::a1 = 0.266515 +/- 0.13378
0x96328c0/V- RooRealVar::bkgfrac = 0.427793 +/- 0.0355511
0xdd934a0/V- RooGaussian::sig1 = 1 [Auto,Dirty]
0xd2790c0/V- RooRealVar::x = 5
0xda376e0/V- RooRealVar::mean = 5
0x1019740/V- RooRealVar::sigma1 = 0.5
0xeffbc80/V- RooRecursiveFraction::model_recursive_fraction_sig1_2 = 0.366918 [Auto,Clean]
0xe100e20/V- RooRealVar::sig1frac = 0.641234 +/- 0.0967457
0x96328c0/V- RooRealVar::bkgfrac = 0.427793 +/- 0.0355511
0xdddef00/V- RooGaussian::sig2 = 1 [Auto,Dirty]
0xd2790c0/V- RooRealVar::x = 5
0xda376e0/V- RooRealVar::mean = 5
0x10112f0/V- RooRealVar::sigma2 = 1
0xf04d530/V- RooRecursiveFraction::model_recursive_fraction_sig2_3 = 0.205289 [Auto,Clean]
0xee26590/V- RooConstVar::1 = 1
0xe100e20/V- RooRealVar::sig1frac = 0.641234 +/- 0.0967457
0x96328c0/V- RooRealVar::bkgfrac = 0.427793 +/- 0.0355511
Draw the frame on the canvas
In [15]:
c = ROOT.TCanvas("rf201_composite", "rf201_composite", 600, 600)
ROOT.gPad.SetLeftMargin(0.15)
xframe.GetYaxis().SetTitleOffset(1.4)
xframe.Draw()
c.SaveAs("rf201_composite.png")
Info in <TCanvas::Print>: png file rf201_composite.png has been created
Draw all canvases
In [16]:
from ROOT import gROOT
gROOT.GetListOfCanvases().Draw()
