Using the RooCustomizer to create multiple PDFs that share a lot of properties, but have unique parameters for each category. As an extra complication, some of the new parameters need to be functions of a mass parameter.
Author: Harshal Shende, Stephan Hageboeck (C++ version)
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Wednesday, April 17, 2024 at 11:19 AM.
import ROOT
E = ROOT.RooRealVar("Energy", "Energy", 0, 3000)
meanG = ROOT.RooRealVar("meanG", "meanG", 100.0, 0.0, 3000.0)
sigmaG = ROOT.RooRealVar("sigmaG", "sigmaG", 3.0)
gauss = ROOT.RooGaussian("gauss", "gauss", E, meanG, sigmaG)
pol1 = ROOT.RooRealVar("pol1", "Constant of the polynomial", 1, -10, 10)
linear = ROOT.RooPolynomial("linear", "linear", E, pol1)
yieldSig = ROOT.RooRealVar("yieldSig", "yieldSig", 1, 0, 1.0e4)
yieldBkg = ROOT.RooRealVar("yieldBkg", "yieldBkg", 1, 0, 1.0e4)
model = ROOT.RooAddPdf("model", "S + B model", [gauss, linear], [yieldSig, yieldBkg])
print("The proto model before customisation:\n")
model.Print("T") # "T" prints the model as a tree
The proto model before customisation:
Build the categories
sample = ROOT.RooCategory("sample", "sample", {"Sample1": 1, "Sample2": 2, "Sample3": 3})
We need two sets for bookkeeping of PDF nodes:
newLeafs = ROOT.RooArgSet()
allCustomiserNodes = ROOT.RooArgSet()
The customiser will make copies of meanG
for each category.
These will all appear in the set newLeafs
, which will own the new nodes.
cust = ROOT.RooCustomizer(model, sample, newLeafs, allCustomiserNodes)
cust.splitArg(meanG, sample)
We need the yields 1 and 2 to be a function of the variable "mass".
For this, we pre-define nodes with exactly the names that the customiser would have created automatically,
that is, "yieldSig
.
mass = ROOT.RooRealVar("M", "M", 1, 0, 12000)
yield1 = ROOT.RooFormulaVar("yieldSig_Sample1", "Signal yield in the first sample", "M/3.360779", mass)
yield2 = ROOT.RooFormulaVar("yieldSig_Sample2", "Signal yield in the second sample", "M/2", mass)
allCustomiserNodes.add(yield1)
allCustomiserNodes.add(yield2)
True
Instruct the customiser to replace all yieldSig nodes for each sample:
cust.splitArg(yieldSig, sample)
Now we can start building the PDFs for all categories:
pdf1 = cust.build("Sample1")
pdf2 = cust.build("Sample2")
pdf3 = cust.build("Sample3")
And we inspect the two PDFs
print("\nPDF 1 with a yield depending on M:\n")
pdf1.Print("T")
print("\nPDF 2 with a yield depending on M:\n")
pdf2.Print("T")
print("\nPDF 3 with a free yield:\n")
pdf3.Print("T")
print("\nThe following leafs have been created automatically while customising:\n")
newLeafs.Print("V")
PDF 1 with a yield depending on M: PDF 2 with a yield depending on M: PDF 3 with a free yield: The following leafs have been created automatically while customising: 0xa7e6440 RooAddPdf::model_Sample1 = 1156.8/1 [Auto,Clean] 0xa7ee440/V- RooGaussian::gauss_Sample1 = 0 [Auto,Dirty] 0x90a0790/V- RooRealVar::Energy = 1500 0xa75bb40/V- RooRealVar::meanG_Sample1 = 100 0x8f8ac30/V- RooRealVar::sigmaG = 3 0xa7c2700/V- RooFormulaVar::yieldSig_Sample1 = 0.29755 [Auto,Clean] 0xa2cd510/V- RooRealVar::M = 1 0x94bdfe0/V- RooPolynomial::linear = 1501 [Auto,Dirty] 0x90a0790/V- RooRealVar::Energy = 1500 0x8bdac60/V- RooRealVar::pol1 = 1 0x64e44a0/V- RooRealVar::yieldBkg = 1 0xa76dd60 RooAddPdf::model_Sample2 = 1000.67/1 [Auto,Clean] 0xa73dee0/V- RooGaussian::gauss_Sample2 = 0 [Auto,Dirty] 0x90a0790/V- RooRealVar::Energy = 1500 0xa7d3e90/V- RooRealVar::meanG_Sample2 = 100 0x8f8ac30/V- RooRealVar::sigmaG = 3 0xa72bcb0/V- RooFormulaVar::yieldSig_Sample2 = 0.5 [Auto,Clean] 0xa2cd510/V- RooRealVar::M = 1 0x94bdfe0/V- RooPolynomial::linear = 1501 [Auto,Dirty] 0x90a0790/V- RooRealVar::Energy = 1500 0x8bdac60/V- RooRealVar::pol1 = 1 0x64e44a0/V- RooRealVar::yieldBkg = 1 0xa83df10 RooAddPdf::model_Sample3 = 750.5/1 [Auto,Clean] 0x8765be0/V- RooGaussian::gauss_Sample3 = 0 [Auto,Dirty] 0x90a0790/V- RooRealVar::Energy = 1500 0xa7d4cb0/V- RooRealVar::meanG_Sample3 = 100 0x8f8ac30/V- RooRealVar::sigmaG = 3 0xa823fc0/V- RooRealVar::yieldSig_Sample3 = 1 0x94bdfe0/V- RooPolynomial::linear = 1501 [Auto,Dirty] 0x90a0790/V- RooRealVar::Energy = 1500 0x8bdac60/V- RooRealVar::pol1 = 1 0x64e44a0/V- RooRealVar::yieldBkg = 1 1) RooRealVar:: meanG_Sample1 = 100 2) RooRealVar:: meanG_Sample2 = 100 3) RooRealVar:: meanG_Sample3 = 100 4) RooRealVar:: yieldSig_Sample3 = 1
If we needed to set reasonable values for the means of the gaussians, this could be done as follows:
meanG1 = allCustomiserNodes["meanG_Sample1"]
meanG1.setVal(200)
meanG2 = allCustomiserNodes["meanG_Sample2"]
meanG2.setVal(300)
print(
"\nThe following leafs 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.):"
)
allCustomiserNodes.Print("V")
The following leafs have been used while customising (partial overlap with the set of automatically created leaves. a new customiser for a different PDF could reuse them if necessary.): 1) RooFormulaVar:: yieldSig_Sample1 = 0.29755 2) RooFormulaVar:: yieldSig_Sample2 = 0.5 3) RooRealVar:: meanG_Sample1 = 200 4) RooRealVar:: meanG_Sample2 = 300 5) RooRealVar:: meanG_Sample3 = 100 6) RooRealVar:: yieldSig_Sample3 = 1