rf207_comptools

'ADDITION AND CONVOLUTION' RooFit tutorial macro #207 Tools and utilities for manipulation of composite objects

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 Wednesday, November 30, 2022 at 11:22 AM.

In [1]:
import ROOT
Welcome to JupyROOT 6.27/01

Set up composite pdf dataset

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)
sigma = ROOT.RooRealVar("sigma", "width of gaussians", 0.5)
sig = ROOT.RooGaussian("sig", "Signal component 1", x, mean, sigma)
[#0] WARNING:InputArguments -- The parameter 'sigma' with range [-1e+30, 1e+30] of the RooGaussian 'sig' exceeds the safe range of (0, inf). Advise to limit its range.

Build Chebychev polynomial p.d.f.

In [4]:
a0 = ROOT.RooRealVar("a0", "a0", 0.5, 0.0, 1.0)
a1 = ROOT.RooRealVar("a1", "a1", 0.2, 0.0, 1.0)
bkg1 = ROOT.RooChebychev("bkg1", "Background 1", x, [a0, a1])

Build expontential pdf

In [5]:
alpha = ROOT.RooRealVar("alpha", "alpha", -1)
bkg2 = ROOT.RooExponential("bkg2", "Background 2", x, alpha)

Sum the background components into a composite background p.d.f.

In [6]:
bkg1frac = ROOT.RooRealVar("bkg1frac", "fraction of component 1 in background", 0.2, 0.0, 1.0)
bkg = ROOT.RooAddPdf("bkg", "Signal", [bkg1, bkg2], [bkg1frac])

Sum the composite signal and background

In [7]:
bkgfrac = ROOT.RooRealVar("bkgfrac", "fraction of background", 0.5, 0.0, 1.0)
model = ROOT.RooAddPdf("model", "g1+g2+a", [bkg, sig], [bkgfrac])

Create dummy dataset that has more observables than the above pdf

In [8]:
y = ROOT.RooRealVar("y", "y", -10, 10)
data = ROOT.RooDataSet("data", "data", {x, y})

Basic information requests

Get list of observables

Get list of observables of pdf in context of a dataset

Observables are define each context as the variables shared between a model and a dataset. In self case that is the variable 'x'

In [9]:
model_obs = model.getObservables(data)
model_obs.Print("v")
  1) 0x8268420 RooRealVar:: x = 5  L(0 - 10)  "x"

Get list of parameters

Get list of parameters, list of observables

In [10]:
model_params = model.getParameters({x})
model_params.Print("v")
  1) 0x818e4d0 RooRealVar::       a0 = 0.5  L(0 - 1)  "a0"
  2) 0x7d8f5a0 RooRealVar::       a1 = 0.2  L(0 - 1)  "a1"
  3) 0x85815e0 RooRealVar::    alpha = -1 C  L(-INF - +INF)  "alpha"
  4) 0x7cede10 RooRealVar:: bkg1frac = 0.2  L(0 - 1)  "fraction of component 1 in background"
  5) 0x17cd9a0 RooRealVar::  bkgfrac = 0.5  L(0 - 1)  "fraction of background"
  6) 0x815adb0 RooRealVar::     mean = 5 C  L(-INF - +INF)  "mean of gaussians"
  7) 0x81b0ca0 RooRealVar::    sigma = 0.5 C  L(-INF - +INF)  "width of gaussians"

Get list of parameters, a dataset (Gives identical results to operation above)

In [11]:
model_params2 = model.getParameters(data)
model_params2.Print()
RooArgSet::parameters = (a0,a1,alpha,bkg1frac,bkgfrac,mean,sigma)

Get list of components

Get list of component objects, top-level node

