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