Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf514_RooCustomizer.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook -nodraw
4## Using the RooCustomizer to create multiple PDFs that share a lot of properties, but have unique parameters for each category.
5## As an extra complication, some of the new parameters need to be functions of a mass parameter.
6##
7## \macro_code
8## \macro_output
9##
10## \date June 2021
11## \author Harshal Shende, Stephan Hageboeck (C++ version)
12
13import ROOT
14
15E = ROOT.RooRealVar("Energy", "Energy", 0, 3000)
16
17meanG = ROOT.RooRealVar("meanG", "meanG", 100.0, 0.0, 3000.0)
18sigmaG = ROOT.RooRealVar("sigmaG", "sigmaG", 3.0)
19gauss = ROOT.RooGaussian("gauss", "gauss", E, meanG, sigmaG)
20
21pol1 = ROOT.RooRealVar("pol1", "Constant of the polynomial", 1, -10, 10)
22linear = ROOT.RooPolynomial("linear", "linear", E, pol1)
23
24yieldSig = ROOT.RooRealVar("yieldSig", "yieldSig", 1, 0, 1.0e4)
25yieldBkg = ROOT.RooRealVar("yieldBkg", "yieldBkg", 1, 0, 1.0e4)
26
27model = ROOT.RooAddPdf("model", "S + B model", [gauss, linear], [yieldSig, yieldBkg])
28
29print("The proto model before customisation:\n")
30model.Print("T") # "T" prints the model as a tree
31
32
33# Build the categories
34sample = ROOT.RooCategory("sample", "sample", {"Sample1": 1, "Sample2": 2, "Sample3": 3})
35
36
37# Start to customise the proto model that was defined above.
38# ---------------------------------------------------------------------------
39
40# We need two sets for bookkeeping of PDF nodes:
41newLeaves = ROOT.RooArgSet()
42allCustomiserNodes = ROOT.RooArgSet()
43
44
45# 1. Each sample should have its own mean for the gaussian
46# The customiser will make copies of `meanG` for each category.
47# These will all appear in the set `newLeaves`, which will own the new nodes.
48cust = ROOT.RooCustomizer(model, sample, newLeaves, allCustomiserNodes)
49cust.splitArg(meanG, sample)
50
51
52# 2. Each sample should have its own signal yield, but there is an extra complication:
53# We need the yields 1 and 2 to be a function of the variable "mass".
54# For this, we pre-define nodes with exactly the names that the customiser would have created automatically,
55# that is, "<nodeName>_<categoryName>", and we register them in the set of customiser nodes.
56# The customiser will pick them up instead of creating new ones.
57# If we don't provide one (e.g. for "yieldSig_Sample3"), it will be created automatically by cloning `yieldSig`.
58mass = ROOT.RooRealVar("M", "M", 1, 0, 12000)
59yield1 = ROOT.RooFormulaVar("yieldSig_Sample1", "Signal yield in the first sample", "M/3.360779", mass)
60yield2 = ROOT.RooFormulaVar("yieldSig_Sample2", "Signal yield in the second sample", "M/2", mass)
61allCustomiserNodes.add(yield1)
62allCustomiserNodes.add(yield2)
63
64# Instruct the customiser to replace all yieldSig nodes for each sample:
65cust.splitArg(yieldSig, sample)
66
67
68# Now we can start building the PDFs for all categories:
69pdf1 = cust.build("Sample1")
70pdf2 = cust.build("Sample2")
71pdf3 = cust.build("Sample3")
72
73# And we inspect the two PDFs
74print("\nPDF 1 with a yield depending on M:\n")
75pdf1.Print("T")
76print("\nPDF 2 with a yield depending on M:\n")
77pdf2.Print("T")
78print("\nPDF 3 with a free yield:\n")
79pdf3.Print("T")
80
81print("\nThe following leaves have been created automatically while customising:\n")
82newLeaves.Print("V")
83
84# If we needed to set reasonable values for the means of the gaussians, this could be done as follows:
85meanG1 = allCustomiserNodes["meanG_Sample1"]
86meanG1.setVal(200)
87meanG2 = allCustomiserNodes["meanG_Sample2"]
88meanG2.setVal(300)
89
90print(
91 "\nThe following leaves have been used while customising\n\t(partial overlap with the set of automatically created leaves.\n\ta new customiser for a different PDF could reuse them if necessary.):"
92)
93allCustomiserNodes.Print("V")