Logo ROOT   6.08/07
Reference Guide
rf208_convolution.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook -js
4 /// 'ADDITION AND CONVOLUTION' RooFit tutorial macro #208
5 ///
6 /// One-dimensional numeric convolution
7 /// (require ROOT to be compiled with --enable-fftw3)
8 ///
9 /// pdf = landau(t) (x) gauss(t)
10 ///
11 /// \macro_image
12 /// \macro_output
13 /// \macro_code
14 /// \author 07/2008 - Wouter Verkerke
15 
16 
17 #include "RooRealVar.h"
18 #include "RooDataSet.h"
19 #include "RooGaussian.h"
20 #include "RooLandau.h"
21 #include "RooFFTConvPdf.h"
22 #include "RooPlot.h"
23 #include "TCanvas.h"
24 #include "TAxis.h"
25 #include "TH1.h"
26 using namespace RooFit ;
27 
28 
29 
30 void rf208_convolution()
31 {
32  // S e t u p c o m p o n e n t p d f s
33  // ---------------------------------------
34 
35  // Construct observable
36  RooRealVar t("t","t",-10,30) ;
37 
38  // Construct landau(t,ml,sl) ;
39  RooRealVar ml("ml","mean landau",5.,-20,20) ;
40  RooRealVar sl("sl","sigma landau",1,0.1,10) ;
41  RooLandau landau("lx","lx",t,ml,sl) ;
42 
43  // Construct gauss(t,mg,sg)
44  RooRealVar mg("mg","mg",0) ;
45  RooRealVar sg("sg","sg",2,0.1,10) ;
46  RooGaussian gauss("gauss","gauss",t,mg,sg) ;
47 
48 
49  // C o n s t r u c t c o n v o l u t i o n p d f
50  // ---------------------------------------
51 
52  // Set #bins to be used for FFT sampling to 10000
53  t.setBins(10000,"cache") ;
54 
55  // Construct landau (x) gauss
56  RooFFTConvPdf lxg("lxg","landau (X) gauss",t,landau,gauss) ;
57 
58 
59 
60  // S a m p l e , f i t a n d p l o t c o n v o l u t e d p d f
61  // ----------------------------------------------------------------------
62 
63  // Sample 1000 events in x from gxlx
64  RooDataSet* data = lxg.generate(t,10000) ;
65 
66  // Fit gxlx to data
67  lxg.fitTo(*data) ;
68 
69  // Plot data, landau pdf, landau (X) gauss pdf
70  RooPlot* frame = t.frame(Title("landau (x) gauss convolution")) ;
71  data->plotOn(frame) ;
72  lxg.plotOn(frame) ;
73  landau.plotOn(frame,LineStyle(kDashed)) ;
74 
75 
76  // Draw frame on canvas
77  new TCanvas("rf208_convolution","rf208_convolution",600,600) ;
78  gPad->SetLeftMargin(0.15) ; frame->GetYaxis()->SetTitleOffset(1.4) ; frame->Draw() ;
79 
80 }
81 
82 
83 
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title Offset is a correction factor with respect to the "s...
Definition: TAttAxis.cxx:262
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1118
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
Plot dataset on specified frame.
Definition: RooAbsData.cxx:552
Landau Distribution p.d.f.
Definition: RooLandau.h:24
RooCmdArg Title(const char *name)
RooCmdArg LineStyle(Style_t style)
Plain Gaussian p.d.f.
Definition: RooGaussian.h:25
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:37
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:29
A RooPlot is a plot frame and a container for graphics objects within that frame. ...
Definition: RooPlot.h:41
The Canvas class.
Definition: TCanvas.h:41
#define gPad
Definition: TVirtualPad.h:289
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:559