From $ROOTSYS/tutorials/roofit/rf309_ndimplot.C

//////////////////////////////////////////////////////////////////////////
//
// 'MULTIDIMENSIONAL MODELS' RooFit tutorial macro #308
// 
// Making 2/3 dimensional plots of p.d.f.s and datasets
//
//
//
// 07/2008 - Wouter Verkerke 
// 
/////////////////////////////////////////////////////////////////////////

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


void rf309_ndimplot()
{

  // C r e a t e   2 D   m o d e l   a n d   d a t a s e t
  // -----------------------------------------------------

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

  // Create parameters
  RooRealVar a0("a0","a0",-3.5,-5,5) ;
  RooRealVar a1("a1","a1",-1.5,-1,1) ;
  RooRealVar sigma("sigma","width of gaussian",1.5) ;

  // Create interpreted function f(y) = a0 - a1*sqrt(10*abs(y))
  RooFormulaVar fy("fy","a0-a1*sqrt(10*abs(y))",RooArgSet(y,a0,a1)) ;

  // Create gauss(x,f(y),s)
  RooGaussian model("model","Gaussian with shifting mean",x,fy,sigma) ;

  // Sample dataset from gauss(x,y)
  RooDataSet* data = model.generate(RooArgSet(x,y),10000) ;


  // M a k e   2 D   p l o t s   o f   d a t a   a n d   m o d e l
  // -------------------------------------------------------------

  // Create and fill ROOT 2D histogram (20x20 bins) with contents of dataset
  //TH2D* hh_data = data->createHistogram("hh_data",x,Binning(20),YVar(y,Binning(20))) ;
  TH1* hh_data = data->createHistogram("x,y",20,20) ;

  // Create and fill ROOT 2D histogram (50x50 bins) with sampling of pdf
  //TH2D* hh_pdf = model.createHistogram("hh_model",x,Binning(50),YVar(y,Binning(50))) ;
  TH1* hh_pdf = model.createHistogram("x,y",50,50) ;
  hh_pdf->SetLineColor(kBlue) ;


  // C r e a t e   3 D   m o d e l   a n d   d a t a s e t
  // -----------------------------------------------------

  // Create observables
  RooRealVar z("z","z",-5,5) ;

  RooGaussian gz("gz","gz",z,RooConst(0),RooConst(2)) ;
  RooProdPdf model3("model3","model3",RooArgSet(model,gz)) ;

  RooDataSet* data3 = model3.generate(RooArgSet(x,y,z),10000) ;

  
  // M a k e   3 D   p l o t s   o f   d a t a   a n d   m o d e l
  // -------------------------------------------------------------

  // Create and fill ROOT 2D histogram (8x8x8 bins) with contents of dataset
  TH1* hh_data3 = data3->createHistogram("hh_data3",x,Binning(8),YVar(y,Binning(8)),ZVar(z,Binning(8))) ;

  // Create and fill ROOT 2D histogram (20x20x20 bins) with sampling of pdf
  TH1* hh_pdf3 = model3.createHistogram("hh_model3",x,Binning(20),YVar(y,Binning(20)),ZVar(z,Binning(20))) ;
  hh_pdf3->SetFillColor(kBlue) ;



  TCanvas* c1 = new TCanvas("rf309_2dimplot","rf309_2dimplot",800,800) ;
  c1->Divide(2,2) ;
  c1->cd(1) ; gPad->SetLeftMargin(0.15) ; hh_data->GetZaxis()->SetTitleOffset(1.4) ; hh_data->Draw("lego") ; 
  c1->cd(2) ; gPad->SetLeftMargin(0.20) ; hh_pdf->GetZaxis()->SetTitleOffset(2.5) ; hh_pdf->Draw("surf") ;
  c1->cd(3) ; gPad->SetLeftMargin(0.15) ; hh_data->GetZaxis()->SetTitleOffset(1.4) ; hh_data->Draw("box") ; 
  c1->cd(4) ; gPad->SetLeftMargin(0.15) ; hh_pdf->GetZaxis()->SetTitleOffset(2.5) ; hh_pdf->Draw("cont3") ;
  
  TCanvas* c2 = new TCanvas("rf309_3dimplot","rf309_3dimplot",800,400) ;
  c2->Divide(2) ;
  c2->cd(1) ; gPad->SetLeftMargin(0.15) ; hh_data3->GetZaxis()->SetTitleOffset(1.4) ; hh_data3->Draw("lego") ;
  c2->cd(2) ; gPad->SetLeftMargin(0.15) ; hh_pdf3->GetZaxis()->SetTitleOffset(1.4) ; hh_pdf3->Draw("iso") ;
  
}
 rf309_ndimplot.C:1
 rf309_ndimplot.C:2
 rf309_ndimplot.C:3
 rf309_ndimplot.C:4
 rf309_ndimplot.C:5
 rf309_ndimplot.C:6
 rf309_ndimplot.C:7
 rf309_ndimplot.C:8
 rf309_ndimplot.C:9
 rf309_ndimplot.C:10
 rf309_ndimplot.C:11
 rf309_ndimplot.C:12
 rf309_ndimplot.C:13
 rf309_ndimplot.C:14
 rf309_ndimplot.C:15
 rf309_ndimplot.C:16
 rf309_ndimplot.C:17
 rf309_ndimplot.C:18
 rf309_ndimplot.C:19
 rf309_ndimplot.C:20
 rf309_ndimplot.C:21
 rf309_ndimplot.C:22
 rf309_ndimplot.C:23
 rf309_ndimplot.C:24
 rf309_ndimplot.C:25
 rf309_ndimplot.C:26
 rf309_ndimplot.C:27
 rf309_ndimplot.C:28
 rf309_ndimplot.C:29
 rf309_ndimplot.C:30
 rf309_ndimplot.C:31
 rf309_ndimplot.C:32
 rf309_ndimplot.C:33
 rf309_ndimplot.C:34
 rf309_ndimplot.C:35
 rf309_ndimplot.C:36
 rf309_ndimplot.C:37
 rf309_ndimplot.C:38
 rf309_ndimplot.C:39
 rf309_ndimplot.C:40
 rf309_ndimplot.C:41
 rf309_ndimplot.C:42
 rf309_ndimplot.C:43
 rf309_ndimplot.C:44
 rf309_ndimplot.C:45
 rf309_ndimplot.C:46
 rf309_ndimplot.C:47
 rf309_ndimplot.C:48
 rf309_ndimplot.C:49
 rf309_ndimplot.C:50
 rf309_ndimplot.C:51
 rf309_ndimplot.C:52
 rf309_ndimplot.C:53
 rf309_ndimplot.C:54
 rf309_ndimplot.C:55
 rf309_ndimplot.C:56
 rf309_ndimplot.C:57
 rf309_ndimplot.C:58
 rf309_ndimplot.C:59
 rf309_ndimplot.C:60
 rf309_ndimplot.C:61
 rf309_ndimplot.C:62
 rf309_ndimplot.C:63
 rf309_ndimplot.C:64
 rf309_ndimplot.C:65
 rf309_ndimplot.C:66
 rf309_ndimplot.C:67
 rf309_ndimplot.C:68
 rf309_ndimplot.C:69
 rf309_ndimplot.C:70
 rf309_ndimplot.C:71
 rf309_ndimplot.C:72
 rf309_ndimplot.C:73
 rf309_ndimplot.C:74
 rf309_ndimplot.C:75
 rf309_ndimplot.C:76
 rf309_ndimplot.C:77
 rf309_ndimplot.C:78
 rf309_ndimplot.C:79
 rf309_ndimplot.C:80
 rf309_ndimplot.C:81
 rf309_ndimplot.C:82
 rf309_ndimplot.C:83
 rf309_ndimplot.C:84
 rf309_ndimplot.C:85
 rf309_ndimplot.C:86
 rf309_ndimplot.C:87
 rf309_ndimplot.C:88
 rf309_ndimplot.C:89
 rf309_ndimplot.C:90
 rf309_ndimplot.C:91
 rf309_ndimplot.C:92
 rf309_ndimplot.C:93
 rf309_ndimplot.C:94
 rf309_ndimplot.C:95
 rf309_ndimplot.C:96
 rf309_ndimplot.C:97
 rf309_ndimplot.C:98
 rf309_ndimplot.C:99
 rf309_ndimplot.C:100
 rf309_ndimplot.C:101
 rf309_ndimplot.C:102
 rf309_ndimplot.C:103