Logo ROOT  
Reference Guide
rf312_multirangefit.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook -nodraw
4##
5## Multidimensional models: performing fits in multiple (disjoint) ranges in one or more dimensions
6##
7## \macro_code
8##
9## \date February 2018
10## \author Clemens Lange, Wouter Verkerke (C++ version)
11
12import ROOT
13
14
15# Create 2D pdf and data
16# -------------------------------------------
17
18# Define observables x,y
19x = ROOT.RooRealVar("x", "x", -10, 10)
20y = ROOT.RooRealVar("y", "y", -10, 10)
21
22# Construct the signal pdf gauss(x)*gauss(y)
23mx = ROOT.RooRealVar("mx", "mx", 1, -10, 10)
24my = ROOT.RooRealVar("my", "my", 1, -10, 10)
25
26gx = ROOT.RooGaussian("gx", "gx", x, mx, ROOT.RooFit.RooConst(1))
27gy = ROOT.RooGaussian("gy", "gy", y, my, ROOT.RooFit.RooConst(1))
28
29sig = ROOT.RooProdPdf("sig", "sig", gx, gy)
30
31# Construct the background pdf (flat in x,y)
32px = ROOT.RooPolynomial("px", "px", x)
33py = ROOT.RooPolynomial("py", "py", y)
34bkg = ROOT.RooProdPdf("bkg", "bkg", px, py)
35
36# Construct the composite model sig+bkg
37f = ROOT.RooRealVar("f", "f", 0., 1.)
38model = ROOT.RooAddPdf(
39 "model", "model", ROOT.RooArgList(
40 sig, bkg), ROOT.RooArgList(f))
41
42# Sample 10000 events in (x,y) from the model
43modelData = model.generate(ROOT.RooArgSet(x, y), 10000)
44
45# Define signal and sideband regions
46# -------------------------------------------------------------------
47
48# Construct the SideBand1,SideBand2, regions
49#
50# |
51# +-------------+-----------+
52# | | |
53# | Side | Sig |
54# | Band1 | nal |
55# | | |
56# --+-------------+-----------+--
57# | |
58# | Side |
59# | Band2 |
60# | |
61# +-------------+-----------+
62# |
63
64x.setRange("SB1", -10, +10)
65y.setRange("SB1", -10, 0)
66
67x.setRange("SB2", -10, 0)
68y.setRange("SB2", 0, +10)
69
70x.setRange("SIG", 0, +10)
71y.setRange("SIG", 0, +10)
72
73x.setRange("FULL", -10, +10)
74y.setRange("FULL", -10, +10)
75
76# Perform fits in individual sideband regions
77# -------------------------------------------------------------------------------------
78
79# Perform fit in SideBand1 region (ROOT.RooAddPdf coefficients will be
80# interpreted in full range)
81r_sb1 = model.fitTo(modelData, ROOT.RooFit.Range(
82 "SB1"), ROOT.RooFit.Save())
83
84# Perform fit in SideBand2 region (ROOT.RooAddPdf coefficients will be
85# interpreted in full range)
86r_sb2 = model.fitTo(modelData, ROOT.RooFit.Range(
87 "SB2"), ROOT.RooFit.Save())
88
89# Perform fits in joint sideband regions
90# -----------------------------------------------------------------------------
91
92# Now perform fit to joint 'L-shaped' sideband region 'SB1|SB2'
93# (ROOT.RooAddPdf coefficients will be interpreted in full range)
94r_sb12 = model.fitTo(modelData, ROOT.RooFit.Range(
95 "SB1,SB2"), ROOT.RooFit.Save())
96
97# Print results for comparison
98r_sb1.Print()
99r_sb2.Print()
100r_sb12.Print()