Logo ROOT  
Reference Guide
rf101_basics.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4##
5## This tutorial illustrates the basic features of RooFit.
6##
7## \macro_code
8##
9## \date February 2018
10## \authors Clemens Lange, Wouter Verkerke (C++ version)
11
12import ROOT
13
14# Set up model
15# ---------------------
16# Declare variables x,mean,sigma with associated name, title, initial
17# value and allowed range
18x = ROOT.RooRealVar("x", "x", -10, 10)
19mean = ROOT.RooRealVar("mean", "mean of gaussian", 1, -10, 10)
20sigma = ROOT.RooRealVar("sigma", "width of gaussian", 1, 0.1, 10)
21
22# Build gaussian p.d.f in terms of x,mean and sigma
23gauss = ROOT.RooGaussian("gauss", "gaussian PDF", x, mean, sigma)
24
25# Construct plot frame in 'x'
26xframe = x.frame(ROOT.RooFit.Title("Gaussian p.d.f.")) # RooPlot
27
28# Plot model and change parameter values
29# ---------------------------------------------------------------------------
30# Plot gauss in frame (i.e. in x)
31gauss.plotOn(xframe)
32
33# Change the value of sigma to 3
34sigma.setVal(3)
35
36# Plot gauss in frame (i.e. in x) and draw frame on canvas
37gauss.plotOn(xframe, ROOT.RooFit.LineColor(ROOT.kRed))
38
39# Generate events
40# -----------------------------
41# Generate a dataset of 1000 events in x from gauss
42data = gauss.generate(ROOT.RooArgSet(x), 10000) # ROOT.RooDataSet
43
44# Make a second plot frame in x and draw both the
45# data and the p.d.f in the frame
46xframe2 = x.frame(ROOT.RooFit.Title(
47 "Gaussian p.d.f. with data")) # RooPlot
48data.plotOn(xframe2)
49gauss.plotOn(xframe2)
50
51# Fit model to data
52# -----------------------------
53# Fit pdf to data
54gauss.fitTo(data)
55
56# Print values of mean and sigma (that now reflect fitted values and
57# errors)
58mean.Print()
59sigma.Print()
60
61# Draw all frames on a canvas
62c = ROOT.TCanvas("rf101_basics", "rf101_basics", 800, 400)
63c.Divide(2)
64c.cd(1)
65ROOT.gPad.SetLeftMargin(0.15)
66xframe.GetYaxis().SetTitleOffset(1.6)
67xframe.Draw()
68c.cd(2)
69ROOT.gPad.SetLeftMargin(0.15)
70xframe2.GetYaxis().SetTitleOffset(1.6)
71xframe2.Draw()
72
73c.SaveAs("rf101_basics.png")