Logo ROOT  
Reference Guide
rf708_bphysics.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook -js
4 /// Special pdf's: special decay pdf for B physics with mixing and/or CP violation
5 ///
6 /// \macro_image
7 /// \macro_output
8 /// \macro_code
9 ///
10 /// \date July 2008
11 /// \author Wouter Verkerke
12 
13 #include "RooRealVar.h"
14 #include "RooDataSet.h"
15 #include "RooGaussian.h"
16 #include "RooConstVar.h"
17 #include "RooCategory.h"
18 #include "RooBMixDecay.h"
19 #include "RooBCPEffDecay.h"
20 #include "RooBDecay.h"
21 #include "RooFormulaVar.h"
22 #include "RooTruthModel.h"
23 #include "TCanvas.h"
24 #include "TAxis.h"
25 #include "RooPlot.h"
26 using namespace RooFit;
27 
28 void rf708_bphysics()
29 {
30  // -------------------------------------
31  // B - D e c a y w i t h m i x i n g
32  // =====================================
33 
34  // C o n s t r u c t p d f
35  // -------------------------
36 
37  // Observable
38  RooRealVar dt("dt", "dt", -10, 10);
39  dt.setBins(40);
40 
41  // Parameters
42  RooRealVar dm("dm", "delta m(B0)", 0.472);
43  RooRealVar tau("tau", "tau (B0)", 1.547);
44  RooRealVar w("w", "flavour mistag rate", 0.1);
45  RooRealVar dw("dw", "delta mistag rate for B0/B0bar", 0.1);
46 
47  RooCategory mixState("mixState", "B0/B0bar mixing state");
48  mixState.defineType("mixed", -1);
49  mixState.defineType("unmixed", 1);
50 
51  RooCategory tagFlav("tagFlav", "Flavour of the tagged B0");
52  tagFlav.defineType("B0", 1);
53  tagFlav.defineType("B0bar", -1);
54 
55  // Use delta function resolution model
56  RooTruthModel truthModel("tm", "truth model", dt);
57 
58  // Construct Bdecay with mixing
59  RooBMixDecay bmix("bmix", "decay", dt, mixState, tagFlav, tau, dm, w, dw, truthModel, RooBMixDecay::DoubleSided);
60 
61  // P l o t p d f i n v a r i o u s s l i c e s
62  // ---------------------------------------------------
63 
64  // Generate some data
65  RooDataSet *data = bmix.generate(RooArgSet(dt, mixState, tagFlav), 10000);
66 
67  // Plot B0 and B0bar tagged data separately
68  // For all plots below B0 and B0 tagged data will look somewhat differently
69  // if the flavor tagging mistag rate for B0 and B0 is different (i.e. dw!=0)
70  RooPlot *frame1 = dt.frame(Title("B decay distribution with mixing (B0/B0bar)"));
71 
72  data->plotOn(frame1, Cut("tagFlav==tagFlav::B0"));
73  bmix.plotOn(frame1, Slice(tagFlav, "B0"));
74 
75  data->plotOn(frame1, Cut("tagFlav==tagFlav::B0bar"), MarkerColor(kCyan));
76  bmix.plotOn(frame1, Slice(tagFlav, "B0bar"), LineColor(kCyan));
77 
78  // Plot mixed slice for B0 and B0bar tagged data separately
79  RooPlot *frame2 = dt.frame(Title("B decay distribution of mixed events (B0/B0bar)"));
80 
81  data->plotOn(frame2, Cut("mixState==mixState::mixed&&tagFlav==tagFlav::B0"));
82  bmix.plotOn(frame2, Slice(tagFlav, "B0"), Slice(mixState, "mixed"));
83 
84  data->plotOn(frame2, Cut("mixState==mixState::mixed&&tagFlav==tagFlav::B0bar"), MarkerColor(kCyan));
85  bmix.plotOn(frame2, Slice(tagFlav, "B0bar"), Slice(mixState, "mixed"), LineColor(kCyan));
86 
87  // Plot unmixed slice for B0 and B0bar tagged data separately
88  RooPlot *frame3 = dt.frame(Title("B decay distribution of unmixed events (B0/B0bar)"));
89 
90  data->plotOn(frame3, Cut("mixState==mixState::unmixed&&tagFlav==tagFlav::B0"));
91  bmix.plotOn(frame3, Slice(tagFlav, "B0"), Slice(mixState, "unmixed"));
92 
93  data->plotOn(frame3, Cut("mixState==mixState::unmixed&&tagFlav==tagFlav::B0bar"), MarkerColor(kCyan));
94  bmix.plotOn(frame3, Slice(tagFlav, "B0bar"), Slice(mixState, "unmixed"), LineColor(kCyan));
95 
96  // -------------------------------------------------
97  // B - D e c a y w i t h C P v i o l a t i o n
98  // =================================================
99 
100  // C o n s t r u c t p d f
101  // -------------------------
102 
103  // Additional parameters needed for B decay with CPV
104  RooRealVar CPeigen("CPeigen", "CP eigen value", -1);
105  RooRealVar absLambda("absLambda", "|lambda|", 1, 0, 2);
106  RooRealVar argLambda("absLambda", "|lambda|", 0.7, -1, 1);
107  RooRealVar effR("effR", "B0/B0bar reco efficiency ratio", 1);
108 
109  // Construct Bdecay with CP violation
110  RooBCPEffDecay bcp("bcp", "bcp", dt, tagFlav, tau, dm, w, CPeigen, absLambda, argLambda, effR, dw, truthModel,
112 
113  // P l o t s c e n a r i o 1 - s i n ( 2 b ) = 0 . 7 , | l | = 1
114  // ---------------------------------------------------------------------------
115 
116  // Generate some data
117  RooDataSet *data2 = bcp.generate(RooArgSet(dt, tagFlav), 10000);
118 
119  // Plot B0 and B0bar tagged data separately
120  RooPlot *frame4 = dt.frame(Title("B decay distribution with CPV(|l|=1,Im(l)=0.7) (B0/B0bar)"));
121 
122  data2->plotOn(frame4, Cut("tagFlav==tagFlav::B0"));
123  bcp.plotOn(frame4, Slice(tagFlav, "B0"));
124 
125  data2->plotOn(frame4, Cut("tagFlav==tagFlav::B0bar"), MarkerColor(kCyan));
126  bcp.plotOn(frame4, Slice(tagFlav, "B0bar"), LineColor(kCyan));
127 
128  // P l o t s c e n a r i o 2 - s i n ( 2 b ) = 0 . 7 , | l | = 0 . 7
129  // -------------------------------------------------------------------------------
130 
131  absLambda = 0.7;
132 
133  // Generate some data
134  RooDataSet *data3 = bcp.generate(RooArgSet(dt, tagFlav), 10000);
135 
136  // Plot B0 and B0bar tagged data separately (sin2b = 0.7 plus direct CPV |l|=0.5)
137  RooPlot *frame5 = dt.frame(Title("B decay distribution with CPV(|l|=0.7,Im(l)=0.7) (B0/B0bar)"));
138 
139  data3->plotOn(frame5, Cut("tagFlav==tagFlav::B0"));
140  bcp.plotOn(frame5, Slice(tagFlav, "B0"));
141 
142  data3->plotOn(frame5, Cut("tagFlav==tagFlav::B0bar"), MarkerColor(kCyan));
143  bcp.plotOn(frame5, Slice(tagFlav, "B0bar"), LineColor(kCyan));
144 
145  // ---------------------------------------------------------------------------
146  // G e n e r i c B d e c a y w i t h u s e r c o e f f i c i e n t s
147  // ===========================================================================
148 
149  // C o n s t r u c t p d f
150  // -------------------------
151 
152  // Model parameters
153  RooRealVar DGbG("DGbG", "DGamma/GammaAvg", 0.5, -1, 1);
154  RooRealVar Adir("Adir", "-[1-abs(l)**2]/[1+abs(l)**2]", 0);
155  RooRealVar Amix("Amix", "2Im(l)/[1+abs(l)**2]", 0.7);
156  RooRealVar Adel("Adel", "2Re(l)/[1+abs(l)**2]", 0.7);
157 
158  // Derived input parameters for pdf
159  RooFormulaVar DG("DG", "Delta Gamma", "@1/@0", RooArgList(tau, DGbG));
160 
161  // Construct coefficient functions for sin,cos,sinh modulations of decay distribution
162  RooFormulaVar fsin("fsin", "fsin", "@0*@1*(1-2*@2)", RooArgList(Amix, tagFlav, w));
163  RooFormulaVar fcos("fcos", "fcos", "@0*@1*(1-2*@2)", RooArgList(Adir, tagFlav, w));
164  RooFormulaVar fsinh("fsinh", "fsinh", "@0", RooArgList(Adel));
165 
166  // Construct generic B decay pdf using above user coefficients
167  RooBDecay bcpg("bcpg", "bcpg", dt, tau, DG, RooConst(1), fsinh, fcos, fsin, dm, truthModel, RooBDecay::DoubleSided);
168 
169  // P l o t - I m ( l ) = 0 . 7 , R e ( l ) = 0 . 7 | l | = 1, d G / G = 0 . 5
170  // -------------------------------------------------------------------------------------
171 
172  // Generate some data
173  RooDataSet *data4 = bcpg.generate(RooArgSet(dt, tagFlav), 10000);
174 
175  // Plot B0 and B0bar tagged data separately
176  RooPlot *frame6 = dt.frame(Title("B decay distribution with CPV(Im(l)=0.7,Re(l)=0.7,|l|=1,dG/G=0.5) (B0/B0bar)"));
177 
178  data4->plotOn(frame6, Cut("tagFlav==tagFlav::B0"));
179  bcpg.plotOn(frame6, Slice(tagFlav, "B0"));
180 
181  data4->plotOn(frame6, Cut("tagFlav==tagFlav::B0bar"), MarkerColor(kCyan));
182  bcpg.plotOn(frame6, Slice(tagFlav, "B0bar"), LineColor(kCyan));
183 
184  TCanvas *c = new TCanvas("rf708_bphysics", "rf708_bphysics", 1200, 800);
185  c->Divide(3, 2);
186  c->cd(1);
187  gPad->SetLeftMargin(0.15);
188  frame1->GetYaxis()->SetTitleOffset(1.6);
189  frame1->Draw();
190  c->cd(2);
191  gPad->SetLeftMargin(0.15);
192  frame2->GetYaxis()->SetTitleOffset(1.6);
193  frame2->Draw();
194  c->cd(3);
195  gPad->SetLeftMargin(0.15);
196  frame3->GetYaxis()->SetTitleOffset(1.6);
197  frame3->Draw();
198  c->cd(4);
199  gPad->SetLeftMargin(0.15);
200  frame4->GetYaxis()->SetTitleOffset(1.6);
201  frame4->Draw();
202  c->cd(5);
203  gPad->SetLeftMargin(0.15);
204  frame5->GetYaxis()->SetTitleOffset(1.6);
205  frame5->Draw();
206  c->cd(6);
207  gPad->SetLeftMargin(0.15);
208  frame6->GetYaxis()->SetTitleOffset(1.6);
209  frame6->Draw();
210 }
RooFormulaVar.h
c
#define c(i)
Definition: RSha256.hxx:101
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
RooBMixDecay
Class RooBMixDecay is a RooAbsAnaConvPdf implementation that describes the decay of B mesons with the...
Definition: RooBMixDecay.h:23
RooBDecay.h
RooArgList
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooGaussian.h
RooBCPEffDecay.h
TCanvas.h
RooDataSet.h
RooBMixDecay.h
RooFit::Slice
RooCmdArg Slice(const RooArgSet &sliceSet)
Definition: RooGlobalFunc.cxx:40
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
kCyan
@ kCyan
Definition: Rtypes.h:66
RooFit::Cut
RooCmdArg Cut(const char *cutSpec)
Definition: RooGlobalFunc.cxx:81
RooFormulaVar
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Definition: RooFormulaVar.h:30
RooBCPEffDecay::DoubleSided
@ DoubleSided
Definition: RooBCPEffDecay.h:26
RooFit
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition: RooCFunction1Binding.h:29
RooBDecay
Most general description of B decay time distribution with effects of CP violation,...
Definition: RooBDecay.h:25
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
RooCategory.h
RooRealVar.h
RooConstVar.h
rf708_bphysics
Definition: rf708_bphysics.py:1
RooFit::LineColor
RooCmdArg LineColor(Color_t color)
Definition: RooGlobalFunc.cxx:57
RooBCPEffDecay
PDF describing decay time distribution of B meson including effects of standard model CP violation.
Definition: RooBCPEffDecay.h:23
TCanvas
The Canvas class.
Definition: TCanvas.h:23
RooBMixDecay::DoubleSided
@ DoubleSided
Definition: RooBMixDecay.h:26
RooCategory
RooCategory is an object to represent discrete states.
Definition: RooCategory.h:27
TAxis.h
RooBDecay::DoubleSided
@ DoubleSided
Definition: RooBDecay.h:29
gPad
#define gPad
Definition: TVirtualPad.h:287
RooDataSet
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:33
RooFit::MarkerColor
RooCmdArg MarkerColor(Color_t color)
Definition: RooGlobalFunc.cxx:88
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
RooFit::Title
RooCmdArg Title(const char *name)
Definition: RooGlobalFunc.cxx:176
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:29
RooTruthModel
RooTruthModel is an implementation of RooResolution model that provides a delta-function resolution m...
Definition: RooTruthModel.h:21
RooFit::RooConst
RooConstVar & RooConst(Double_t val)
Definition: RooGlobalFunc.cxx:347