In [12]:
model_comps = model.getComponents()
model_comps.Print("v")
  1) 0x86d3860 RooAddPdf:: model[ bkgfrac * bkg + [%] * sig ] = 0.582695/1  "g1+g2+a"
  2) 0x876c6b0 RooAddPdf::   bkg[ bkg1frac * bkg1 + [%] * bkg2 ] = 0.16539/1  "Signal"
  3) 0x85d46f0 RooChebychev::  bkg1[ x=x coefficients=(a0,a1) ] = 0.8  "Background 1"
  4) 0x86c7510 RooExponential::  bkg2[ x=x c=alpha ] = 0.00673795  "Background 2"
  5) 0x83d1dd0 RooGaussian::   sig[ x=x mean=mean sigma=sigma ] = 1  "Signal component 1"

Modifications to structure of composites

Create a second Gaussian

In [13]:
sigma2 = ROOT.RooRealVar("sigma2", "width of gaussians", 1)
sig2 = ROOT.RooGaussian("sig2", "Signal component 1", x, mean, sigma2)
[#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.

Create a sum of the original Gaussian plus the second Gaussian

In [14]:
sig1frac = ROOT.RooRealVar("sig1frac", "fraction of component 1 in signal", 0.8, 0.0, 1.0)
sigsum = ROOT.RooAddPdf("sigsum", "sig+sig2", [sig, sig2], [sig1frac])

Construct a customizer utility to customize model

In [15]:
cust = ROOT.RooCustomizer(model, "cust")

Instruct the customizer to replace node 'sig' with node 'sigsum'

In [16]:
cust.replaceArg(sig, sigsum)

Build a clone of the input pdf according to the above customization instructions. Each node that requires modified is clone so that the original pdf remained untouched. The name of each cloned node is that of the original node suffixed by the name of the customizer object

The returned head node own all nodes that were cloned as part of the build process so when cust_clone is deleted so will all other nodes that were created in the process.

In [17]:
cust_clone = cust.build(ROOT.kTRUE)
[#1] INFO:ObjectHandling -- RooCustomizer::build(model): tree node sig will be replaced by sigsum
[#1] INFO:ObjectHandling -- RooCustomizer::build(model) Branch node RooAddPdf::model cloned: depends on a replaced parameter
[#1] INFO:ObjectHandling -- RooCustomizer::build(model) Branch node sig is already replaced

Print structure of clone of model with sig.sigsum replacement.

In [18]:
cust_clone.Print("t")
0x8c0c300 RooAddPdf::model_cust = 0.582695/1 [Auto,Clean] 
  0x876c6b0/V- RooAddPdf::bkg = 0.16539/1 [Auto,Clean] 
    0x85d46f0/V- RooChebychev::bkg1 = 0.8 [Auto,Dirty] 
      0x8268420/V- RooRealVar::x = 5
      0x818e4d0/V- RooRealVar::a0 = 0.5
      0x7d8f5a0/V- RooRealVar::a1 = 0.2
    0x7cede10/V- RooRealVar::bkg1frac = 0.2
    0x86c7510/V- RooExponential::bkg2 = 0.00673795 [Auto,Dirty] 
      0x8268420/V- RooRealVar::x = 5
      0x85815e0/V- RooRealVar::alpha = -1
  0x17cd9a0/V- RooRealVar::bkgfrac = 0.5
  0x85c4670/V- RooAddPdf::sigsum = 1/1 [Auto,Clean] 
    0x83d1dd0/V- RooGaussian::sig = 1 [Auto,Dirty] 
      0x8268420/V- RooRealVar::x = 5
      0x815adb0/V- RooRealVar::mean = 5
      0x81b0ca0/V- RooRealVar::sigma = 0.5
    0x8aeb7a0/V- RooRealVar::sig1frac = 0.8
    0x8bc6bd0/V- RooGaussian::sig2 = 1 [Auto,Dirty] 
      0x8268420/V- RooRealVar::x = 5
      0x815adb0/V- RooRealVar::mean = 5
      0x8412b10/V- RooRealVar::sigma2 = 1

The RooCustomizer has the be deleted first. Otherwise, it might happen that sig or sigsum are deleted first, in which case the internal TLists in the RooCustomizer will complain about deleted objects.

In [19]:
del cust

Draw all canvases

In [20]:
from ROOT import gROOT 
gROOT.GetListOfCanvases().Draw()