Logo ROOT   6.18/05
Reference Guide
rf310_sliceplot.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_roofit
3/// \notebook
4/// Multidimensional models: projecting p.d.f and data slices in discrete observables
5///
6/// \macro_image
7/// \macro_output
8/// \macro_code
9/// \author 07/2008 - Wouter Verkerke
10
11#include "RooRealVar.h"
12#include "RooDataSet.h"
13#include "RooGaussModel.h"
14#include "RooDecay.h"
15#include "RooBMixDecay.h"
16#include "RooCategory.h"
17#include "TCanvas.h"
18#include "TAxis.h"
19#include "RooPlot.h"
20using namespace RooFit;
21
22void rf310_sliceplot()
23{
24
25 // C r e a t e B d e c a y p d f w it h m i x i n g
26 // ----------------------------------------------------------
27
28 // Decay time observables
29 RooRealVar dt("dt", "dt", -20, 20);
30
31 // Discrete observables mixState (B0tag==B0reco?) and tagFlav (B0tag==B0(bar)?)
32 RooCategory mixState("mixState", "B0/B0bar mixing state");
33 RooCategory tagFlav("tagFlav", "Flavour of the tagged B0");
34
35 // Define state labels of discrete observables
36 mixState.defineType("mixed", -1);
37 mixState.defineType("unmixed", 1);
38 tagFlav.defineType("B0", 1);
39 tagFlav.defineType("B0bar", -1);
40
41 // Model parameters
42 RooRealVar dm("dm", "delta m(B)", 0.472, 0., 1.0);
43 RooRealVar tau("tau", "B0 decay time", 1.547, 1.0, 2.0);
44 RooRealVar w("w", "Flavor Mistag rate", 0.03, 0.0, 1.0);
45 RooRealVar dw("dw", "Flavor Mistag rate difference between B0 and B0bar", 0.01);
46
47 // Build a gaussian resolution model
48 RooRealVar bias1("bias1", "bias1", 0);
49 RooRealVar sigma1("sigma1", "sigma1", 0.01);
50 RooGaussModel gm1("gm1", "gauss model 1", dt, bias1, sigma1);
51
52 // Construct a decay pdf, smeared with single gaussian resolution model
53 RooBMixDecay bmix_gm1("bmix", "decay", dt, mixState, tagFlav, tau, dm, w, dw, gm1, RooBMixDecay::DoubleSided);
54
55 // Generate BMixing data with above set of event errors
56 RooDataSet *data = bmix_gm1.generate(RooArgSet(dt, tagFlav, mixState), 20000);
57
58 // P l o t f u l l d e c a y d i s t r i b u t i o n
59 // ----------------------------------------------------------
60
61 // Create frame, plot data and pdf projection (integrated over tagFlav and mixState)
62 RooPlot *frame = dt.frame(Title("Inclusive decay distribution"));
63 data->plotOn(frame);
64 bmix_gm1.plotOn(frame);
65
66 // P l o t d e c a y d i s t r . f o r m i x e d a n d u n m i x e d s l i c e o f m i x S t a t e
67 // ------------------------------------------------------------------------------------------------------------------
68
69 // Create frame, plot data (mixed only)
70 RooPlot *frame2 = dt.frame(Title("Decay distribution of mixed events"));
71 data->plotOn(frame2, Cut("mixState==mixState::mixed"));
72
73 // Position slice in mixState at "mixed" and plot slice of pdf in mixstate over data (integrated over tagFlav)
74 bmix_gm1.plotOn(frame2, Slice(mixState, "mixed"));
75
76 // Create frame, plot data (unmixed only)
77 RooPlot *frame3 = dt.frame(Title("Decay distribution of unmixed events"));
78 data->plotOn(frame3, Cut("mixState==mixState::unmixed"));
79
80 // Position slice in mixState at "unmixed" and plot slice of pdf in mixstate over data (integrated over tagFlav)
81 bmix_gm1.plotOn(frame3, Slice(mixState, "unmixed"));
82
83 TCanvas *c = new TCanvas("rf310_sliceplot", "rf310_sliceplot", 1200, 400);
84 c->Divide(3);
85 c->cd(1);
86 gPad->SetLeftMargin(0.15);
87 frame->GetYaxis()->SetTitleOffset(1.4);
88 gPad->SetLogy();
89 frame->Draw();
90 c->cd(2);
91 gPad->SetLeftMargin(0.15);
92 frame2->GetYaxis()->SetTitleOffset(1.4);
93 gPad->SetLogy();
94 frame2->Draw();
95 c->cd(3);
96 gPad->SetLeftMargin(0.15);
97 frame3->GetYaxis()->SetTitleOffset(1.4);
98 gPad->SetLogy();
99 frame3->Draw();
100}
#define c(i)
Definition: RSha256.hxx:101
#define gPad
Definition: TVirtualPad.h:286
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
Calls RooPlot* plotOn(RooPlot* frame, const RooLinkedList& cmdList) const ;.
Definition: RooAbsData.cxx:552
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
Class RooBMixDecay is a RooAbsAnaConvPdf implementation that describes the decay of B mesons with the...
Definition: RooBMixDecay.h:23
RooCategory represents a fundamental (non-derived) discrete value object.
Definition: RooCategory.h:24
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:31
Class RooGaussModel implements a RooResolutionModel that models a Gaussian distribution.
Definition: RooGaussModel.h:26
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition: RooPlot.h:41
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1123
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:558
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
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:294
The Canvas class.
Definition: TCanvas.h:31
Template specialisation used in RooAbsArg:
RooCmdArg Slice(const RooArgSet &sliceSet)
RooCmdArg Cut(const char *cutSpec)
const char * Title
Definition: TXMLSetup.cxx:67