ROOT logo

From $ROOTSYS/tutorials/roofit/rf109_chi2residpull.C

//////////////////////////////////////////////////////////////////////////
//
// 'BASIC FUNCTIONALITY' RooFit tutorial macro #109
// 
// Calculating chi^2 from histograms and curves in RooPlots, 
// making histogram of residual and pull distributions
//
//
//
// 07/2008 - Wouter Verkerke 
// 
/////////////////////////////////////////////////////////////////////////

#ifndef __CINT__
#include "RooGlobalFunc.h"
#endif
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooConstVar.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
#include "RooHist.h"
using namespace RooFit ;


void rf109_chi2residpull()
{

  // S e t u p   m o d e l 
  // ---------------------

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

  // Create Gaussian
  RooRealVar sigma("sigma","sigma",3,0.1,10) ;
  RooRealVar mean("mean","mean",0,-10,10) ;
  RooGaussian gauss("gauss","gauss",x,RooConst(0),sigma) ;

  // Generate a sample of 1000 events with sigma=3
  RooDataSet* data = gauss.generate(x,10000) ;

  // Change sigma to 3.15
  sigma=3.15 ;


  // P l o t   d a t a   a n d   s l i g h t l y   d i s t o r t e d   m o d e l
  // ---------------------------------------------------------------------------

  // Overlay projection of gauss with sigma=3.15 on data with sigma=3.0
  RooPlot* frame1 = x.frame(Title("Data with distorted Gaussian pdf"),Bins(40)) ;
  data->plotOn(frame1,DataError(RooAbsData::SumW2)) ;
  gauss.plotOn(frame1) ;


  // C a l c u l a t e   c h i ^ 2 
  // ------------------------------

  // Show the chi^2 of the curve w.r.t. the histogram
  // If multiple curves or datasets live in the frame you can specify
  // the name of the relevant curve and/or dataset in chiSquare()
  cout << "chi^2 = " << frame1->chiSquare() << endl ;


  // S h o w   r e s i d u a l   a n d   p u l l   d i s t s
  // -------------------------------------------------------

  // Construct a histogram with the residuals of the data w.r.t. the curve
  RooHist* hresid = frame1->residHist() ;

  // Construct a histogram with the pulls of the data w.r.t the curve
  RooHist* hpull = frame1->pullHist() ;

  // Create a new frame to draw the residual distribution and add the distribution to the frame
  RooPlot* frame2 = x.frame(Title("Residual Distribution")) ;
  frame2->addPlotable(hresid,"P") ;

  // Create a new frame to draw the pull distribution and add the distribution to the frame
  RooPlot* frame3 = x.frame(Title("Pull Distribution")) ;
  frame3->addPlotable(hpull,"P") ;



  TCanvas* c = new TCanvas("rf109_chi2residpull","rf109_chi2residpull",900,300) ;
  c->Divide(3) ;
  c->cd(1) ; gPad->SetLeftMargin(0.15) ; frame1->GetYaxis()->SetTitleOffset(1.6) ; frame1->Draw() ;
  c->cd(2) ; gPad->SetLeftMargin(0.15) ; frame2->GetYaxis()->SetTitleOffset(1.6) ; frame2->Draw() ;
  c->cd(3) ; gPad->SetLeftMargin(0.15) ; frame3->GetYaxis()->SetTitleOffset(1.6) ; frame3->Draw() ;
  
}
 rf109_chi2residpull.C:1
 rf109_chi2residpull.C:2
 rf109_chi2residpull.C:3
 rf109_chi2residpull.C:4
 rf109_chi2residpull.C:5
 rf109_chi2residpull.C:6
 rf109_chi2residpull.C:7
 rf109_chi2residpull.C:8
 rf109_chi2residpull.C:9
 rf109_chi2residpull.C:10
 rf109_chi2residpull.C:11
 rf109_chi2residpull.C:12
 rf109_chi2residpull.C:13
 rf109_chi2residpull.C:14
 rf109_chi2residpull.C:15
 rf109_chi2residpull.C:16
 rf109_chi2residpull.C:17
 rf109_chi2residpull.C:18
 rf109_chi2residpull.C:19
 rf109_chi2residpull.C:20
 rf109_chi2residpull.C:21
 rf109_chi2residpull.C:22
 rf109_chi2residpull.C:23
 rf109_chi2residpull.C:24
 rf109_chi2residpull.C:25
 rf109_chi2residpull.C:26
 rf109_chi2residpull.C:27
 rf109_chi2residpull.C:28
 rf109_chi2residpull.C:29
 rf109_chi2residpull.C:30
 rf109_chi2residpull.C:31
 rf109_chi2residpull.C:32
 rf109_chi2residpull.C:33
 rf109_chi2residpull.C:34
 rf109_chi2residpull.C:35
 rf109_chi2residpull.C:36
 rf109_chi2residpull.C:37
 rf109_chi2residpull.C:38
 rf109_chi2residpull.C:39
 rf109_chi2residpull.C:40
 rf109_chi2residpull.C:41
 rf109_chi2residpull.C:42
 rf109_chi2residpull.C:43
 rf109_chi2residpull.C:44
 rf109_chi2residpull.C:45
 rf109_chi2residpull.C:46
 rf109_chi2residpull.C:47
 rf109_chi2residpull.C:48
 rf109_chi2residpull.C:49
 rf109_chi2residpull.C:50
 rf109_chi2residpull.C:51
 rf109_chi2residpull.C:52
 rf109_chi2residpull.C:53
 rf109_chi2residpull.C:54
 rf109_chi2residpull.C:55
 rf109_chi2residpull.C:56
 rf109_chi2residpull.C:57
 rf109_chi2residpull.C:58
 rf109_chi2residpull.C:59
 rf109_chi2residpull.C:60
 rf109_chi2residpull.C:61
 rf109_chi2residpull.C:62
 rf109_chi2residpull.C:63
 rf109_chi2residpull.C:64
 rf109_chi2residpull.C:65
 rf109_chi2residpull.C:66
 rf109_chi2residpull.C:67
 rf109_chi2residpull.C:68
 rf109_chi2residpull.C:69
 rf109_chi2residpull.C:70
 rf109_chi2residpull.C:71
 rf109_chi2residpull.C:72
 rf109_chi2residpull.C:73
 rf109_chi2residpull.C:74
 rf109_chi2residpull.C:75
 rf109_chi2residpull.C:76
 rf109_chi2residpull.C:77
 rf109_chi2residpull.C:78
 rf109_chi2residpull.C:79
 rf109_chi2residpull.C:80
 rf109_chi2residpull.C:81
 rf109_chi2residpull.C:82
 rf109_chi2residpull.C:83
 rf109_chi2residpull.C:84
 rf109_chi2residpull.C:85
 rf109_chi2residpull.C:86
 rf109_chi2residpull.C:87
 rf109_chi2residpull.C:88
 rf109_chi2residpull.C:89
 rf109_chi2residpull.C:90
 rf109_chi2residpull.C:91
 rf109_chi2residpull.C:92
 rf109_chi2residpull.C:93