ROOT   Reference Guide
Searching...
No Matches
rf305_condcorrprod.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## Multidimensional models: multi-dimensional pdfs with conditional pdfs in product
5##
6## pdf = gauss(x,f(y),sx | y ) * gauss(y,ms,sx) with f(y) = a0 + a1*y
7##
8## \macro_image
9## \macro_code
10## \macro_output
11##
12## \date February 2018
13## \authors Clemens Lange, Wouter Verkerke (C++ version)
14
15import ROOT
16
17# Create conditional pdf gx(x|y)
18# -----------------------------------------------------------
19
20# Create observables
21x = ROOT.RooRealVar("x", "x", -5, 5)
22y = ROOT.RooRealVar("y", "y", -5, 5)
23
24# Create function f(y) = a0 + a1*y
25a0 = ROOT.RooRealVar("a0", "a0", -0.5, -5, 5)
26a1 = ROOT.RooRealVar("a1", "a1", -0.5, -1, 1)
27fy = ROOT.RooPolyVar("fy", "fy", y, [a0, a1])
28
29# Create gaussx(x,f(y),sx)
30sigmax = ROOT.RooRealVar("sigma", "width of gaussian", 0.5)
31gaussx = ROOT.RooGaussian("gaussx", "Gaussian in x with shifting mean in y", x, fy, sigmax)
32
33# Create pdf gy(y)
34# -----------------------------------------------------------
35
36# Create gaussy(y,0,5)
37gaussy = ROOT.RooGaussian("gaussy", "Gaussian in y", y, 0.0, 3.0)
38
39# Create product gx(x|y)*gy(y)
40# -------------------------------------------------------
41
42# Create gaussx(x,sx|y) * gaussy(y)
43model = ROOT.RooProdPdf("model", "gaussx(x|y)*gaussy(y)", {gaussy}, Conditional=({gaussx}, {x}))
44
45# Sample, fit and plot product pdf
46# ---------------------------------------------------------------
47
48# Generate 1000 events in x and y from model
49data = model.generate({x, y}, 10000)
50
51# Plot x distribution of data and projection of model x = Int(dy)
52# model(x,y)
53xframe = x.frame()
54data.plotOn(xframe)
55model.plotOn(xframe)
56
57# Plot x distribution of data and projection of model y = Int(dx)
58# model(x,y)
59yframe = y.frame()
60data.plotOn(yframe)
61model.plotOn(yframe)
62
63# Make two-dimensional plot in x vs y
64hh_model = model.createHistogram("hh_model", x, ROOT.RooFit.Binning(50), ROOT.RooFit.YVar(y, ROOT.RooFit.Binning(50)))
65hh_model.SetLineColor(ROOT.kBlue)
66
67# Make canvas and draw ROOT.RooPlots
68c = ROOT.TCanvas("rf305_condcorrprod", "rf05_condcorrprod", 1200, 400)
69c.Divide(3)
70c.cd(1)