Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf101_basics.C File Reference

Detailed Description

View in nbviewer Open in SWAN
Basic functionality: fitting, plotting, toy data generation on one-dimensional PDFs.

pdf = gauss(x,m,s)

#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "TCanvas.h"
#include "RooPlot.h"
#include "TAxis.h"
using namespace RooFit;
{
// S e t u p m o d e l
// ---------------------
// Declare variables x,mean,sigma with associated name, title, initial value and allowed range
RooRealVar x("x", "x", -10, 10);
RooRealVar mean("mean", "mean of gaussian", 1, -10, 10);
RooRealVar sigma("sigma", "width of gaussian", 1, 0.1, 10);
// Build gaussian pdf in terms of x,mean and sigma
RooGaussian gauss("gauss", "gaussian PDF", x, mean, sigma);
// Construct plot frame in 'x'
RooPlot *xframe = x.frame(Title("Gaussian pdf."));
// P l o t m o d e l a n d c h a n g e p a r a m e t e r v a l u e s
// ---------------------------------------------------------------------------
// Plot gauss in frame (i.e. in x)
gauss.plotOn(xframe);
// Change the value of sigma to 3
sigma.setVal(3);
// Plot gauss in frame (i.e. in x) and draw frame on canvas
gauss.plotOn(xframe, LineColor(kRed));
// G e n e r a t e e v e n t s
// -----------------------------
// Generate a dataset of 1000 events in x from gauss
std::unique_ptr<RooDataSet> data{gauss.generate(x, 10000)};
// Make a second plot frame in x and draw both the
// data and the pdf in the frame
RooPlot *xframe2 = x.frame(Title("Gaussian pdf with data"));
data->plotOn(xframe2);
gauss.plotOn(xframe2);
// F i t m o d e l t o d a t a
// -----------------------------
// Fit pdf to data
gauss.fitTo(*data, PrintLevel(-1));
// Print values of mean and sigma (that now reflect fitted values and errors)
mean.Print();
sigma.Print();
// Draw all frames on a canvas
TCanvas *c = new TCanvas("rf101_basics", "rf101_basics", 800, 400); c->Divide(2);
c->cd(1);
gPad->SetLeftMargin(0.15);
xframe->GetYaxis()->SetTitleOffset(1.6);
xframe->Draw();
c->cd(2);
gPad->SetLeftMargin(0.15);
xframe2->GetYaxis()->SetTitleOffset(1.6);
xframe2->Draw();
}
#define c(i)
Definition RSha256.hxx:101
@ kRed
Definition Rtypes.h:66
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
#define gPad
Plain Gaussian p.d.f.
Definition RooGaussian.h:24
Plot frame and a container for graphics objects within that frame.
Definition RooPlot.h:45
static RooPlot * frame(const RooAbsRealLValue &var, double xmin, double xmax, Int_t nBins)
Create a new frame for a given variable in x.
Definition RooPlot.cxx:225
TAxis * GetYaxis() const
Definition RooPlot.cxx:1264
void Draw(Option_t *options=nullptr) override
Draw this plot and all of the elements it contains.
Definition RooPlot.cxx:637
Variable that can be changed from the outside.
Definition RooRealVar.h:37
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition TAttAxis.cxx:298
The Canvas class.
Definition TCanvas.h:23
const Double_t sigma
Double_t x[n]
Definition legend1.C:17
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition JSONIO.h:26
[#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
RooRealVar::mean = 1.01746 +/- 0.0300144 L(-10 - 10)
RooRealVar::sigma = 2.9787 +/- 0.0219217 L(0.1 - 10)
Date
July 2008
Author
Wouter Verkerke

Definition in file rf101_basics.C.