Logo ROOT  
Reference Guide
rf105_funcbinding.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook -js
4 ///
5 /// Basic functionality: binding ROOT math functions as RooFit functions and pdfs
6 ///
7 /// \macro_image
8 /// \macro_output
9 /// \macro_code
10 ///
11 /// \date 07/2008
12 /// \author Wouter Verkerke
13 
14 #include "RooRealVar.h"
15 #include "RooDataSet.h"
16 #include "RooGaussian.h"
17 #include "TCanvas.h"
18 #include "TAxis.h"
19 #include "RooPlot.h"
20 #include "TMath.h"
21 #include "TF1.h"
22 #include "Math/DistFunc.h"
23 #include "RooTFnBinding.h"
24 
25 using namespace RooFit;
26 
27 void rf105_funcbinding()
28 {
29 
30  // B i n d T M a t h : : E r f C f u n c t i o n
31  // ---------------------------------------------------
32 
33  // Bind one-dimensional TMath::Erf function as RooAbsReal function
34  RooRealVar x("x", "x", -3, 3);
35  RooAbsReal *errorFunc = bindFunction("erf", TMath::Erf, x);
36 
37  // Print erf definition
38  errorFunc->Print();
39 
40  // Plot erf on frame
41  RooPlot *frame1 = x.frame(Title("TMath::Erf bound as RooFit function"));
42  errorFunc->plotOn(frame1);
43 
44  // 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
45  // -----------------------------------------------------------------------
46 
47  // Bind pdf ROOT::Math::Beta with three variables as RooAbsPdf function
48  RooRealVar x2("x2", "x2", 0, 0.999);
49  RooRealVar a("a", "a", 5, 0, 10);
50  RooRealVar b("b", "b", 2, 0, 10);
52 
53  // Perf beta definition
54  beta->Print();
55 
56  // Generate some events and fit
57  RooDataSet *data = beta->generate(x2, 10000);
58  beta->fitTo(*data);
59 
60  // Plot data and pdf on frame
61  RooPlot *frame2 = x2.frame(Title("ROOT::Math::Beta bound as RooFit pdf"));
62  data->plotOn(frame2);
63  beta->plotOn(frame2);
64 
65  // 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
66  // ---------------------------------------------------------------
67 
68  // Create a ROOT TF1 function
69  TF1 *fa1 = new TF1("fa1", "sin(x)/x", 0, 10);
70 
71  // Create an observable
72  RooRealVar x3("x3", "x3", 0.01, 20);
73 
74  // Create binding of TF1 object to above observable
75  RooAbsReal *rfa1 = bindFunction(fa1, x3);
76 
77  // Print rfa1 definition
78  rfa1->Print();
79 
80  // Make plot frame in observable, plot TF1 binding function
81  RooPlot *frame3 = x3.frame(Title("TF1 bound as RooFit function"));
82  rfa1->plotOn(frame3);
83 
84  TCanvas *c = new TCanvas("rf105_funcbinding", "rf105_funcbinding", 1200, 400);
85  c->Divide(3);
86  c->cd(1);
87  gPad->SetLeftMargin(0.15);
88  frame1->GetYaxis()->SetTitleOffset(1.6);
89  frame1->Draw();
90  c->cd(2);
91  gPad->SetLeftMargin(0.15);
92  frame2->GetYaxis()->SetTitleOffset(1.6);
93  frame2->Draw();
94  c->cd(3);
95  gPad->SetLeftMargin(0.15);
96  frame3->GetYaxis()->SetTitleOffset(1.6);
97  frame3->Draw();
98 }
c
#define c(i)
Definition: RSha256.hxx:119
RooPlot::Draw
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:691
ROOT::Math::beta_pdf
double beta_pdf(double x, double a, double b)
Probability density function of the beta distribution.
Definition: PdfFuncMathCore.h:82
RooAbsArg::Print
virtual void Print(Option_t *options=0) const
Print the object to the defaultPrintStream().
Definition: RooAbsArg.h:320
RooGaussian.h
x
Double_t x[n]
Definition: legend1.C:17
TMath::Erf
Double_t Erf(Double_t x)
Computation of the error function erf(x).
Definition: TMath.cxx:184
RooAbsReal
Definition: RooAbsReal.h:61
TCanvas.h
ROOT::Math::beta
double beta(double x, double y)
Calculates the beta function.
Definition: SpecFuncMathCore.cxx:111
RooDataSet.h
b
#define b(i)
Definition: RSha256.hxx:118
DistFunc.h
RooFit::bindFunction
RooAbsReal * bindFunction(const char *name, CFUNCD1D func, RooAbsReal &x)
Definition: RooCFunction1Binding.cxx:59
x3
static const double x3[11]
Definition: RooGaussKronrodIntegrator1D.cxx:392
RooTFnBinding.h
RooFit
Definition: RooCFunction1Binding.h:29
RooAbsData::plotOn
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Definition: RooAbsData.cxx:547
a
auto * a
Definition: textangle.C:12
RooPlot.h
RooPlot::GetYaxis
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1256
RooPlot
Definition: RooPlot.h:44
RooRealVar.h
RooAbsReal::plotOn
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg(), const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg(), const RooCmdArg &arg9=RooCmdArg(), const RooCmdArg &arg10=RooCmdArg()) const
Plot (project) PDF on specified frame.
Definition: RooAbsReal.cxx:1714
TF1.h
TCanvas
Definition: TCanvas.h:23
TAxis.h
x2
static const double x2[5]
Definition: RooGaussKronrodIntegrator1D.cxx:364
gPad
#define gPad
Definition: TVirtualPad.h:287
RooDataSet
Definition: RooDataSet.h:33
RooAbsPdf
Definition: RooAbsPdf.h:40
TAttAxis::SetTitleOffset
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition: TAttAxis.cxx:293
TF1
1-Dim function class
Definition: TF1.h:212
RooRealVar
Definition: RooRealVar.h:35
RooFit::bindPdf
RooAbsPdf * bindPdf(const char *name, CFUNCD1D func, RooAbsReal &x)
Definition: RooCFunction1Binding.cxx:67
RooFit::Title
RooCmdArg Title(const char *name)
Definition: RooGlobalFunc.cxx:173
TMath.h