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