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