Logo ROOT  
Reference Guide
rf211_paramconv.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook
4 ///
5 /// Addition and convolution: working with a p.d.f. with a convolution operator in terms of a parameter
6 ///
7 /// This tutorial requires FFT3 to be enabled.
8 ///
9 /// \macro_image
10 /// \macro_output
11 /// \macro_code
12 ///
13 /// \date 04/2009
14 /// \author Wouter Verkerke
15 
16 #include "RooRealVar.h"
17 #include "RooDataHist.h"
18 #include "RooGaussian.h"
19 #include "RooGenericPdf.h"
20 #include "RooFormulaVar.h"
21 #include "RooFFTConvPdf.h"
22 #include "RooPlot.h"
23 #include "TCanvas.h"
24 #include "TAxis.h"
25 #include "TH2.h"
26 using namespace RooFit;
27 
28 void rf211_paramconv()
29 {
30  // S e t u p c o m p o n e n t p d f s
31  // ---------------------------------------
32 
33  // Gaussian g(x ; mean,sigma)
34  RooRealVar x("x", "x", -10, 10);
35  RooRealVar mean("mean", "mean", -3, 3);
36  RooRealVar sigma("sigma", "sigma", 0.5, 0.1, 10);
37  RooGaussian modelx("gx", "gx", x, mean, sigma);
38 
39  // Block function in mean
40  RooRealVar a("a", "a", 2, 1, 10);
41  RooGenericPdf model_mean("model_mean", "abs(mean)<a", RooArgList(mean, a));
42 
43  // Convolution in mean parameter model = g(x,mean,sigma) (x) block(mean)
44  x.setBins(1000, "cache");
45  mean.setBins(50, "cache");
46  RooFFTConvPdf model("model", "model", mean, modelx, model_mean);
47 
48  // Configure convolution to construct a 2-D cache in (x,mean)
49  // rather than a 1-d cache in mean that needs to be recalculated
50  // for each value of x
51  model.setCacheObservables(x);
52  model.setBufferFraction(1.0);
53 
54  // Integrate model over mean projModel = Int model dmean
55  RooAbsPdf *projModel = model.createProjection(mean);
56 
57  // Generate 1000 toy events
58  RooDataHist *d = projModel->generateBinned(x, 1000);
59 
60  // Fit p.d.f. to toy data
61  projModel->fitTo(*d, Verbose());
62 
63  // Plot data and fitted p.d.f.
64  RooPlot *frame = x.frame(Bins(25));
65  d->plotOn(frame);
66  projModel->plotOn(frame);
67 
68  // Make 2d histogram of model(x;mean)
69  TH1 *hh = model.createHistogram("hh", x, Binning(50), YVar(mean, Binning(50)), ConditionalObservables(mean));
70  hh->SetTitle("histogram of model(x|mean)");
71  hh->SetLineColor(kBlue);
72 
73  // Draw frame on canvas
74  TCanvas *c = new TCanvas("rf211_paramconv", "rf211_paramconv", 800, 400);
75  c->Divide(2);
76  c->cd(1);
77  gPad->SetLeftMargin(0.15);
78  frame->GetYaxis()->SetTitleOffset(1.4);
79  frame->Draw();
80  c->cd(2);
81  gPad->SetLeftMargin(0.20);
82  hh->GetZaxis()->SetTitleOffset(2.5);
83  hh->Draw("surf");
84 }
RooFormulaVar.h
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
RooFit::Bins
RooCmdArg Bins(Int_t nbin)
Definition: RooGlobalFunc.cxx:174
RooArgList
Definition: RooArgList.h:21
RooGaussian.h
x
Double_t x[n]
Definition: legend1.C:17
RooGaussian
Definition: RooGaussian.h:25
TCanvas.h
RooFit::YVar
RooCmdArg YVar(const RooAbsRealLValue &var, const RooCmdArg &arg=RooCmdArg::none())
Definition: RooGlobalFunc.cxx:240
RooFit::Binning
RooCmdArg Binning(const RooAbsBinning &binning)
Definition: RooGlobalFunc.cxx:82
RooGenericPdf.h
TH1::GetZaxis
TAxis * GetZaxis()
Definition: TH1.h:319
RooAbsPdf::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 RooCmdArg &arg9=RooCmdArg::none(), const RooCmdArg &arg10=RooCmdArg::none()) const
Helper calling plotOn(RooPlot*, RooLinkedList&) const.
Definition: RooAbsPdf.h:118
TH1::SetTitle
virtual void SetTitle(const char *title)
See GetStatOverflows for more information.
Definition: TH1.cxx:6344
RooDataHist
Definition: RooDataHist.h:39
RooFit
Definition: RooCFunction1Binding.h:29
a
auto * a
Definition: textangle.C:12
RooFit::Verbose
RooCmdArg Verbose(Bool_t flag=kTRUE)
Definition: RooGlobalFunc.cxx:186
RooDataHist.h
RooPlot.h
RooPlot::GetYaxis
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1256
RooPlot
Definition: RooPlot.h:44
RooRealVar.h
RooFit::ConditionalObservables
RooCmdArg ConditionalObservables(const RooArgSet &set)
Definition: RooGlobalFunc.cxx:196
TH2.h
sigma
const Double_t sigma
Definition: h1analysisProxy.h:11
RooFFTConvPdf
PDF for the numerical (FFT) convolution of two PDFs.
Definition: RooFFTConvPdf.h:25
TCanvas
Definition: TCanvas.h:23
TAxis.h
TH1
Definition: TH1.h:57
RooGenericPdf
Definition: RooGenericPdf.h:25
kBlue
@ kBlue
Definition: Rtypes.h:66
RooFFTConvPdf.h
d
#define d(i)
Definition: RSha256.hxx:120
gPad
#define gPad
Definition: TVirtualPad.h:287
RooAbsPdf::generateBinned
virtual RooDataHist * generateBinned(const RooArgSet &whatVars, Double_t nEvents, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none()) const
As RooAbsPdf::generateBinned(const RooArgSet&, const RooCmdArg&,const RooCmdArg&, const RooCmdArg&,...
Definition: RooAbsPdf.h:104
make_cnn_model.model
model
Definition: make_cnn_model.py:6
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
RooRealVar
Definition: RooRealVar.h:35
TH1::Draw
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2997
RooAbsPdf::fitTo
virtual RooFitResult * fitTo(RooAbsData &data, 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())
Fit PDF to given dataset.
Definition: RooAbsPdf.cxx:1261