ROOT logo

From $ROOTSYS/tutorials/roofit/rf308_normintegration2d.C

/////////////////////////////////////////////////////////////////////////
//
// 'MULTIDIMENSIONAL MODELS' RooFit tutorial macro #308
// 
// Examples on normalization of p.d.f.s,
// integration of p.d.fs, construction
// of cumulative distribution functions from p.d.f.s
// in two dimensions
//
// 07/2008 - Wouter Verkerke 
//
/////////////////////////////////////////////////////////////////////////

#ifndef __CINT__
#include "RooGlobalFunc.h"
#endif
#include "RooRealVar.h"
#include "RooGaussian.h"
#include "RooConstVar.h"
#include "RooProdPdf.h"
#include "RooAbsReal.h"
#include "RooPlot.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "TH1.h"
using namespace RooFit ;


void rf308_normintegration2d()
{
  // S e t u p   m o d e l 
  // ---------------------

  // Create observables x,y
  RooRealVar x("x","x",-10,10) ;
  RooRealVar y("y","y",-10,10) ;

  // Create p.d.f. gaussx(x,-2,3), gaussy(y,2,2) 
  RooGaussian gx("gx","gx",x,RooConst(-2),RooConst(3)) ;
  RooGaussian gy("gy","gy",y,RooConst(+2),RooConst(2)) ;

  // Create gxy = gx(x)*gy(y)
  RooProdPdf gxy("gxy","gxy",RooArgSet(gx,gy)) ;



  // R e t r i e v e   r a w  &   n o r m a l i z e d   v a l u e s   o f   R o o F i t   p . d . f . s
  // --------------------------------------------------------------------------------------------------

  // Return 'raw' unnormalized value of gx
  cout << "gxy = " << gxy.getVal() << endl ;
  
  // Return value of gxy normalized over x _and_ y in range [-10,10]
  RooArgSet nset_xy(x,y) ;
  cout << "gx_Norm[x,y] = " << gxy.getVal(&nset_xy) << endl ;

  // Create object representing integral over gx
  // which is used to calculate  gx_Norm[x,y] == gx / gx_Int[x,y]
  RooAbsReal* igxy = gxy.createIntegral(RooArgSet(x,y)) ;
  cout << "gx_Int[x,y] = " << igxy->getVal() << endl ;

  // NB: it is also possible to do the following

  // Return value of gxy normalized over x in range [-10,10] (i.e. treating y as parameter)
  RooArgSet nset_x(x) ;
  cout << "gx_Norm[x] = " << gxy.getVal(&nset_x) << endl ;

  // Return value of gxy normalized over y in range [-10,10] (i.e. treating x as parameter)
  RooArgSet nset_y(y) ;
  cout << "gx_Norm[y] = " << gxy.getVal(&nset_y) << endl ;



  // I n t e g r a t e   n o r m a l i z e d   p d f   o v e r   s u b r a n g e
  // ----------------------------------------------------------------------------

  // Define a range named "signal" in x from -5,5
  x.setRange("signal",-5,5) ;
  y.setRange("signal",-3,3) ;
  
  // Create an integral of gxy_Norm[x,y] over x and y in range "signal"
  // This is the fraction of of p.d.f. gxy_Norm[x,y] which is in the
  // range named "signal"
  RooAbsReal* igxy_sig = gxy.createIntegral(RooArgSet(x,y),NormSet(RooArgSet(x,y)),Range("signal")) ;
  cout << "gx_Int[x,y|signal]_Norm[x,y] = " << igxy_sig->getVal() << endl ;




  // C o n s t r u c t   c u m u l a t i v e   d i s t r i b u t i o n   f u n c t i o n   f r o m   p d f
  // -----------------------------------------------------------------------------------------------------

  // Create the cumulative distribution function of gx
  // i.e. calculate Int[-10,x] gx(x') dx'
  RooAbsReal* gxy_cdf = gxy.createCdf(RooArgSet(x,y)) ;
  
  // Plot cdf of gx versus x
  TH1* hh_cdf = gxy_cdf->createHistogram("hh_cdf",x,Binning(40),YVar(y,Binning(40))) ;
  hh_cdf->SetLineColor(kBlue) ;

  new TCanvas("rf308_normintegration2d","rf308_normintegration2d",600,600) ;
  gPad->SetLeftMargin(0.15) ; hh_cdf->GetZaxis()->SetTitleOffset(1.8) ; 
  hh_cdf->Draw("surf") ;

}
 rf308_normintegration2d.C:1
 rf308_normintegration2d.C:2
 rf308_normintegration2d.C:3
 rf308_normintegration2d.C:4
 rf308_normintegration2d.C:5
 rf308_normintegration2d.C:6
 rf308_normintegration2d.C:7
 rf308_normintegration2d.C:8
 rf308_normintegration2d.C:9
 rf308_normintegration2d.C:10
 rf308_normintegration2d.C:11
 rf308_normintegration2d.C:12
 rf308_normintegration2d.C:13
 rf308_normintegration2d.C:14
 rf308_normintegration2d.C:15
 rf308_normintegration2d.C:16
 rf308_normintegration2d.C:17
 rf308_normintegration2d.C:18
 rf308_normintegration2d.C:19
 rf308_normintegration2d.C:20
 rf308_normintegration2d.C:21
 rf308_normintegration2d.C:22
 rf308_normintegration2d.C:23
 rf308_normintegration2d.C:24
 rf308_normintegration2d.C:25
 rf308_normintegration2d.C:26
 rf308_normintegration2d.C:27
 rf308_normintegration2d.C:28
 rf308_normintegration2d.C:29
 rf308_normintegration2d.C:30
 rf308_normintegration2d.C:31
 rf308_normintegration2d.C:32
 rf308_normintegration2d.C:33
 rf308_normintegration2d.C:34
 rf308_normintegration2d.C:35
 rf308_normintegration2d.C:36
 rf308_normintegration2d.C:37
 rf308_normintegration2d.C:38
 rf308_normintegration2d.C:39
 rf308_normintegration2d.C:40
 rf308_normintegration2d.C:41
 rf308_normintegration2d.C:42
 rf308_normintegration2d.C:43
 rf308_normintegration2d.C:44
 rf308_normintegration2d.C:45
 rf308_normintegration2d.C:46
 rf308_normintegration2d.C:47
 rf308_normintegration2d.C:48
 rf308_normintegration2d.C:49
 rf308_normintegration2d.C:50
 rf308_normintegration2d.C:51
 rf308_normintegration2d.C:52
 rf308_normintegration2d.C:53
 rf308_normintegration2d.C:54
 rf308_normintegration2d.C:55
 rf308_normintegration2d.C:56
 rf308_normintegration2d.C:57
 rf308_normintegration2d.C:58
 rf308_normintegration2d.C:59
 rf308_normintegration2d.C:60
 rf308_normintegration2d.C:61
 rf308_normintegration2d.C:62
 rf308_normintegration2d.C:63
 rf308_normintegration2d.C:64
 rf308_normintegration2d.C:65
 rf308_normintegration2d.C:66
 rf308_normintegration2d.C:67
 rf308_normintegration2d.C:68
 rf308_normintegration2d.C:69
 rf308_normintegration2d.C:70
 rf308_normintegration2d.C:71
 rf308_normintegration2d.C:72
 rf308_normintegration2d.C:73
 rf308_normintegration2d.C:74
 rf308_normintegration2d.C:75
 rf308_normintegration2d.C:76
 rf308_normintegration2d.C:77
 rf308_normintegration2d.C:78
 rf308_normintegration2d.C:79
 rf308_normintegration2d.C:80
 rf308_normintegration2d.C:81
 rf308_normintegration2d.C:82
 rf308_normintegration2d.C:83
 rf308_normintegration2d.C:84
 rf308_normintegration2d.C:85
 rf308_normintegration2d.C:86
 rf308_normintegration2d.C:87
 rf308_normintegration2d.C:88
 rf308_normintegration2d.C:89
 rf308_normintegration2d.C:90
 rf308_normintegration2d.C:91
 rf308_normintegration2d.C:92
 rf308_normintegration2d.C:93
 rf308_normintegration2d.C:94
 rf308_normintegration2d.C:95
 rf308_normintegration2d.C:96
 rf308_normintegration2d.C:97
 rf308_normintegration2d.C:98
 rf308_normintegration2d.C:99
 rf308_normintegration2d.C:100
 rf308_normintegration2d.C:101
 rf308_normintegration2d.C:102
 rf308_normintegration2d.C:103
 rf308_normintegration2d.C:104
 rf308_normintegration2d.C:105
 rf308_normintegration2d.C:106