From $ROOTSYS/tutorials/roofit/rf708_bphysics.C

//////////////////////////////////////////////////////////////////////////
//
// 'SPECIAL PDFS' RooFit tutorial macro #708
// 
// Special decay pdf for B physics with mixing and/or CP violation
//
//
//
// 07/2008 - Wouter Verkerke 
// 
/////////////////////////////////////////////////////////////////////////

#ifndef __CINT__
#include "RooGlobalFunc.h"
#endif
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooConstVar.h"
#include "RooCategory.h"
#include "RooBMixDecay.h"
#include "RooBCPEffDecay.h"
#include "RooBDecay.h"
#include "RooFormulaVar.h"
#include "RooTruthModel.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
using namespace RooFit ;

void rf708_bphysics()
{
  ////////////////////////////////////////////////////
  // B - D e c a y   w i t h   m i x i n g          //
  ////////////////////////////////////////////////////

  // C o n s t r u c t   p d f 
  // -------------------------
  
  // Observable
  RooRealVar dt("dt","dt",-10,10) ;
  dt.setBins(40) ;

  // Parameters
  RooRealVar dm("dm","delta m(B0)",0.472) ;
  RooRealVar tau("tau","tau (B0)",1.547) ;
  RooRealVar w("w","flavour mistag rate",0.1) ;
  RooRealVar dw("dw","delta mistag rate for B0/B0bar",0.1) ;

  RooCategory mixState("mixState","B0/B0bar mixing state") ;
  mixState.defineType("mixed",-1) ;
  mixState.defineType("unmixed",1) ;

  RooCategory tagFlav("tagFlav","Flavour of the tagged B0") ;
  tagFlav.defineType("B0",1) ;
  tagFlav.defineType("B0bar",-1) ;

  // Use delta function resolution model
  RooTruthModel tm("tm","truth model",dt) ;

  // Construct Bdecay with mixing
  RooBMixDecay bmix("bmix","decay",dt,mixState,tagFlav,tau,dm,w,dw,tm,RooBMixDecay::DoubleSided) ;



  // P l o t   p d f   i n   v a r i o u s   s l i c e s
  // ---------------------------------------------------

  // Generate some data
  RooDataSet* data = bmix.generate(RooArgSet(dt,mixState,tagFlav),10000) ;

  // Plot B0 and B0bar tagged data separately
  // For all plots below B0 and B0 tagged data will look somewhat differently
  // if the flavor tagging mistag rate for B0 and B0 is different (i.e. dw!=0)
  RooPlot* frame1 = dt.frame(Title("B decay distribution with mixing (B0/B0bar)")) ;

  data->plotOn(frame1,Cut("tagFlav==tagFlav::B0")) ;
  bmix.plotOn(frame1,Slice(tagFlav,"B0")) ;

  data->plotOn(frame1,Cut("tagFlav==tagFlav::B0bar"),MarkerColor(kCyan)) ;
  bmix.plotOn(frame1,Slice(tagFlav,"B0bar"),LineColor(kCyan)) ;


  // Plot mixed slice for B0 and B0bar tagged data separately
  RooPlot* frame2 = dt.frame(Title("B decay distribution of mixed events (B0/B0bar)")) ;

  data->plotOn(frame2,Cut("mixState==mixState::mixed&&tagFlav==tagFlav::B0")) ;
  bmix.plotOn(frame2,Slice(tagFlav,"B0"),Slice(mixState,"mixed")) ;

  data->plotOn(frame2,Cut("mixState==mixState::mixed&&tagFlav==tagFlav::B0bar"),MarkerColor(kCyan)) ;
  bmix.plotOn(frame2,Slice(tagFlav,"B0bar"),Slice(mixState,"mixed"),LineColor(kCyan)) ;


  // Plot unmixed slice for B0 and B0bar tagged data separately
  RooPlot* frame3 = dt.frame(Title("B decay distribution of unmixed events (B0/B0bar)")) ;

  data->plotOn(frame3,Cut("mixState==mixState::unmixed&&tagFlav==tagFlav::B0")) ;
  bmix.plotOn(frame3,Slice(tagFlav,"B0"),Slice(mixState,"unmixed")) ;

  data->plotOn(frame3,Cut("mixState==mixState::unmixed&&tagFlav==tagFlav::B0bar"),MarkerColor(kCyan)) ;
  bmix.plotOn(frame3,Slice(tagFlav,"B0bar"),Slice(mixState,"unmixed"),LineColor(kCyan)) ;





  ///////////////////////////////////////////////////////
  // B - D e c a y   w i t h   C P   v i o l a t i o n //
  ///////////////////////////////////////////////////////

  // C o n s t r u c t   p d f 
  // -------------------------
  
  // Additional parameters needed for B decay with CPV
  RooRealVar CPeigen("CPeigen","CP eigen value",-1) ;
  RooRealVar absLambda("absLambda","|lambda|",1,0,2) ;
  RooRealVar argLambda("absLambda","|lambda|",0.7,-1,1) ;
  RooRealVar effR("effR","B0/B0bar reco efficiency ratio",1) ;

  // Construct Bdecay with CP violation
  RooBCPEffDecay bcp("bcp","bcp", dt, tagFlav, tau, dm, w, CPeigen, absLambda, argLambda, effR, dw, tm, RooBCPEffDecay::DoubleSided) ;
      


  // P l o t   s c e n a r i o   1   -   s i n ( 2 b )   =   0 . 7 ,   | l | = 1 
  // ---------------------------------------------------------------------------

  // Generate some data
  RooDataSet* data2 = bcp.generate(RooArgSet(dt,tagFlav),10000) ;

  // Plot B0 and B0bar tagged data separately
  RooPlot* frame4 = dt.frame(Title("B decay distribution with CPV(|l|=1,Im(l)=0.7) (B0/B0bar)")) ;

  data2->plotOn(frame4,Cut("tagFlav==tagFlav::B0")) ;
  bcp.plotOn(frame4,Slice(tagFlav,"B0")) ;

  data2->plotOn(frame4,Cut("tagFlav==tagFlav::B0bar"),MarkerColor(kCyan)) ;
  bcp.plotOn(frame4,Slice(tagFlav,"B0bar"),LineColor(kCyan)) ;



  // P l o t   s c e n a r i o   2   -   s i n ( 2 b )   =   0 . 7 ,   | l | = 0 . 7 
  // -------------------------------------------------------------------------------

  absLambda=0.7 ;

  // Generate some data
  RooDataSet* data3 = bcp.generate(RooArgSet(dt,tagFlav),10000) ;

  // Plot B0 and B0bar tagged data separately (sin2b = 0.7 plus direct CPV |l|=0.5)
  RooPlot* frame5 = dt.frame(Title("B decay distribution with CPV(|l|=0.7,Im(l)=0.7) (B0/B0bar)")) ;  

  data3->plotOn(frame5,Cut("tagFlav==tagFlav::B0")) ;
  bcp.plotOn(frame5,Slice(tagFlav,"B0")) ;

  data3->plotOn(frame5,Cut("tagFlav==tagFlav::B0bar"),MarkerColor(kCyan)) ;
  bcp.plotOn(frame5,Slice(tagFlav,"B0bar"),LineColor(kCyan)) ;



  //////////////////////////////////////////////////////////////////////////////////
  // 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  //
  //////////////////////////////////////////////////////////////////////////////////

  // C o n s t r u c t   p d f 
  // -------------------------
  
  // Model parameters
  RooRealVar DGbG("DGbG","DGamma/GammaAvg",0.5,-1,1);
  RooRealVar Adir("Adir","-[1-abs(l)**2]/[1+abs(l)**2]",0);
  RooRealVar Amix("Amix","2Im(l)/[1+abs(l)**2]",0.7);
  RooRealVar Adel("Adel","2Re(l)/[1+abs(l)**2]",0.7);
  
  // Derived input parameters for pdf
  RooFormulaVar DG("DG","Delta Gamma","@1/@0",RooArgList(tau,DGbG));
  
  // Construct coefficient functions for sin,cos,sinh modulations of decay distribution
  RooFormulaVar fsin("fsin","fsin","@0*@1*(1-2*@2)",RooArgList(Amix,tagFlav,w));
  RooFormulaVar fcos("fcos","fcos","@0*@1*(1-2*@2)",RooArgList(Adir,tagFlav,w));
  RooFormulaVar fsinh("fsinh","fsinh","@0",RooArgList(Adel));
  
  // Construct generic B decay pdf using above user coefficients
  RooBDecay bcpg("bcpg","bcpg",dt,tau,DG,RooConst(1),fsinh,fcos,fsin,dm,tm,RooBDecay::DoubleSided);
  
  
  
  // P l o t   -   I m ( l ) = 0 . 7 ,   R e ( l ) = 0 . 7   | l | = 1,   d G / G = 0 . 5 
  // -------------------------------------------------------------------------------------
  
  // Generate some data
  RooDataSet* data4 = bcpg.generate(RooArgSet(dt,tagFlav),10000) ;
  
  // Plot B0 and B0bar tagged data separately 
  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)")) ;  
  
  data4->plotOn(frame6,Cut("tagFlav==tagFlav::B0")) ;
  bcpg.plotOn(frame6,Slice(tagFlav,"B0")) ;
  
  data4->plotOn(frame6,Cut("tagFlav==tagFlav::B0bar"),MarkerColor(kCyan)) ;
  bcpg.plotOn(frame6,Slice(tagFlav,"B0bar"),LineColor(kCyan)) ;
  
 
 

  TCanvas* c = new TCanvas("rf708_bphysics","rf708_bphysics",1200,800) ;
  c->Divide(3,2) ;
  c->cd(1) ; gPad->SetLeftMargin(0.15) ; frame1->GetYaxis()->SetTitleOffset(1.6) ; frame1->Draw() ;
  c->cd(2) ; gPad->SetLeftMargin(0.15) ; frame2->GetYaxis()->SetTitleOffset(1.6) ; frame2->Draw() ;
  c->cd(3) ; gPad->SetLeftMargin(0.15) ; frame3->GetYaxis()->SetTitleOffset(1.6) ; frame3->Draw() ;
  c->cd(4) ; gPad->SetLeftMargin(0.15) ; frame4->GetYaxis()->SetTitleOffset(1.6) ; frame4->Draw() ;
  c->cd(5) ; gPad->SetLeftMargin(0.15) ; frame5->GetYaxis()->SetTitleOffset(1.6) ; frame5->Draw() ;
  c->cd(6) ; gPad->SetLeftMargin(0.15) ; frame6->GetYaxis()->SetTitleOffset(1.6) ; frame6->Draw() ;
  
}
 rf708_bphysics.C:1
 rf708_bphysics.C:2
 rf708_bphysics.C:3
 rf708_bphysics.C:4
 rf708_bphysics.C:5
 rf708_bphysics.C:6
 rf708_bphysics.C:7
 rf708_bphysics.C:8
 rf708_bphysics.C:9
 rf708_bphysics.C:10
 rf708_bphysics.C:11
 rf708_bphysics.C:12
 rf708_bphysics.C:13
 rf708_bphysics.C:14
 rf708_bphysics.C:15
 rf708_bphysics.C:16
 rf708_bphysics.C:17
 rf708_bphysics.C:18
 rf708_bphysics.C:19
 rf708_bphysics.C:20
 rf708_bphysics.C:21
 rf708_bphysics.C:22
 rf708_bphysics.C:23
 rf708_bphysics.C:24
 rf708_bphysics.C:25
 rf708_bphysics.C:26
 rf708_bphysics.C:27
 rf708_bphysics.C:28
 rf708_bphysics.C:29
 rf708_bphysics.C:30
 rf708_bphysics.C:31
 rf708_bphysics.C:32
 rf708_bphysics.C:33
 rf708_bphysics.C:34
 rf708_bphysics.C:35
 rf708_bphysics.C:36
 rf708_bphysics.C:37
 rf708_bphysics.C:38
 rf708_bphysics.C:39
 rf708_bphysics.C:40
 rf708_bphysics.C:41
 rf708_bphysics.C:42
 rf708_bphysics.C:43
 rf708_bphysics.C:44
 rf708_bphysics.C:45
 rf708_bphysics.C:46
 rf708_bphysics.C:47
 rf708_bphysics.C:48
 rf708_bphysics.C:49
 rf708_bphysics.C:50
 rf708_bphysics.C:51
 rf708_bphysics.C:52
 rf708_bphysics.C:53
 rf708_bphysics.C:54
 rf708_bphysics.C:55
 rf708_bphysics.C:56
 rf708_bphysics.C:57
 rf708_bphysics.C:58
 rf708_bphysics.C:59
 rf708_bphysics.C:60
 rf708_bphysics.C:61
 rf708_bphysics.C:62
 rf708_bphysics.C:63
 rf708_bphysics.C:64
 rf708_bphysics.C:65
 rf708_bphysics.C:66
 rf708_bphysics.C:67
 rf708_bphysics.C:68
 rf708_bphysics.C:69
 rf708_bphysics.C:70
 rf708_bphysics.C:71
 rf708_bphysics.C:72
 rf708_bphysics.C:73
 rf708_bphysics.C:74
 rf708_bphysics.C:75
 rf708_bphysics.C:76
 rf708_bphysics.C:77
 rf708_bphysics.C:78
 rf708_bphysics.C:79
 rf708_bphysics.C:80
 rf708_bphysics.C:81
 rf708_bphysics.C:82
 rf708_bphysics.C:83
 rf708_bphysics.C:84
 rf708_bphysics.C:85
 rf708_bphysics.C:86
 rf708_bphysics.C:87
 rf708_bphysics.C:88
 rf708_bphysics.C:89
 rf708_bphysics.C:90
 rf708_bphysics.C:91
 rf708_bphysics.C:92
 rf708_bphysics.C:93
 rf708_bphysics.C:94
 rf708_bphysics.C:95
 rf708_bphysics.C:96
 rf708_bphysics.C:97
 rf708_bphysics.C:98
 rf708_bphysics.C:99
 rf708_bphysics.C:100
 rf708_bphysics.C:101
 rf708_bphysics.C:102
 rf708_bphysics.C:103
 rf708_bphysics.C:104
 rf708_bphysics.C:105
 rf708_bphysics.C:106
 rf708_bphysics.C:107
 rf708_bphysics.C:108
 rf708_bphysics.C:109
 rf708_bphysics.C:110
 rf708_bphysics.C:111
 rf708_bphysics.C:112
 rf708_bphysics.C:113
 rf708_bphysics.C:114
 rf708_bphysics.C:115
 rf708_bphysics.C:116
 rf708_bphysics.C:117
 rf708_bphysics.C:118
 rf708_bphysics.C:119
 rf708_bphysics.C:120
 rf708_bphysics.C:121
 rf708_bphysics.C:122
 rf708_bphysics.C:123
 rf708_bphysics.C:124
 rf708_bphysics.C:125
 rf708_bphysics.C:126
 rf708_bphysics.C:127
 rf708_bphysics.C:128
 rf708_bphysics.C:129
 rf708_bphysics.C:130
 rf708_bphysics.C:131
 rf708_bphysics.C:132
 rf708_bphysics.C:133
 rf708_bphysics.C:134
 rf708_bphysics.C:135
 rf708_bphysics.C:136
 rf708_bphysics.C:137
 rf708_bphysics.C:138
 rf708_bphysics.C:139
 rf708_bphysics.C:140
 rf708_bphysics.C:141
 rf708_bphysics.C:142
 rf708_bphysics.C:143
 rf708_bphysics.C:144
 rf708_bphysics.C:145
 rf708_bphysics.C:146
 rf708_bphysics.C:147
 rf708_bphysics.C:148
 rf708_bphysics.C:149
 rf708_bphysics.C:150
 rf708_bphysics.C:151
 rf708_bphysics.C:152
 rf708_bphysics.C:153
 rf708_bphysics.C:154
 rf708_bphysics.C:155
 rf708_bphysics.C:156
 rf708_bphysics.C:157
 rf708_bphysics.C:158
 rf708_bphysics.C:159
 rf708_bphysics.C:160
 rf708_bphysics.C:161
 rf708_bphysics.C:162
 rf708_bphysics.C:163
 rf708_bphysics.C:164
 rf708_bphysics.C:165
 rf708_bphysics.C:166
 rf708_bphysics.C:167
 rf708_bphysics.C:168
 rf708_bphysics.C:169
 rf708_bphysics.C:170
 rf708_bphysics.C:171
 rf708_bphysics.C:172
 rf708_bphysics.C:173
 rf708_bphysics.C:174
 rf708_bphysics.C:175
 rf708_bphysics.C:176
 rf708_bphysics.C:177
 rf708_bphysics.C:178
 rf708_bphysics.C:179
 rf708_bphysics.C:180
 rf708_bphysics.C:181
 rf708_bphysics.C:182
 rf708_bphysics.C:183
 rf708_bphysics.C:184
 rf708_bphysics.C:185
 rf708_bphysics.C:186
 rf708_bphysics.C:187
 rf708_bphysics.C:188
 rf708_bphysics.C:189
 rf708_bphysics.C:190
 rf708_bphysics.C:191
 rf708_bphysics.C:192
 rf708_bphysics.C:193
 rf708_bphysics.C:194
 rf708_bphysics.C:195
 rf708_bphysics.C:196
 rf708_bphysics.C:197
 rf708_bphysics.C:198
 rf708_bphysics.C:199
 rf708_bphysics.C:200
 rf708_bphysics.C:201
 rf708_bphysics.C:202
 rf708_bphysics.C:203
 rf708_bphysics.C:204
 rf708_bphysics.C:205
 rf708_bphysics.C:206
 rf708_bphysics.C:207
 rf708_bphysics.C:208
 rf708_bphysics.C:209
 rf708_bphysics.C:210
 rf708_bphysics.C:211
 rf708_bphysics.C:212
 rf708_bphysics.C:213
 rf708_bphysics.C:214
 rf708_bphysics.C:215