Logo ROOT  
Reference Guide
RooBCPEffDecay.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitModels *
4 * @(#)root/roofit:$Id$
5 * Authors: *
6 * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7 * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8 * *
9 * Copyright (c) 2000-2005, Regents of the University of California *
10 * and Stanford University. All rights reserved. *
11 * *
12 * Redistribution and use in source and binary forms, *
13 * with or without modification, are permitted according to the terms *
14 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15 *****************************************************************************/
16
17/** \class RooBCPEffDecay
18 \ingroup Roofit
19
20PDF describing decay time distribution of B meson including effects of standard model CP violation.
21This function can be analytically convolved with any RooResolutionModel implementation.
22*/
23
24
25#include "RooFit.h"
26
27#include "Riostream.h"
28#include "RooRealVar.h"
29#include "RooRandom.h"
30#include "RooBCPEffDecay.h"
31#include "RooRealIntegral.h"
32
33using namespace std;
34
36
37////////////////////////////////////////////////////////////////////////////////
38/// Constructor.
39
40RooBCPEffDecay::RooBCPEffDecay(const char *name, const char *title,
42 RooAbsReal& tau, RooAbsReal& dm,
43 RooAbsReal& avgMistag, RooAbsReal& CPeigenval,
45 RooAbsReal& effRatio, RooAbsReal& delMistag,
46 const RooResolutionModel& model, DecayType type) :
47 RooAbsAnaConvPdf(name,title,model,t),
48 _absLambda("absLambda","Absolute value of lambda",this,a),
49 _argLambda("argLambda","Arg(Lambda)",this,b),
50 _effRatio("effRatio","B0/B0bar efficiency ratio",this,effRatio),
51 _CPeigenval("CPeigenval","CP eigen value",this,CPeigenval),
52 _avgMistag("avgMistag","Average mistag rate",this,avgMistag),
53 _delMistag("delMistag","Delta mistag rate",this,delMistag),
54 _t("t","time",this,t),
55 _tau("tau","decay time",this,tau),
56 _dm("dm","mixing frequency",this,dm),
57 _tag("tag","CP state",this,tag),
58 _genB0Frac(0),
59 _type(type)
60{
61 switch(type) {
62 case SingleSided:
63 _basisExp = declareBasis("exp(-@0/@1)",RooArgList(tau,dm)) ;
64 _basisSin = declareBasis("exp(-@0/@1)*sin(@0*@2)",RooArgList(tau,dm)) ;
65 _basisCos = declareBasis("exp(-@0/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
66 break ;
67 case Flipped:
68 _basisExp = declareBasis("exp(@0)/@1)",RooArgList(tau,dm)) ;
69 _basisSin = declareBasis("exp(@0/@1)*sin(@0*@2)",RooArgList(tau,dm)) ;
70 _basisCos = declareBasis("exp(@0/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
71 break ;
72 case DoubleSided:
73 _basisExp = declareBasis("exp(-abs(@0)/@1)",RooArgList(tau,dm)) ;
74 _basisSin = declareBasis("exp(-abs(@0)/@1)*sin(@0*@2)",RooArgList(tau,dm)) ;
75 _basisCos = declareBasis("exp(-abs(@0)/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
76 break ;
77 }
78}
79
80////////////////////////////////////////////////////////////////////////////////
81/// Copy constructor.
82
85 _absLambda("absLambda",this,other._absLambda),
86 _argLambda("argLambda",this,other._argLambda),
87 _effRatio("effRatio",this,other._effRatio),
88 _CPeigenval("CPeigenval",this,other._CPeigenval),
89 _avgMistag("avgMistag",this,other._avgMistag),
90 _delMistag("delMistag",this,other._delMistag),
91 _t("t",this,other._t),
92 _tau("tau",this,other._tau),
93 _dm("dm",this,other._dm),
94 _tag("tag",this,other._tag),
95 _genB0Frac(other._genB0Frac),
96 _type(other._type),
97 _basisExp(other._basisExp),
98 _basisSin(other._basisSin),
99 _basisCos(other._basisCos)
100{
101}
102
103////////////////////////////////////////////////////////////////////////////////
104/// Destructor.
105
107{
108}
109
110////////////////////////////////////////////////////////////////////////////////
111/// B0 : _tag = +1
112///
113/// B0bar : _tag = -1
114/// \param[in] basisIndex
115
117{
118 if (basisIndex==_basisExp) {
119 //exp term: (1 -/+ dw)(1+a^2)/2
120 return (1 - _tag*_delMistag)*(1+_absLambda*_absLambda)/2 ;
121 // = 1 + _tag*deltaDil/2
122 }
123
124 if (basisIndex==_basisSin) {
125 //sin term: +/- (1-2w)*ImLambda(= -etaCP * |l| * sin2b)
127 // = _tag*avgDil * ...
128 }
129
130 if (basisIndex==_basisCos) {
131 //cos term: +/- (1-2w)*(1-a^2)/2
132 return -1*_tag*(1-2*_avgMistag)*(1-_absLambda*_absLambda)/2 ;
133 // = -_tag*avgDil * ...
134 }
135
136 return 0 ;
137}
138
139////////////////////////////////////////////////////////////////////////////////
140
141Int_t RooBCPEffDecay::getCoefAnalyticalIntegral(Int_t /*code*/, RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
142{
143 if (rangeName) return 0 ;
144
145 if (matchArgs(allVars,analVars,_tag)) return 1 ;
146 return 0 ;
147}
148
149////////////////////////////////////////////////////////////////////////////////
150
151Double_t RooBCPEffDecay::coefAnalyticalIntegral(Int_t basisIndex, Int_t code, const char* /*rangeName*/) const
152{
153 switch(code) {
154 // No integration
155 case 0: return coefficient(basisIndex) ;
156
157 // Integration over 'tag'
158 case 1:
159 if (basisIndex==_basisExp) {
160 return (1+_absLambda*_absLambda) ;
161 }
162
163 if (basisIndex==_basisSin || basisIndex==_basisCos) {
164 return 0 ;
165 }
166 break ;
167
168 default:
169 assert(0) ;
170 }
171
172 return 0 ;
173}
174
175////////////////////////////////////////////////////////////////////////////////
176
177Int_t RooBCPEffDecay::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t staticInitOK) const
178{
179 if (staticInitOK) {
180 if (matchArgs(directVars,generateVars,_t,_tag)) return 2 ;
181 }
182 if (matchArgs(directVars,generateVars,_t)) return 1 ;
183 return 0 ;
184}
185
186////////////////////////////////////////////////////////////////////////////////
187
189{
190 if (code==2) {
191 // Calculate the fraction of mixed events to generate
192 Double_t sumInt = RooRealIntegral("sumInt","sum integral",*this,RooArgSet(_t.arg(),_tag.arg())).getVal() ;
193 _tag = 1 ;
194 Double_t b0Int = RooRealIntegral("mixInt","mix integral",*this,RooArgSet(_t.arg())).getVal() ;
195 _genB0Frac = b0Int/sumInt ;
196 }
197}
198
199////////////////////////////////////////////////////////////////////////////////
200/// Generates mix-state dependent.
201/// \param[in] code
202
204{
205 if (code==2) {
207 _tag = (rand<=_genB0Frac) ? 1 : -1 ;
208 }
209
210 // Generate delta-t dependent
211 while(1) {
213 Double_t tval(0) ;
214
215 switch(_type) {
216 case SingleSided:
217 tval = -_tau*log(rand);
218 break ;
219 case Flipped:
220 tval= +_tau*log(rand);
221 break ;
222 case DoubleSided:
223 tval = (rand<=0.5) ? -_tau*log(2*rand) : +_tau*log(2*(rand-0.5)) ;
224 break ;
225 }
226
227 // Accept event if T is in generated range
228 Double_t maxDil = 1.0 ;
230 Double_t maxAcceptProb = (1+al2) + fabs(maxDil*_CPeigenval*_absLambda*_argLambda) + fabs(maxDil*(1-al2)/2);
231 Double_t acceptProb = (1+al2)/2*(1-_tag*_delMistag)
232 - (_tag*(1-2*_avgMistag))*(_CPeigenval*_absLambda*_argLambda)*sin(_dm*tval)
233 - (_tag*(1-2*_avgMistag))*(1-al2)/2*cos(_dm*tval);
234
235 Bool_t accept = maxAcceptProb*RooRandom::uniform() < acceptProb ? kTRUE : kFALSE ;
236
237 if (tval<_t.max() && tval>_t.min() && accept) {
238 _t = tval ;
239 break ;
240 }
241 }
242
243}
#define b(i)
Definition: RSha256.hxx:100
int Int_t
Definition: RtypesCore.h:45
const Bool_t kFALSE
Definition: RtypesCore.h:101
bool Bool_t
Definition: RtypesCore.h:63
double Double_t
Definition: RtypesCore.h:59
const Bool_t kTRUE
Definition: RtypesCore.h:100
#define ClassImp(name)
Definition: Rtypes.h:364
char name[80]
Definition: TGX11.cxx:110
int type
Definition: TGX11.cxx:121
RooAbsAnaConvPdf is the base class for PDFs that represent a physics model that can be analytically c...
Int_t declareBasis(const char *expression, const RooArgList &params)
Declare a basis function for use in this physics model.
friend class RooArgSet
Definition: RooAbsArg.h:600
RooAbsCategory is the base class for objects that represent a discrete value with a finite number of ...
friend class RooRealIntegral
Definition: RooAbsPdf.h:345
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:22
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:35
PDF describing decay time distribution of B meson including effects of standard model CP violation.
DecayType _type
Double_t _genB0Frac
virtual Double_t coefAnalyticalIntegral(Int_t coef, Int_t code, const char *rangeName=0) const
Default implementation of function implementing advertised integrals.
RooRealProxy _t
RooRealProxy _avgMistag
virtual Double_t coefficient(Int_t basisIndex) const
B0 : _tag = +1.
virtual Int_t getCoefAnalyticalIntegral(Int_t coef, RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Default implementation of function advertising integration capabilities.
RooRealProxy _delMistag
void initGenerator(Int_t code)
Interface for one-time initialization to setup the generator for the specified code.
RooRealProxy _dm
RooRealProxy _absLambda
void generateEvent(Int_t code)
Generates mix-state dependent.
RooRealProxy _argLambda
virtual ~RooBCPEffDecay()
Destructor.
RooCategoryProxy _tag
RooRealProxy _CPeigenval
Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
RooRealProxy _tau
static Double_t uniform(TRandom *generator=randomGenerator())
Return a number uniformly distributed from (0,1)
Definition: RooRandom.cxx:83
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:39
RooResolutionModel is the base class for PDFs that represent a resolution model that can be convolute...
double min(const char *rname=0) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
double max(const char *rname=0) const
Query upper limit of range. This requires the payload to be RooAbsRealLValue or derived.
const T & arg() const
Return reference to object held in proxy.
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
auto * a
Definition: textangle.C:12