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