Logo ROOT  
Reference Guide
rf704_amplitudefit.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook
4 ///
5 /// Special p.d.f.'s: using a p.d.f defined by a sum of real-valued amplitude components
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 "RooConstVar.h"
18 #include "RooTruthModel.h"
19 #include "RooFormulaVar.h"
20 #include "RooRealSumPdf.h"
21 #include "RooPolyVar.h"
22 #include "RooProduct.h"
23 #include "TH1.h"
24 #include "TCanvas.h"
25 #include "TAxis.h"
26 #include "RooPlot.h"
27 using namespace RooFit;
28 
29 void rf704_amplitudefit()
30 {
31  // S e t u p 2 D a m p l i t u d e f u n c t i o n s
32  // -------------------------------------------------------
33 
34  // Observables
35  RooRealVar t("t", "time", -1., 15.);
36  RooRealVar cosa("cosa", "cos(alpha)", -1., 1.);
37 
38  // Use RooTruthModel to obtain compiled implementation of sinh/cosh modulated decay functions
39  RooRealVar tau("tau", "#tau", 1.5);
40  RooRealVar deltaGamma("deltaGamma", "deltaGamma", 0.3);
41  RooTruthModel truthModel("tm", "tm", t);
42  RooFormulaVar coshGBasis("coshGBasis", "exp(-@0/ @1)*cosh(@0*@2/2)", RooArgList(t, tau, deltaGamma));
43  RooFormulaVar sinhGBasis("sinhGBasis", "exp(-@0/ @1)*sinh(@0*@2/2)", RooArgList(t, tau, deltaGamma));
44  RooAbsReal *coshGConv = truthModel.convolution(&coshGBasis, &t);
45  RooAbsReal *sinhGConv = truthModel.convolution(&sinhGBasis, &t);
46 
47  // Construct polynomial amplitudes in cos(a)
48  RooPolyVar poly1("poly1", "poly1", cosa, RooArgList(RooConst(0.5), RooConst(0.2), RooConst(0.2)), 0);
49  RooPolyVar poly2("poly2", "poly2", cosa, RooArgList(RooConst(1), RooConst(-0.2), RooConst(3)), 0);
50 
51  // Construct 2D amplitude as uncorrelated product of amp(t)*amp(cosa)
52  RooProduct ampl1("ampl1", "amplitude 1", RooArgSet(poly1, *coshGConv));
53  RooProduct ampl2("ampl2", "amplitude 2", RooArgSet(poly2, *sinhGConv));
54 
55  // C o n s t r u c t a m p l i t u d e s u m p d f
56  // -----------------------------------------------------
57 
58  // Amplitude strengths
59  RooRealVar f1("f1", "f1", 1, 0, 2);
60  RooRealVar f2("f2", "f2", 0.5, 0, 2);
61 
62  // Construct pdf
63  RooRealSumPdf pdf("pdf", "pdf", RooArgList(ampl1, ampl2), RooArgList(f1, f2));
64 
65  // Generate some toy data from pdf
66  RooDataSet *data = pdf.generate(RooArgSet(t, cosa), 10000);
67 
68  // Fit pdf to toy data with only amplitude strength floating
69  pdf.fitTo(*data);
70 
71  // P l o t a m p l i t u d e s u m p d f
72  // -------------------------------------------
73 
74  // Make 2D plots of amplitudes
75  TH1 *hh_cos = ampl1.createHistogram("hh_cos", t, Binning(50), YVar(cosa, Binning(50)));
76  TH1 *hh_sin = ampl2.createHistogram("hh_sin", t, Binning(50), YVar(cosa, Binning(50)));
77  hh_cos->SetLineColor(kBlue);
78  hh_sin->SetLineColor(kRed);
79 
80  // Make projection on t, plot data, pdf and its components
81  // Note component projections may be larger than sum because amplitudes can be negative
82  RooPlot *frame1 = t.frame();
83  data->plotOn(frame1);
84  pdf.plotOn(frame1);
85  pdf.plotOn(frame1, Components(ampl1), LineStyle(kDashed));
86  pdf.plotOn(frame1, Components(ampl2), LineStyle(kDashed), LineColor(kRed));
87 
88  // Make projection on cosa, plot data, pdf and its components
89  // Note that components projection may be larger than sum because amplitudes can be negative
90  RooPlot *frame2 = cosa.frame();
91  data->plotOn(frame2);
92  pdf.plotOn(frame2);
93  pdf.plotOn(frame2, Components(ampl1), LineStyle(kDashed));
94  pdf.plotOn(frame2, Components(ampl2), LineStyle(kDashed), LineColor(kRed));
95 
96  TCanvas *c = new TCanvas("rf704_amplitudefit", "rf704_amplitudefit", 800, 800);
97  c->Divide(2, 2);
98  c->cd(1);
99  gPad->SetLeftMargin(0.15);
100  frame1->GetYaxis()->SetTitleOffset(1.4);
101  frame1->Draw();
102  c->cd(2);
103  gPad->SetLeftMargin(0.15);
104  frame2->GetYaxis()->SetTitleOffset(1.4);
105  frame2->Draw();
106  c->cd(3);
107  gPad->SetLeftMargin(0.20);
108  hh_cos->GetZaxis()->SetTitleOffset(2.3);
109  hh_cos->Draw("surf");
110  c->cd(4);
111  gPad->SetLeftMargin(0.20);
112  hh_sin->GetZaxis()->SetTitleOffset(2.3);
113  hh_sin->Draw("surf");
114 }
RooFormulaVar.h
c
#define c(i)
Definition: RSha256.hxx:119
RooTruthModel.h
RooPlot::Draw
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:691
rf704_amplitudefit
Definition: rf704_amplitudefit.py:1
RooArgList
Definition: RooArgList.h:21
RooGaussian.h
RooAbsReal
Definition: RooAbsReal.h:61
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
TH1::GetZaxis
TAxis * GetZaxis()
Definition: TH1.h:319
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
RooPolyVar.h
RooProduct
Definition: RooProduct.h:30
RooFormulaVar
Definition: RooFormulaVar.h:30
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
RooPlot.h
RooPlot::GetYaxis
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1256
RooPlot
Definition: RooPlot.h:44
RooRealVar.h
kRed
@ kRed
Definition: Rtypes.h:66
RooPolyVar
Definition: RooPolyVar.h:28
RooProduct.h
RooConstVar.h
RooFit::LineColor
RooCmdArg LineColor(Color_t color)
Definition: RooGlobalFunc.cxx:56
f1
TF1 * f1
Definition: legend1.C:11
TCanvas
Definition: TCanvas.h:23
TAxis.h
TH1
Definition: TH1.h:57
kBlue
@ kBlue
Definition: Rtypes.h:66
kDashed
@ kDashed
Definition: TAttLine.h:48
gPad
#define gPad
Definition: TVirtualPad.h:287
RooDataSet
Definition: RooDataSet.h:33
RooFit::LineStyle
RooCmdArg LineStyle(Style_t style)
Definition: RooGlobalFunc.cxx:57
TAttAxis::SetTitleOffset
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition: TAttAxis.cxx:293
RooRealSumPdf.h
RooRealVar
Definition: RooRealVar.h:35
RooFit::Components
RooCmdArg Components(const RooArgSet &compSet)
Definition: RooGlobalFunc.cxx:74
TH1.h
RooRealSumPdf
Definition: RooRealSumPdf.h:25
RooArgSet
Definition: RooArgSet.h:28
RooTruthModel
Definition: RooTruthModel.h:21
TH1::Draw
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2997
RooFit::RooConst
RooConstVar & RooConst(Double_t val)
Definition: RooGlobalFunc.cxx:341