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