Logo ROOT  
Reference Guide
RooBMixDecay.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 RooBMixDecay
18 \ingroup Roofit
19
20Class RooBMixDecay is a RooAbsAnaConvPdf implementation that describes
21the decay of B mesons with the effects of B0/B0bar mixing.
22This function can be analytically convolved with any RooResolutionModel implementation
23**/
24
25#include "RooFit.h"
26
27#include "Riostream.h"
28#include "TMath.h"
29#include "RooRealVar.h"
30#include "RooBMixDecay.h"
31#include "RooRealIntegral.h"
32#include "RooRandom.h"
33
34using namespace std;
35
37
38////////////////////////////////////////////////////////////////////////////////
39/// Constructor
40
41RooBMixDecay::RooBMixDecay(const char *name, const char *title,
42 RooRealVar& t, RooAbsCategory& mixState,
43 RooAbsCategory& tagFlav,
44 RooAbsReal& tau, RooAbsReal& dm,
45 RooAbsReal& mistag, RooAbsReal& delMistag,
46 const RooResolutionModel& model,
48 RooAbsAnaConvPdf(name,title,model,t),
49 _type(type),
50 _mistag("mistag","Mistag rate",this,mistag),
51 _delMistag("delMistag","Delta mistag rate",this,delMistag),
52 _mixState("mixState","Mixing state",this,mixState),
53 _tagFlav("tagFlav","Flavour of tagged B0",this,tagFlav),
54 _tau("tau","Mixing life time",this,tau),
55 _dm("dm","Mixing frequency",this,dm),
56 _t("_t","time",this,t), _genMixFrac(0)
57{
58 switch(type) {
59 case SingleSided:
60 _basisExp = declareBasis("exp(-@0/@1)",RooArgList(tau)) ;
61 _basisCos = declareBasis("exp(-@0/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
62 break ;
63 case Flipped:
64 _basisExp = declareBasis("exp(@0/@1)",RooArgList(tau)) ;
65 _basisCos = declareBasis("exp(@0/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
66 break ;
67 case DoubleSided:
68 _basisExp = declareBasis("exp(-abs(@0)/@1)",RooArgList(tau)) ;
69 _basisCos = declareBasis("exp(-abs(@0)/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
70 break ;
71 }
72}
73
74////////////////////////////////////////////////////////////////////////////////
75/// Copy constructor
76
77RooBMixDecay::RooBMixDecay(const RooBMixDecay& other, const char* name) :
79 _type(other._type),
80 _mistag("mistag",this,other._mistag),
81 _delMistag("delMistag",this,other._delMistag),
82 _mixState("mixState",this,other._mixState),
83 _tagFlav("tagFlav",this,other._tagFlav),
84 _tau("tau",this,other._tau),
85 _dm("dm",this,other._dm),
86 _t("t",this,other._t),
87 _basisExp(other._basisExp),
88 _basisCos(other._basisCos),
89 _genMixFrac(other._genMixFrac),
90 _genFlavFrac(other._genFlavFrac),
91 _genFlavFracMix(other._genFlavFracMix),
92 _genFlavFracUnmix(other._genFlavFracUnmix)
93{
94}
95
96////////////////////////////////////////////////////////////////////////////////
97/// Destructor
98
100{
101}
102
103////////////////////////////////////////////////////////////////////////////////
104/// Comp with tFit MC: must be (1 - tagFlav*...)
105
107{
108 if (basisIndex==_basisExp) {
109 return (1 - _tagFlav*_delMistag) ;
110 }
111
112 if (basisIndex==_basisCos) {
113 return _mixState*(1-2*_mistag) ;
114 }
115
116 return 0 ;
117}
118
119////////////////////////////////////////////////////////////////////////////////
120/// cout << "RooBMixDecay::getCoefAI " ; allVars.Print("1") ;
121
122Int_t RooBMixDecay::getCoefAnalyticalIntegral(Int_t /*code*/, RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
123{
124 if (rangeName) {
125 return 0 ;
126 }
127
128 if (matchArgs(allVars,analVars,_mixState,_tagFlav)) return 3 ;
129 if (matchArgs(allVars,analVars,_mixState)) return 2 ;
130 if (matchArgs(allVars,analVars,_tagFlav)) return 1 ;
131 return 0 ;
132}
133
134////////////////////////////////////////////////////////////////////////////////
135
136Double_t RooBMixDecay::coefAnalyticalIntegral(Int_t basisIndex, Int_t code, const char* /*rangeName*/) const
137{
138 switch(code) {
139 // No integration
140 case 0: return coefficient(basisIndex) ;
141
142 // Integration over 'mixState' and 'tagFlav'
143 case 3:
144 if (basisIndex==_basisExp) {
145 return 4.0 ;
146 }
147 if (basisIndex==_basisCos) {
148 return 0.0 ;
149 }
150 break ;
151
152 // Integration over 'mixState'
153 case 2:
154 if (basisIndex==_basisExp) {
155 return 2.0*coefficient(basisIndex) ;
156 }
157 if (basisIndex==_basisCos) {
158 return 0.0 ;
159 }
160 break ;
161
162 // Integration over 'tagFlav'
163 case 1:
164 if (basisIndex==_basisExp) {
165 return 2.0 ;
166 }
167 if (basisIndex==_basisCos) {
168 return 2.0*coefficient(basisIndex) ;
169 }
170 break ;
171
172 default:
173 assert(0) ;
174 }
175
176 return 0 ;
177}
178
179////////////////////////////////////////////////////////////////////////////////
180
181Int_t RooBMixDecay::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t staticInitOK) const
182{
183 if (staticInitOK) {
184 if (matchArgs(directVars,generateVars,_t,_mixState,_tagFlav)) return 4 ;
185 if (matchArgs(directVars,generateVars,_t,_mixState)) return 3 ;
186 if (matchArgs(directVars,generateVars,_t,_tagFlav)) return 2 ;
187 }
188
189 if (matchArgs(directVars,generateVars,_t)) return 1 ;
190 return 0 ;
191}
192
193////////////////////////////////////////////////////////////////////////////////
194
196{
197 switch (code) {
198 case 2:
199 {
200 // Calculate the fraction of B0bar events to generate
201 Double_t sumInt = RooRealIntegral("sumInt","sum integral",*this,RooArgSet(_t.arg(),_tagFlav.arg())).getVal() ;
202 _tagFlav = 1 ; // B0
203 Double_t flavInt = RooRealIntegral("flavInt","flav integral",*this,RooArgSet(_t.arg())).getVal() ;
204 _genFlavFrac = flavInt/sumInt ;
205 break ;
206 }
207 case 3:
208 {
209 // Calculate the fraction of mixed events to generate
210 Double_t sumInt = RooRealIntegral("sumInt","sum integral",*this,RooArgSet(_t.arg(),_mixState.arg())).getVal() ;
211 _mixState = -1 ; // mixed
212 Double_t mixInt = RooRealIntegral("mixInt","mix integral",*this,RooArgSet(_t.arg())).getVal() ;
213 _genMixFrac = mixInt/sumInt ;
214 break ;
215 }
216 case 4:
217 {
218 // Calculate the fraction of mixed events to generate
219 Double_t sumInt = RooRealIntegral("sumInt","sum integral",*this,RooArgSet(_t.arg(),_mixState.arg(),_tagFlav.arg())).getVal() ;
220 _mixState = -1 ; // mixed
221 Double_t mixInt = RooRealIntegral("mixInt","mix integral",*this,RooArgSet(_t.arg(),_tagFlav.arg())).getVal() ;
222 _genMixFrac = mixInt/sumInt ;
223
224 // Calculate the fraction of B0bar tags for mixed and unmixed
225 RooRealIntegral dtInt("mixInt","mix integral",*this,RooArgSet(_t.arg())) ;
226 _mixState = -1 ; // Mixed
227 _tagFlav = 1 ; // B0
228 _genFlavFracMix = dtInt.getVal() / mixInt ;
229 _mixState = 1 ; // Unmixed
230 _tagFlav = 1 ; // B0
231 _genFlavFracUnmix = dtInt.getVal() / (sumInt - mixInt) ;
232 break ;
233 }
234 }
235}
236
237////////////////////////////////////////////////////////////////////////////////
238/// Generate mix-state dependent
239
241{
242 switch(code) {
243 case 2:
244 {
246 _tagFlav = (Int_t) ((rand<=_genFlavFrac) ? 1 : -1) ;
247 break ;
248 }
249 case 3:
250 {
252 _mixState = (Int_t) ((rand<=_genMixFrac) ? -1 : 1) ;
253 break ;
254 }
255 case 4:
256 {
258 _mixState = (Int_t) ((rand<=_genMixFrac) ? -1 : 1) ;
259
260 rand = RooRandom::uniform() ;
261 Double_t genFlavFrac = (_mixState==-1) ? _genFlavFracMix : _genFlavFracUnmix ;
262 _tagFlav = (Int_t) ((rand<=genFlavFrac) ? 1 : -1) ;
263 break ;
264 }
265 }
266
267 // Generate delta-t dependent
268 while(1) {
270 Double_t tval(0) ;
271
272 switch(_type) {
273 case SingleSided:
274 tval = -_tau*log(rand);
275 break ;
276 case Flipped:
277 tval= +_tau*log(rand);
278 break ;
279 case DoubleSided:
280 tval = (rand<=0.5) ? -_tau*log(2*rand) : +_tau*log(2*(rand-0.5)) ;
281 break ;
282 }
283
284 // Accept event if T is in generated range
285 Double_t dil = 1-2.*_mistag ;
286 Double_t maxAcceptProb = 1 + TMath::Abs(_delMistag) + TMath::Abs(dil) ;
287 Double_t acceptProb = (1-_tagFlav*_delMistag) + _mixState*dil*cos(_dm*tval);
288 Bool_t mixAccept = maxAcceptProb*RooRandom::uniform() < acceptProb ? kTRUE : kFALSE ;
289
290 if (tval<_t.max() && tval>_t.min() && mixAccept) {
291 _t = tval ;
292 break ;
293 }
294 }
295}
int Int_t
Definition: RtypesCore.h:45
const Bool_t kFALSE
Definition: RtypesCore.h:101
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:624
RooAbsCategory is the base class for objects that represent a discrete value with a finite number of ...
friend class RooRealIntegral
Definition: RooAbsPdf.h:357
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:63
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:93
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
Class RooBMixDecay is a RooAbsAnaConvPdf implementation that describes the decay of B mesons with the...
Definition: RooBMixDecay.h:23
Double_t _genFlavFrac
do not persist
Definition: RooBMixDecay.h:63
Int_t _basisExp
Definition: RooBMixDecay.h:59
RooRealProxy _tau
Definition: RooBMixDecay.h:56
Int_t _basisCos
Definition: RooBMixDecay.h:60
void generateEvent(Int_t code)
Generate mix-state dependent.
virtual Double_t coefficient(Int_t basisIndex) const
Comp with tFit MC: must be (1 - tagFlav*...)
DecayType _type
Definition: RooBMixDecay.h:51
Double_t _genFlavFracUnmix
Definition: RooBMixDecay.h:65
void initGenerator(Int_t code)
Interface for one-time initialization to setup the generator for the specified code.
RooRealProxy _mistag
Definition: RooBMixDecay.h:52
RooRealProxy _t
Definition: RooBMixDecay.h:58
Double_t _genFlavFracMix
Definition: RooBMixDecay.h:64
RooCategoryProxy _mixState
Definition: RooBMixDecay.h:54
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...
Double_t _genMixFrac
Definition: RooBMixDecay.h:62
RooRealProxy _dm
Definition: RooBMixDecay.h:57
virtual Double_t coefAnalyticalIntegral(Int_t coef, Int_t code, const char *rangeName=0) const
Default implementation of function implementing advertised integrals.
RooCategoryProxy _tagFlav
Definition: RooBMixDecay.h:55
RooRealProxy _delMistag
Definition: RooBMixDecay.h:53
virtual ~RooBMixDecay()
Destructor.
virtual Int_t getCoefAnalyticalIntegral(Int_t coef, RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
cout << "RooBMixDecay::getCoefAI " ; allVars.Print("1") ;
static Double_t uniform(TRandom *generator=randomGenerator())
Return a number uniformly distributed from (0,1)
Definition: RooRandom.cxx:83
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
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.
Short_t Abs(Short_t d)
Definition: TMathBase.h:120