Logo ROOT  
Reference Guide
rf310_sliceplot.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook
4 ///
5 /// Multidimensional models: projecting p.d.f and data slices in discrete observables
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 "RooGaussModel.h"
17 #include "RooDecay.h"
18 #include "RooBMixDecay.h"
19 #include "RooCategory.h"
20 #include "TCanvas.h"
21 #include "TAxis.h"
22 #include "RooPlot.h"
23 using namespace RooFit;
24 
25 void rf310_sliceplot()
26 {
27 
28  // C r e a t e B d e c a y p d f w it h m i x i n g
29  // ----------------------------------------------------------
30 
31  // Decay time observables
32  RooRealVar dt("dt", "dt", -20, 20);
33 
34  // Discrete observables mixState (B0tag==B0reco?) and tagFlav (B0tag==B0(bar)?)
35  RooCategory mixState("mixState", "B0/B0bar mixing state");
36  RooCategory tagFlav("tagFlav", "Flavour of the tagged B0");
37 
38  // Define state labels of discrete observables
39  mixState.defineType("mixed", -1);
40  mixState.defineType("unmixed", 1);
41  tagFlav.defineType("B0", 1);
42  tagFlav.defineType("B0bar", -1);
43 
44  // Model parameters
45  RooRealVar dm("dm", "delta m(B)", 0.472, 0., 1.0);
46  RooRealVar tau("tau", "B0 decay time", 1.547, 1.0, 2.0);
47  RooRealVar w("w", "Flavor Mistag rate", 0.03, 0.0, 1.0);
48  RooRealVar dw("dw", "Flavor Mistag rate difference between B0 and B0bar", 0.01);
49 
50  // Build a gaussian resolution model
51  RooRealVar bias1("bias1", "bias1", 0);
52  RooRealVar sigma1("sigma1", "sigma1", 0.01);
53  RooGaussModel gm1("gm1", "gauss model 1", dt, bias1, sigma1);
54 
55  // Construct a decay pdf, smeared with single gaussian resolution model
56  RooBMixDecay bmix_gm1("bmix", "decay", dt, mixState, tagFlav, tau, dm, w, dw, gm1, RooBMixDecay::DoubleSided);
57 
58  // Generate BMixing data with above set of event errors
59  RooDataSet *data = bmix_gm1.generate(RooArgSet(dt, tagFlav, mixState), 20000);
60 
61  // 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
62  // ----------------------------------------------------------
63 
64  // Create frame, plot data and pdf projection (integrated over tagFlav and mixState)
65  RooPlot *frame = dt.frame(Title("Inclusive decay distribution"));
66  data->plotOn(frame);
67  bmix_gm1.plotOn(frame);
68 
69  // 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
70  // ------------------------------------------------------------------------------------------------------------------
71 
72  // Create frame, plot data (mixed only)
73  RooPlot *frame2 = dt.frame(Title("Decay distribution of mixed events"));
74  data->plotOn(frame2, Cut("mixState==mixState::mixed"));
75 
76  // Position slice in mixState at "mixed" and plot slice of pdf in mixstate over data (integrated over tagFlav)
77  bmix_gm1.plotOn(frame2, Slice(mixState, "mixed"));
78 
79  // Create frame, plot data (unmixed only)
80  RooPlot *frame3 = dt.frame(Title("Decay distribution of unmixed events"));
81  data->plotOn(frame3, Cut("mixState==mixState::unmixed"));
82 
83  // Position slice in mixState at "unmixed" and plot slice of pdf in mixstate over data (integrated over tagFlav)
84  bmix_gm1.plotOn(frame3, Slice(mixState, "unmixed"));
85 
86  TCanvas *c = new TCanvas("rf310_sliceplot", "rf310_sliceplot", 1200, 400);
87  c->Divide(3);
88  c->cd(1);
89  gPad->SetLeftMargin(0.15);
90  frame->GetYaxis()->SetTitleOffset(1.4);
91  gPad->SetLogy();
92  frame->Draw();
93  c->cd(2);
94  gPad->SetLeftMargin(0.15);
95  frame2->GetYaxis()->SetTitleOffset(1.4);
96  gPad->SetLogy();
97  frame2->Draw();
98  c->cd(3);
99  gPad->SetLeftMargin(0.15);
100  frame3->GetYaxis()->SetTitleOffset(1.4);
101  gPad->SetLogy();
102  frame3->Draw();
103 }
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
RooGaussModel.h
RooBMixDecay
Definition: RooBMixDecay.h:23
RooDecay.h
TCanvas.h
RooDataSet.h
RooBMixDecay.h
RooFit::Slice
RooCmdArg Slice(const RooArgSet &sliceSet)
Definition: RooGlobalFunc.cxx:39
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
RooFit::Cut
RooCmdArg Cut(const char *cutSpec)
Definition: RooGlobalFunc.cxx:80
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
rf310_sliceplot
Definition: rf310_sliceplot.py:1
RooPlot.h
RooPlot::GetYaxis
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1256
RooPlot
Definition: RooPlot.h:44
RooCategory.h
RooRealVar.h
TCanvas
Definition: TCanvas.h:23
RooBMixDecay::DoubleSided
@ DoubleSided
Definition: RooBMixDecay.h:54
RooCategory
Definition: RooCategory.h:27
TAxis.h
gPad
#define gPad
Definition: TVirtualPad.h:287
RooDataSet
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
Definition: RooRealVar.h:35
RooGaussModel
Definition: RooGaussModel.h:26
RooFit::Title
RooCmdArg Title(const char *name)
Definition: RooGlobalFunc.cxx:173
RooArgSet
Definition: RooArgSet.h:28