Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf101_basics.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## This tutorial illustrates the basic features of RooFit.
5##
6## \macro_image
7## \macro_code
8## \macro_output
9##
10## \date February 2018
11## \authors Clemens Lange, Wouter Verkerke (C++ version)
12
13import ROOT
14
15# Set up model
16# ---------------------
17# Declare variables x,mean,sigma with associated name, title, initial
18# value and allowed range
19x = ROOT.RooRealVar("x", "x", -10, 10)
20mean = ROOT.RooRealVar("mean", "mean of gaussian", 1, -10, 10)
21sigma = ROOT.RooRealVar("sigma", "width of gaussian", 1, 0.1, 10)
22
23# Build gaussian pdf in terms of x,mean and sigma
24gauss = ROOT.RooGaussian("gauss", "gaussian PDF", x, mean, sigma)
25
26# Construct plot frame in 'x'
27xframe = x.frame(Title="Gaussian pdf") # RooPlot
28
29# Plot model and change parameter values
30# ---------------------------------------------------------------------------
31# Plot gauss in frame (i.e. in x)
32gauss.plotOn(xframe)
33
34# Change the value of sigma to 3
35sigma.setVal(3)
36
37# Plot gauss in frame (i.e. in x) and draw frame on canvas
38gauss.plotOn(xframe, LineColor="r")
39
40# Generate events
41# -----------------------------
42# Generate a dataset of 1000 events in x from gauss
43data = gauss.generate({x}, 10000) # ROOT.RooDataSet
44
45# Make a second plot frame in x and draw both the
46# data and the pdf in the frame
47xframe2 = x.frame(Title="Gaussian pdf with data") # RooPlot
48data.plotOn(xframe2)
49gauss.plotOn(xframe2)
50
51# Fit model to data
52# -----------------------------
53# Fit pdf to data
54gauss.fitTo(data, PrintLevel=-1)
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")