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 "Riostream.h"
26#include "TMath.h"
27#include "RooRealVar.h"
28#include "RooBMixDecay.h"
29#include "RooRealIntegral.h"
30#include "RooRandom.h"
31
32using namespace std;
33
35
36////////////////////////////////////////////////////////////////////////////////
37/// Constructor
38
39RooBMixDecay::RooBMixDecay(const char *name, const char *title,
40 RooRealVar& t, RooAbsCategory& mixState,
41 RooAbsCategory& tagFlav,
42 RooAbsReal& tau, RooAbsReal& dm,
43 RooAbsReal& mistag, RooAbsReal& delMistag,
44 const RooResolutionModel& model,
46 RooAbsAnaConvPdf(name,title,model,t),
47 _type(type),
48 _mistag("mistag","Mistag rate",this,mistag),
49 _delMistag("delMistag","Delta mistag rate",this,delMistag),
50 _mixState("mixState","Mixing state",this,mixState),
51 _tagFlav("tagFlav","Flavour of tagged B0",this,tagFlav),
52 _tau("tau","Mixing life time",this,tau),
53 _dm("dm","Mixing frequency",this,dm),
54 _t("_t","time",this,t), _genMixFrac(0)
55{
56 switch(type) {
57 case SingleSided:
58 _basisExp = declareBasis("exp(-@0/@1)",RooArgList(tau)) ;
59 _basisCos = declareBasis("exp(-@0/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
60 break ;
61 case Flipped:
62 _basisExp = declareBasis("exp(@0/@1)",RooArgList(tau)) ;
63 _basisCos = declareBasis("exp(@0/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
64 break ;
65 case DoubleSided:
66 _basisExp = declareBasis("exp(-abs(@0)/@1)",RooArgList(tau)) ;
67 _basisCos = declareBasis("exp(-abs(@0)/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
68 break ;
69 }
70}
71
72////////////////////////////////////////////////////////////////////////////////
73/// Copy constructor
74
75RooBMixDecay::RooBMixDecay(const RooBMixDecay& other, const char* name) :
77 _type(other._type),
78 _mistag("mistag",this,other._mistag),
79 _delMistag("delMistag",this,other._delMistag),
80 _mixState("mixState",this,other._mixState),
81 _tagFlav("tagFlav",this,other._tagFlav),
82 _tau("tau",this,other._tau),
83 _dm("dm",this,other._dm),
84 _t("t",this,other._t),
85 _basisExp(other._basisExp),
86 _basisCos(other._basisCos),
87 _genMixFrac(other._genMixFrac),
88 _genFlavFrac(other._genFlavFrac),
89 _genFlavFracMix(other._genFlavFracMix),
90 _genFlavFracUnmix(other._genFlavFracUnmix)
91{
92}
93
94////////////////////////////////////////////////////////////////////////////////
95/// Destructor
96
98{
99}
100
101////////////////////////////////////////////////////////////////////////////////
102/// Comp with tFit MC: must be (1 - tagFlav*...)
103
104double RooBMixDecay::coefficient(Int_t basisIndex) const
105{
106 if (basisIndex==_basisExp) {
107 return (1 - _tagFlav*_delMistag) ;
108 }
109
110 if (basisIndex==_basisCos) {
111 return _mixState*(1-2*_mistag) ;
112 }
113
114 return 0 ;
115}
116
117////////////////////////////////////////////////////////////////////////////////
118/// cout << "RooBMixDecay::getCoefAI " ; allVars.Print("1") ;
119
120Int_t RooBMixDecay::getCoefAnalyticalIntegral(Int_t /*code*/, RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
121{
122 if (rangeName) {
123 return 0 ;
124 }
125
126 if (matchArgs(allVars,analVars,_mixState,_tagFlav)) return 3 ;
127 if (matchArgs(allVars,analVars,_mixState)) return 2 ;
128 if (matchArgs(allVars,analVars,_tagFlav)) return 1 ;
129 return 0 ;
130}
131
132////////////////////////////////////////////////////////////////////////////////
133
134double RooBMixDecay::coefAnalyticalIntegral(Int_t basisIndex, Int_t code, const char* /*rangeName*/) const
135{
136 switch(code) {
137 // No integration
138 case 0: return coefficient(basisIndex) ;
139
140 // Integration over 'mixState' and 'tagFlav'
141 case 3:
142 if (basisIndex==_basisExp) {
143 return 4.0 ;
144 }
145 if (basisIndex==_basisCos) {
146 return 0.0 ;
147 }
148 break ;
149
150 // Integration over 'mixState'
151 case 2:
152 if (basisIndex==_basisExp) {
153 return 2.0*coefficient(basisIndex) ;
154 }
155 if (basisIndex==_basisCos) {
156 return 0.0 ;
157 }
158 break ;
159
160 // Integration over 'tagFlav'
161 case 1:
162 if (basisIndex==_basisExp) {
163 return 2.0 ;
164 }
165 if (basisIndex==_basisCos) {
166 return 2.0*coefficient(basisIndex) ;
167 }
168 break ;
169
170 default:
171 assert(0) ;
172 }
173
174 return 0 ;
175}
176
177////////////////////////////////////////////////////////////////////////////////
178
179Int_t RooBMixDecay::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, bool staticInitOK) const
180{
181 if (staticInitOK) {
182 if (matchArgs(directVars,generateVars,_t,_mixState,_tagFlav)) return 4 ;
183 if (matchArgs(directVars,generateVars,_t,_mixState)) return 3 ;
184 if (matchArgs(directVars,generateVars,_t,_tagFlav)) return 2 ;
185 }
186
187 if (matchArgs(directVars,generateVars,_t)) return 1 ;
188 return 0 ;
189}
190
191////////////////////////////////////////////////////////////////////////////////
192
194{
195 switch (code) {
196 case 2:
197 {
198 // Calculate the fraction of B0bar events to generate
199 double sumInt = RooRealIntegral("sumInt","sum integral",*this,RooArgSet(_t.arg(),_tagFlav.arg())).getVal() ;
200 _tagFlav = 1 ; // B0
201 double flavInt = RooRealIntegral("flavInt","flav integral",*this,RooArgSet(_t.arg())).getVal() ;
202 _genFlavFrac = flavInt/sumInt ;
203 break ;
204 }
205 case 3:
206 {
207 // Calculate the fraction of mixed events to generate
208 double sumInt = RooRealIntegral("sumInt","sum integral",*this,RooArgSet(_t.arg(),_mixState.arg())).getVal() ;
209 _mixState = -1 ; // mixed
210 double mixInt = RooRealIntegral("mixInt","mix integral",*this,RooArgSet(_t.arg())).getVal() ;
211 _genMixFrac = mixInt/sumInt ;
212 break ;
213 }
214 case 4:
215 {
216 // Calculate the fraction of mixed events to generate
217 double sumInt = RooRealIntegral("sumInt","sum integral",*this,RooArgSet(_t.arg(),_mixState.arg(),_tagFlav.arg())).getVal() ;
218 _mixState = -1 ; // mixed
219 double mixInt = RooRealIntegral("mixInt","mix integral",*this,RooArgSet(_t.arg(),_tagFlav.arg())).getVal() ;
220 _genMixFrac = mixInt/sumInt ;
221
222 // Calculate the fraction of B0bar tags for mixed and unmixed
223 RooRealIntegral dtInt("mixInt","mix integral",*this,RooArgSet(_t.arg())) ;
224 _mixState = -1 ; // Mixed
225 _tagFlav = 1 ; // B0
226 _genFlavFracMix = dtInt.getVal() / mixInt ;
227 _mixState = 1 ; // Unmixed
228 _tagFlav = 1 ; // B0
229 _genFlavFracUnmix = dtInt.getVal() / (sumInt - mixInt) ;
230 break ;
231 }
232 }
233}
234
235////////////////////////////////////////////////////////////////////////////////
236/// Generate mix-state dependent
237
239{
240 switch(code) {
241 case 2:
242 {
243 double rand = RooRandom::uniform() ;
244 _tagFlav = (Int_t) ((rand<=_genFlavFrac) ? 1 : -1) ;
245 break ;
246 }
247 case 3:
248 {
249 double rand = RooRandom::uniform() ;
250 _mixState = (Int_t) ((rand<=_genMixFrac) ? -1 : 1) ;
251 break ;
252 }
253 case 4:
254 {
255 double rand = RooRandom::uniform() ;
256 _mixState = (Int_t) ((rand<=_genMixFrac) ? -1 : 1) ;
257
258 rand = RooRandom::uniform() ;
259 double genFlavFrac = (_mixState==-1) ? _genFlavFracMix : _genFlavFracUnmix ;
260 _tagFlav = (Int_t) ((rand<=genFlavFrac) ? 1 : -1) ;
261 break ;
262 }
263 }
264
265 // Generate delta-t dependent
266 while(1) {
267 double rand = RooRandom::uniform() ;
268 double tval(0) ;
269
270 switch(_type) {
271 case SingleSided:
272 tval = -_tau*log(rand);
273 break ;
274 case Flipped:
275 tval= +_tau*log(rand);
276 break ;
277 case DoubleSided:
278 tval = (rand<=0.5) ? -_tau*log(2*rand) : +_tau*log(2*(rand-0.5)) ;
279 break ;
280 }
281
282 // Accept event if T is in generated range
283 double dil = 1-2.*_mistag ;
284 double maxAcceptProb = 1 + TMath::Abs(_delMistag) + TMath::Abs(dil) ;
285 double acceptProb = (1-_tagFlav*_delMistag) + _mixState*dil*cos(_dm*tval);
286 bool mixAccept = maxAcceptProb*RooRandom::uniform() < acceptProb ? true : false ;
287
288 if (tval<_t.max() && tval>_t.min() && mixAccept) {
289 _t = tval ;
290 break ;
291 }
292 }
293}
int Int_t
Definition: RtypesCore.h:45
#define ClassImp(name)
Definition: Rtypes.h:375
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
char name[80]
Definition: TGX11.cxx:110
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 RooRealIntegral
Definition: RooAbsArg.h:644
A space to attach TBranches.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:62
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
bool 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:56
Class RooBMixDecay is a RooAbsAnaConvPdf implementation that describes the decay of B mesons with the...
Definition: RooBMixDecay.h:23
Int_t _basisExp
Definition: RooBMixDecay.h:59
Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, bool staticInitOK=true) const override
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
RooRealProxy _tau
Definition: RooBMixDecay.h:56
double _genMixFrac
Definition: RooBMixDecay.h:62
double coefAnalyticalIntegral(Int_t coef, Int_t code, const char *rangeName=nullptr) const override
Default implementation of function implementing advertised integrals.
Int_t _basisCos
Definition: RooBMixDecay.h:60
DecayType _type
Definition: RooBMixDecay.h:51
double _genFlavFracMix
Definition: RooBMixDecay.h:64
RooRealProxy _mistag
Definition: RooBMixDecay.h:52
RooRealProxy _t
Definition: RooBMixDecay.h:58
void initGenerator(Int_t code) override
Interface for one-time initialization to setup the generator for the specified code.
void generateEvent(Int_t code) override
Generate mix-state dependent.
RooCategoryProxy _mixState
Definition: RooBMixDecay.h:54
Int_t getCoefAnalyticalIntegral(Int_t coef, RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=nullptr) const override
cout << "RooBMixDecay::getCoefAI " ; allVars.Print("1") ;
double coefficient(Int_t basisIndex) const override
Comp with tFit MC: must be (1 - tagFlav*...)
double _genFlavFracUnmix
Definition: RooBMixDecay.h:65
RooRealProxy _dm
Definition: RooBMixDecay.h:57
~RooBMixDecay() override
Destructor.
double _genFlavFrac
do not persist
Definition: RooBMixDecay.h:63
RooCategoryProxy _tagFlav
Definition: RooBMixDecay.h:55
RooRealProxy _delMistag
Definition: RooBMixDecay.h:53
static double uniform(TRandom *generator=randomGenerator())
Return a number uniformly distributed from (0,1)
Definition: RooRandom.cxx:81
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
RooResolutionModel is the base class for PDFs that represent a resolution model that can be convolute...
double max(const char *rname=nullptr) 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.
double min(const char *rname=nullptr) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
RVec< PromoteType< T > > cos(const RVec< T > &v)
Definition: RVec.hxx:1776
RVec< PromoteType< T > > log(const RVec< T > &v)
Definition: RVec.hxx:1765
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition: TMathBase.h:123