ROOT logo

From $ROOTSYS/tutorials/roofit/rf105_funcbinding.C

//////////////////////////////////////////////////////////////////////////
//
// 'BASIC FUNCTIONALITY' RooFit tutorial macro #105
// 
//  Demonstration of binding ROOT Math functions as RooFit functions
//  and pdfs
//
// 07/2008 - Wouter Verkerke 
// 
/////////////////////////////////////////////////////////////////////////

#ifndef __CINT__
#include "RooGlobalFunc.h"
#else
// Refer to a class implemented in libRooFit to force its loading
// via the autoloader.
class Roo2DKeysPdf;
#endif
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
#include "TMath.h"
#include "TF1.h"
#include "Math/DistFunc.h"
#ifndef __CINT__
#include "RooCFunction1Binding.h" 
#include "RooCFunction3Binding.h"
#endif
#include "RooTFnBinding.h" 

using namespace RooFit;


void rf105_funcbinding()
{

  // B i n d   T M a t h : : E r f   C   f u n c t i o n
  // ---------------------------------------------------

  // Bind one-dimensional TMath::Erf function as RooAbsReal function
  RooRealVar x("x","x",-3,3) ;
  RooAbsReal* erf = bindFunction("erf",TMath::Erf,x) ;

  // Print erf definition
  erf->Print() ;

  // Plot erf on frame 
  RooPlot* frame1 = x.frame(Title("TMath::Erf bound as RooFit function")) ;
  erf->plotOn(frame1) ;


  // B i n d   R O O T : : M a t h : : b e t a _ p d f   C   f u n c t i o n
  // -----------------------------------------------------------------------

  // Bind pdf ROOT::Math::Beta with three variables as RooAbsPdf function
  RooRealVar x2("x2","x2",0,0.999) ;
  RooRealVar a("a","a",5,0,10) ;
  RooRealVar b("b","b",2,0,10) ;
  RooAbsPdf* beta = bindPdf("beta",ROOT::Math::beta_pdf,x2,a,b) ;

  // Perf beta definition
  beta->Print() ;

  // Generate some events and fit
  RooDataSet* data = beta->generate(x2,10000) ;
  beta->fitTo(*data) ;

  // Plot data and pdf on frame
  RooPlot* frame2 = x2.frame(Title("ROOT::Math::Beta bound as RooFit pdf")) ;
  data->plotOn(frame2) ;
  beta->plotOn(frame2) ;



  // B i n d   R O O T   T F 1   a s   R o o F i t   f u n c t i o n
  // ---------------------------------------------------------------

  // Create a ROOT TF1 function
  TF1 *fa1 = new TF1("fa1","sin(x)/x",0,10);
  
  // Create an observable 
  RooRealVar x3("x3","x3",0.01,20) ;

  // Create binding of TF1 object to above observable
  RooAbsReal* rfa1 = bindFunction(fa1,x3) ;

  // Print rfa1 definition
  rfa1->Print() ;

  // Make plot frame in observable, plot TF1 binding function
  RooPlot* frame3 = x3.frame(Title("TF1 bound as RooFit function")) ;
  rfa1->plotOn(frame3) ;



  TCanvas* c = new TCanvas("rf105_funcbinding","rf105_funcbinding",1200,400) ;
  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() ;

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