Logo ROOT  
Reference Guide
rf603_multicpu.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4##
5## Likelihood and minimization: setting up a multi-core parallelized unbinned maximum likelihood fit
6##
7## \macro_code
8##
9## \date February 2018
10## \author Clemens Lange, Wouter Verkerke (C++ version)
11
12import ROOT
13
14
15# Create 3D pdf and data
16# -------------------------------------------
17
18# Create observables
19x = ROOT.RooRealVar("x", "x", -5, 5)
20y = ROOT.RooRealVar("y", "y", -5, 5)
21z = ROOT.RooRealVar("z", "z", -5, 5)
22
23# Create signal pdf gauss(x)*gauss(y)*gauss(z)
24gx = ROOT.RooGaussian(
25 "gx", "gx", x, ROOT.RooFit.RooConst(0), ROOT.RooFit.RooConst(1))
26gy = ROOT.RooGaussian(
27 "gy", "gy", y, ROOT.RooFit.RooConst(0), ROOT.RooFit.RooConst(1))
28gz = ROOT.RooGaussian(
29 "gz", "gz", z, ROOT.RooFit.RooConst(0), ROOT.RooFit.RooConst(1))
30sig = ROOT.RooProdPdf("sig", "sig", ROOT.RooArgList(gx, gy, gz))
31
32# Create background pdf poly(x)*poly(y)*poly(z)
33px = ROOT.RooPolynomial("px", "px", x, ROOT.RooArgList(
34 ROOT.RooFit.RooConst(-0.1), ROOT.RooFit.RooConst(0.004)))
35py = ROOT.RooPolynomial("py", "py", y, ROOT.RooArgList(
36 ROOT.RooFit.RooConst(0.1), ROOT.RooFit.RooConst(-0.004)))
37pz = ROOT.RooPolynomial("pz", "pz", z)
38bkg = ROOT.RooProdPdf("bkg", "bkg", ROOT.RooArgList(px, py, pz))
39
40# Create composite pdf sig+bkg
41fsig = ROOT.RooRealVar("fsig", "signal fraction", 0.1, 0., 1.)
42model = ROOT.RooAddPdf("model", "model", ROOT.RooArgList(
43 sig, bkg), ROOT.RooArgList(fsig))
44
45# Generate large dataset
46data = model.generate(ROOT.RooArgSet(x, y, z), 200000)
47
48# Parallel fitting
49# -------------------------------
50
51# In parallel mode the likelihood calculation is split in N pieces,
52# that are calculated in parallel and added a posteriori before passing
53# it back to MINUIT.
54
55# Use four processes and time results both in wall time and CPU time
56model.fitTo(data, ROOT.RooFit.NumCPU(4), ROOT.RooFit.Timer(ROOT.kTRUE))
57
58# Parallel MC projections
59# ----------------------------------------------
60
61# Construct signal, likelihood projection on (y,z) observables and
62# likelihood ratio
63sigyz = sig.createProjection(ROOT.RooArgSet(x))
64totyz = model.createProjection(ROOT.RooArgSet(x))
65llratio_func = ROOT.RooFormulaVar(
66 "llratio", "log10(@0)-log10(@1)", ROOT.RooArgList(sigyz, totyz))
67
68# Calculate likelihood ratio for each event, subset of events with high
69# signal likelihood
70data.addColumn(llratio_func)
71dataSel = data.reduce(ROOT.RooFit.Cut("llratio>0.7"))
72
73# Make plot frame and plot data
74frame = x.frame(ROOT.RooFit.Title(
75 "Projection on X with LLratio(y,z)>0.7"), ROOT.RooFit.Bins(40))
76dataSel.plotOn(frame)
77
78# Perform parallel projection using MC integration of pdf using given input dataSet.
79# In self mode the data-weighted average of the pdf is calculated by splitting the
80# input dataset in N equal pieces and calculating in parallel the weighted average
81# one each subset. The N results of those calculations are then weighted into the
82# final result
83
84# Use four processes
85model.plotOn(frame, ROOT.RooFit.ProjWData(dataSel), ROOT.RooFit.NumCPU(4))
86
87c = ROOT.TCanvas("rf603_multicpu", "rf603_multicpu", 600, 600)
88ROOT.gPad.SetLeftMargin(0.15)
89frame.GetYaxis().SetTitleOffset(1.6)
90frame.Draw()
91
92c.SaveAs("rf603_multicpu.png")