Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf107_plotstyles.py File Reference

Detailed Description

View in nbviewer Open in SWAN
Basic functionality: demonstration of various plotting styles of data, functions in a RooPlot

import ROOT
# Set up model
# ---------------------
# Create observables
x = ROOT.RooRealVar("x", "x", -10, 10)
# Create Gaussian
sigma = ROOT.RooRealVar("sigma", "sigma", 3, 0.1, 10)
mean = ROOT.RooRealVar("mean", "mean", -3, -10, 10)
gauss = ROOT.RooGaussian("gauss", "gauss", x, mean, sigma)
# Generate a sample of 100 events with sigma=3
data = gauss.generate({x}, 100)
# Fit pdf to data
gauss.fitTo(data, PrintLevel=-1)
# Make plot frames
# -------------------------------
# Make four plot frames to demonstrate various plotting features
frame1 = x.frame(Name="xframe", Title="Red Curve / SumW2 Histo errors", Bins=20)
frame2 = x.frame(Name="xframe", Title="Dashed Curve / No XError bars", Bins=20)
frame3 = x.frame(Name="xframe", Title="Filled Curve / Blue Histo", Bins=20)
frame4 = x.frame(Name="xframe", Title="Partial Range / Filled Bar chart", Bins=20)
# Data plotting styles
# ---------------------------------------
# Use sqrt(sum(weights^2)) error instead of Poisson errors
data.plotOn(frame1, DataError="SumW2")
# Remove horizontal error bars
data.plotOn(frame2, XErrorSize=0)
# Blue markers and error bors
data.plotOn(frame3, MarkerColor="b", LineColor="b")
# Filled bar chart
data.plotOn(frame4, DrawOption="B", DataError=None, XErrorSize=0, FillColor="kGray")
# Function plotting styles
# -----------------------------------------------
# Change line color to red
gauss.plotOn(frame1, LineColor="r")
# Change line style to dashed
gauss.plotOn(frame2, LineStyle="--")
# Filled shapes in green color
gauss.plotOn(frame3, MoveToBack=True, DrawOption="F", FillColor="kOrange")
#
gauss.plotOn(frame4, Range=(-8, 3), LineColor="m")
c = ROOT.TCanvas("rf107_plotstyles", "rf107_plotstyles", 800, 800)
c.Divide(2, 2)
c.cd(1)
ROOT.gPad.SetLeftMargin(0.15)
frame1.GetYaxis().SetTitleOffset(1.6)
frame1.Draw()
c.cd(2)
ROOT.gPad.SetLeftMargin(0.15)
frame2.GetYaxis().SetTitleOffset(1.6)
frame2.Draw()
c.cd(3)
ROOT.gPad.SetLeftMargin(0.15)
frame3.GetYaxis().SetTitleOffset(1.6)
frame3.Draw()
c.cd(4)
ROOT.gPad.SetLeftMargin(0.15)
frame4.GetYaxis().SetTitleOffset(1.6)
frame4.Draw()
c.SaveAs("rf107_plotstyles.png")
[#1] INFO:Fitting -- RooAbsPdf::fitTo(gauss_over_gauss_Int[x]) fixing normalization set for coefficient determination to observables in data
[#1] INFO:Fitting -- using CPU computation library compiled with -mavx2
[#1] INFO:Fitting -- RooAddition::defaultErrorLevel(nll_gauss_over_gauss_Int[x]_gaussData) Summation contains a RooNLLVar, using its error level
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
[#1] INFO:Plotting -- RooAbsPdf::plotOn(gauss) only plotting range [-8,3], curve is normalized to data in given range
[#1] INFO:Eval -- RooRealVar::setRange(x) new range named 'plotRange' created with bounds [-8,3]
Date
February 2018
Authors
Clemens Lange, Wouter Verkerke (C++ version)

Definition in file rf107_plotstyles.py.