Logo ROOT  
Reference Guide
rf210_angularconv.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook -js
4 /// Addition and convolution: convolution in cyclical angular observables theta
5 ///
6 /// and construction of pdf in terms of transformed angular coordinates, e.g. cos(theta),
7 /// where the convolution is performed in theta rather than cos(theta)
8 ///
9 /// ```
10 /// pdf(theta) = T(theta) (x) gauss(theta)
11 /// pdf(cosTheta) = T(acos(cosTheta)) (x) gauss(acos(cosTheta))
12 /// ```
13 ///
14 /// This tutorial requires FFT3 to be enabled.
15 ///
16 /// \macro_image
17 /// \macro_output
18 /// \macro_code
19 ///
20 /// \date April 2009
21 /// \author Wouter Verkerke
22 
23 #include "RooRealVar.h"
24 #include "RooDataSet.h"
25 #include "RooGaussian.h"
26 #include "RooGenericPdf.h"
27 #include "RooFormulaVar.h"
28 #include "RooFFTConvPdf.h"
29 #include "RooPlot.h"
30 #include "TCanvas.h"
31 #include "TAxis.h"
32 #include "TH1.h"
33 using namespace RooFit;
34 
35 void rf210_angularconv()
36 {
37  // S e t u p c o m p o n e n t p d f s
38  // ---------------------------------------
39 
40  // Define angle psi
41  RooRealVar psi("psi", "psi", 0, 3.14159268);
42 
43  // Define physics pdf T(psi)
44  RooGenericPdf Tpsi("Tpsi", "1+sin(2*@0)", psi);
45 
46  // Define resolution R(psi)
47  RooRealVar gbias("gbias", "gbias", 0.2, 0., 1);
48  RooRealVar greso("greso", "greso", 0.3, 0.1, 1.0);
49  RooGaussian Rpsi("Rpsi", "Rpsi", psi, gbias, greso);
50 
51  // Define cos(psi) and function psif that calculates psi from cos(psi)
52  RooRealVar cpsi("cpsi", "cos(psi)", -1, 1);
53  RooFormulaVar psif("psif", "acos(cpsi)", cpsi);
54 
55  // Define physics pdf also as function of cos(psi): T(psif(cpsi)) = T(cpsi) ;
56  RooGenericPdf Tcpsi("T", "1+sin(2*@0)", psif);
57 
58  // C o n s t r u c t c o n v o l u t i o n p d f i n p s i
59  // --------------------------------------------------------------
60 
61  // Define convoluted pdf as function of psi: M=[T(x)R](psi) = M(psi)
62  RooFFTConvPdf Mpsi("Mf", "Mf", psi, Tpsi, Rpsi);
63 
64  // Set the buffer fraction to zero to obtain a true cyclical convolution
65  Mpsi.setBufferFraction(0);
66 
67  // 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 ( p s i )
68  // --------------------------------------------------------------------------------
69 
70  // Generate some events in observable psi
71  RooDataSet *data_psi = Mpsi.generate(psi, 10000);
72 
73  // Fit convoluted model as function of angle psi
74  Mpsi.fitTo(*data_psi);
75 
76  // Plot cos(psi) frame with Mf(cpsi)
77  RooPlot *frame1 = psi.frame(Title("Cyclical convolution in angle psi"));
78  data_psi->plotOn(frame1);
79  Mpsi.plotOn(frame1);
80 
81  // Overlay comparison to unsmeared physics pdf T(psi)
82  Tpsi.plotOn(frame1, LineColor(kRed));
83 
84  // C o n s t r u c t c o n v o l u t i o n p d f i n c o s ( p s i )
85  // --------------------------------------------------------------------------
86 
87  // Define convoluted pdf as function of cos(psi): M=[T(x)R](psif(cpsi)) = M(cpsi)
88  //
89  // Need to give both observable psi here (for definition of convolution)
90  // and function psif here (for definition of observables, ultimately in cpsi)
91  RooFFTConvPdf Mcpsi("Mf", "Mf", psif, psi, Tpsi, Rpsi);
92 
93  // Set the buffer fraction to zero to obtain a true cyclical convolution
94  Mcpsi.setBufferFraction(0);
95 
96  // 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 ( c o s p s i )
97  // --------------------------------------------------------------------------------
98 
99  // Generate some events
100  RooDataSet *data_cpsi = Mcpsi.generate(cpsi, 10000);
101 
102  // set psi constant to exclude to be a parameter of the fit
103  psi.setConstant(true);
104 
105  // Fit convoluted model as function of cos(psi)
106  Mcpsi.fitTo(*data_cpsi);
107 
108  // Plot cos(psi) frame with Mf(cpsi)
109  RooPlot *frame2 = cpsi.frame(Title("Same convolution in psi, expressed in cos(psi)"));
110  data_cpsi->plotOn(frame2);
111  Mcpsi.plotOn(frame2);
112 
113  // Overlay comparison to unsmeared physics pdf Tf(cpsi)
114  Tcpsi.plotOn(frame2, LineColor(kRed));
115 
116  // Draw frame on canvas
117  TCanvas *c = new TCanvas("rf210_angularconv", "rf210_angularconv", 800, 400);
118  c->Divide(2);
119  c->cd(1);
120  gPad->SetLeftMargin(0.15);
121  frame1->GetYaxis()->SetTitleOffset(1.4);
122  frame1->Draw();
123  c->cd(2);
124  gPad->SetLeftMargin(0.15);
125  frame2->GetYaxis()->SetTitleOffset(1.4);
126  frame2->Draw();
127 }
RooFormulaVar.h
c
#define c(i)
Definition: RSha256.hxx:101
RooPlot::Draw
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:691
RooGaussian.h
RooGaussian
Plain Gaussian p.d.f.
Definition: RooGaussian.h:24
TCanvas.h
RooGenericPdf.h
RooDataSet.h
RooPlot::frame
static RooPlot * frame(const RooAbsRealLValue &var, Double_t xmin, Double_t xmax, Int_t nBins)
Create a new frame for a given variable in x.
Definition: RooPlot.cxx:249
RooFormulaVar
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Definition: RooFormulaVar.h:30
RooFit
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
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
RooPlot.h
RooPlot::GetYaxis
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1258
RooPlot
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition: RooPlot.h:44
RooRealVar.h
kRed
@ kRed
Definition: Rtypes.h:66
RooFit::LineColor
RooCmdArg LineColor(Color_t color)
Definition: RooGlobalFunc.cxx:57
RooFFTConvPdf
PDF for the numerical (FFT) convolution of two PDFs.
Definition: RooFFTConvPdf.h:25
TCanvas
The Canvas class.
Definition: TCanvas.h:23
TAxis.h
RooGenericPdf
RooGenericPdf is a concrete implementation of a probability density function, which takes a RooArgLis...
Definition: RooGenericPdf.h:25
RooFFTConvPdf.h
gPad
#define gPad
Definition: TVirtualPad.h:287
RooDataSet
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:33
TAttAxis::SetTitleOffset
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition: TAttAxis.cxx:293
RooRealVar
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:37
TH1.h
RooFit::Title
RooCmdArg Title(const char *name)
Definition: RooGlobalFunc.cxx:176