Logo ROOT  
Reference Guide
rf303_conditional.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## 'MULTIDIMENSIONAL MODELS' RooFit tutorial macro #303
5## Use of tailored p.d.f as conditional p.d.fs.s
6##
7## pdf = gauss(x,f(y),sx | y ) with f(y) = a0 + a1*y
8##
9## \macro_code
10##
11## \date February 2018
12## \authors Clemens Lange, Wouter Verkerke (C version)
13
14import ROOT
15
16
17def makeFakeDataXY():
18 x = ROOT.RooRealVar("x", "x", -10, 10)
19 y = ROOT.RooRealVar("y", "y", -10, 10)
20 coord = {x, y}
21
22 d = ROOT.RooDataSet("d", "d", coord)
23
24 for i in range(10000):
25 tmpy = ROOT.gRandom.Gaus(0, 10)
26 tmpx = ROOT.gRandom.Gaus(0.5 * tmpy, 1)
27 if (abs(tmpy) < 10) and (abs(tmpx) < 10):
28 x.setVal(tmpx)
29 y.setVal(tmpy)
30 d.add(coord)
31
32 return d
33
34
35# Set up composed model gauss(x, m(y), s)
36# -----------------------------------------------------------------------
37
38# Create observables
39x = ROOT.RooRealVar("x", "x", -10, 10)
40y = ROOT.RooRealVar("y", "y", -10, 10)
41
42# Create function f(y) = a0 + a1*y
43a0 = ROOT.RooRealVar("a0", "a0", -0.5, -5, 5)
44a1 = ROOT.RooRealVar("a1", "a1", -0.5, -1, 1)
45fy = ROOT.RooPolyVar("fy", "fy", y, [a0, a1])
46
47# Creat gauss(x,f(y),s)
48sigma = ROOT.RooRealVar("sigma", "width of gaussian", 0.5, 0.1, 2.0)
49model = ROOT.RooGaussian("model", "Gaussian with shifting mean", x, fy, sigma)
50
51# Obtain fake external experimental dataset with values for x and y
52expDataXY = makeFakeDataXY()
53
54# Generate data from conditional p.d.f. model(x|y)
55# ---------------------------------------------------------------------------------------------
56
57# Make subset of experimental data with only y values
58expDataY = expDataXY.reduce({y})
59
60# Generate 10000 events in x obtained from _conditional_ model(x|y) with y
61# values taken from experimental data
62data = model.generate({x}, ProtoData=expDataY)
63data.Print()
64
65# Fit conditional p.d.f model(x|y) to data
66# ---------------------------------------------------------------------------------------------
67
68model.fitTo(expDataXY, ConditionalObservables={y})
69
70# Project conditional p.d.f on x and y dimensions
71# ---------------------------------------------------------------------------------------------
72
73# Plot x distribution of data and projection of model x = 1/Ndata
74# sum(data(y_i)) model(x;y_i)
75xframe = x.frame()
76expDataXY.plotOn(xframe)
77model.plotOn(xframe, ProjWData=expDataY)
78
79# Speed up (and approximate) projection by using binned clone of data for
80# projection
81binnedDataY = expDataY.binnedClone()
82model.plotOn(xframe, ProjWData=binnedDataY, LineColor="c", LineStyle=":")
83
84# Show effect of projection with too coarse binning
85(expDataY.get().find("y")).setBins(5)
86binnedDataY2 = expDataY.binnedClone()
87model.plotOn(xframe, ProjWData=binnedDataY2, LineColor="r")
88
89# Make canvas and draw ROOT.RooPlots
90c = ROOT.TCanvas("rf303_conditional", "rf303_conditional", 600, 460)
91ROOT.gPad.SetLeftMargin(0.15)
92xframe.GetYaxis().SetTitleOffset(1.2)
93xframe.Draw()
94
95c.SaveAs("rf303_conditional.png")
RVec< PromoteType< T > > abs(const RVec< T > &v)
Definition: RVec.hxx:1756