Rf 1 1 0_Normintegration

Basic functionality: normalization and integration of pdfs, construction of cumulative distribution monodimensional functions

Author: Wouter Verkerke
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Wednesday, January 19, 2022 at 10:11 AM.

In [1]:
%%cpp -d
#include "RooRealVar.h"
#include "RooGaussian.h"
#include "RooConstVar.h"
#include "RooAbsReal.h"
#include "RooPlot.h"
#include "TCanvas.h"
#include "TAxis.h"
In [2]:
%%cpp -d
// This is a workaround to make sure the namespace is used inside functions
using namespace RooFit;

Setup model

Create observables x,y

In [3]:
RooRealVar x("x", "x", -10, 10);
RooFit v3.60 -- Developed by Wouter Verkerke and David Kirkby 
                Copyright (C) 2000-2013 NIKHEF, University of California & Stanford University
                All rights reserved, please read http://roofit.sourceforge.net/license.txt

Create pdf gaussx(x,-2,3)

In [4]:
RooGaussian gx("gx", "gx", x, RooConst(-2), RooConst(3));

Retrieve raw & normalized values of roofit p.d.f.s

Return 'raw' unnormalized value of gx

In [5]:
cout << "gx = " << gx.getVal() << endl;
gx = 0.800737

Return value of gx normalized over x in range [-10,10]

In [6]:
RooArgSet nset(x);
cout << "gx_Norm[x] = " << gx.getVal(&nset) << endl;
gx_Norm[x] = 0.106896

Create object representing integral over gx which is used to calculate gx_Norm[x] == gx / gx_Int[x]

In [7]:
RooAbsReal *igx = gx.createIntegral(x);
cout << "gx_Int[x] = " << igx->getVal() << endl;
gx_Int[x] = 7.49084

Integrate normalized pdf over subrange

Define a range named "signal" in x from -5,5

In [8]:
x.setRange("signal", -5, 5);
[#1] INFO:Eval -- RooRealVar::setRange(x) new range named 'signal' created with bounds [-5,5]

Create an integral of gx_norm[x] over x in range "signal" This is the fraction of of pdf gx_Norm[x] which is in the range named "signal"

In [9]:
RooAbsReal *igx_sig = gx.createIntegral(x, NormSet(x), Range("signal"));
cout << "gx_Int[x|signal]_Norm[x] = " << igx_sig->getVal() << endl;
gx_Int[x|signal]_Norm[x] = 0.834753

Construct cumulative distribution function from pdf

Create the cumulative distribution function of gx i.e. calculate Int[-10,x] gx(x') dx'

In [10]:
RooAbsReal *gx_cdf = gx.createCdf(x);

Plot cdf of gx versus x

In [11]:
RooPlot *frame = x.frame(Title("cdf of Gaussian pdf"));
gx_cdf->plotOn(frame);

Draw plot on canvas

In [12]:
new TCanvas("rf110_normintegration", "rf110_normintegration", 600, 600);
gPad->SetLeftMargin(0.15);
frame->GetYaxis()->SetTitleOffset(1.6);
frame->Draw();

Draw all canvases

In [13]:
%jsroot on
gROOT->GetListOfCanvases()->Draw()