Logo ROOT   6.14/05
Reference Guide
RooBDecay.cxx
Go to the documentation of this file.
1 /*****************************************************************************
2  * Project: RooFit *
3  * Package: RooFitModels *
4  * @(#)root/roofit:$Id$
5  * Authors: *
6  * PL, Parker C Lund, UC Irvine *
7  * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8  * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
9  * *
10  * Copyright (c) 2000-2005, Regents of the University of California *
11  * and Stanford University. All rights reserved. *
12  * *
13  * Redistribution and use in source and binary forms, *
14  * with or without modification, are permitted according to the terms *
15  * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
16  *****************************************************************************/
17 
18 
19 /** \class RooBDecay
20  \ingroup Roofit
21 
22 Most general description of B decay time distribution with effects
23 of CP violation, mixing and life time differences. This function can
24 be analytically convolved with any RooResolutionModel implementation
25 **/
26 
27 #include "RooFit.h"
28 
29 #include "Riostream.h"
30 #include "TMath.h"
31 #include "RooBDecay.h"
32 #include "RooRealVar.h"
33 #include "RooRandom.h"
34 
35 #include "TError.h"
36 
37 using namespace std;
38 
40 
41 ////////////////////////////////////////////////////////////////////////////////
42 
43 RooBDecay::RooBDecay(const char *name, const char* title,
44  RooRealVar& t, RooAbsReal& tau, RooAbsReal& dgamma,
47  RooAbsAnaConvPdf(name, title, model, t),
48  _t("t", "time", this, t),
49  _tau("tau", "Average Decay Time", this, tau),
50  _dgamma("dgamma", "Delta Gamma", this, dgamma),
51  _f0("f0", "Cosh Coefficient", this, f0),
52  _f1("f1", "Sinh Coefficient", this, f1),
53  _f2("f2", "Cos Coefficient", this, f2),
54  _f3("f3", "Sin Coefficient", this, f3),
55  _dm("dm", "Delta Mass", this, dm),
56  _type(type)
57 
58 {
59  //Constructor
60  switch(type)
61  {
62  case SingleSided:
63  _basisCosh = declareBasis("exp(-@0/@1)*cosh(@0*@2/2)", RooArgList(tau,dgamma));
64  _basisSinh = declareBasis("exp(-@0/@1)*sinh(@0*@2/2)", RooArgList(tau,dgamma));
65  _basisCos = declareBasis("exp(-@0/@1)*cos(@0*@2)",RooArgList(tau, dm));
66  _basisSin = declareBasis("exp(-@0/@1)*sin(@0*@2)",RooArgList(tau, dm));
67  break;
68  case Flipped:
69  _basisCosh = declareBasis("exp(@0/@1)*cosh(@0*@2/2)", RooArgList(tau,dgamma));
70  _basisSinh = declareBasis("exp(@0/@1)*sinh(@0*@2/2)", RooArgList(tau,dgamma));
71  _basisCos = declareBasis("exp(@0/@1)*cos(@0*@2)",RooArgList(tau, dm));
72  _basisSin = declareBasis("exp(@0/@1)*sin(@0*@2)",RooArgList(tau, dm));
73  break;
74  case DoubleSided:
75  _basisCosh = declareBasis("exp(-abs(@0)/@1)*cosh(@0*@2/2)", RooArgList(tau,dgamma));
76  _basisSinh = declareBasis("exp(-abs(@0)/@1)*sinh(@0*@2/2)", RooArgList(tau,dgamma));
77  _basisCos = declareBasis("exp(-abs(@0)/@1)*cos(@0*@2)",RooArgList(tau, dm));
78  _basisSin = declareBasis("exp(-abs(@0)/@1)*sin(@0*@2)",RooArgList(tau, dm));
79  break;
80  }
81 }
82 
83 ////////////////////////////////////////////////////////////////////////////////
84 ///Copy constructor
85 
86 RooBDecay::RooBDecay(const RooBDecay& other, const char* name) :
87  RooAbsAnaConvPdf(other, name),
88  _t("t", this, other._t),
89  _tau("tau", this, other._tau),
90  _dgamma("dgamma", this, other._dgamma),
91  _f0("f0", this, other._f0),
92  _f1("f1", this, other._f1),
93  _f2("f2", this, other._f2),
94  _f3("f3", this, other._f3),
95  _dm("dm", this, other._dm),
96  _basisCosh(other._basisCosh),
97  _basisSinh(other._basisSinh),
98  _basisCos(other._basisCos),
99  _basisSin(other._basisSin),
100  _type(other._type)
101 {
102 }
103 
104 ////////////////////////////////////////////////////////////////////////////////
105 ///Destructor
106 
108 {
109 }
110 
111 ////////////////////////////////////////////////////////////////////////////////
112 
114 {
115  if(basisIndex == _basisCosh)
116  {
117  return _f0;
118  }
119  if(basisIndex == _basisSinh)
120  {
121  return _f1;
122  }
123  if(basisIndex == _basisCos)
124  {
125  return _f2;
126  }
127  if(basisIndex == _basisSin)
128  {
129  return _f3;
130  }
131 
132  return 0 ;
133 }
134 
135 ////////////////////////////////////////////////////////////////////////////////
136 
138 {
139  if(basisIndex == _basisCosh)
140  {
141  return _f0.arg().getVariables();
142  }
143  if(basisIndex == _basisSinh)
144  {
145  return _f1.arg().getVariables();
146  }
147  if(basisIndex == _basisCos)
148  {
149  return _f2.arg().getVariables();
150  }
151  if(basisIndex == _basisSin)
152  {
153  return _f3.arg().getVariables();
154  }
155 
156  return 0 ;
157 }
158 
159 ////////////////////////////////////////////////////////////////////////////////
160 
161 Int_t RooBDecay::getCoefAnalyticalIntegral(Int_t coef, RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
162 {
163  if(coef == _basisCosh)
164  {
165  return _f0.arg().getAnalyticalIntegral(allVars,analVars,rangeName) ;
166  }
167  if(coef == _basisSinh)
168  {
169  return _f1.arg().getAnalyticalIntegral(allVars,analVars,rangeName) ;
170  }
171  if(coef == _basisCos)
172  {
173  return _f2.arg().getAnalyticalIntegral(allVars,analVars,rangeName) ;
174  }
175  if(coef == _basisSin)
176  {
177  return _f3.arg().getAnalyticalIntegral(allVars,analVars,rangeName) ;
178  }
179 
180  return 0 ;
181 }
182 
183 ////////////////////////////////////////////////////////////////////////////////
184 
185 Double_t RooBDecay::coefAnalyticalIntegral(Int_t coef, Int_t code, const char* rangeName) const
186 {
187  if(coef == _basisCosh)
188  {
189  return _f0.arg().analyticalIntegral(code,rangeName) ;
190  }
191  if(coef == _basisSinh)
192  {
193  return _f1.arg().analyticalIntegral(code,rangeName) ;
194  }
195  if(coef == _basisCos)
196  {
197  return _f2.arg().analyticalIntegral(code,rangeName) ;
198  }
199  if(coef == _basisSin)
200  {
201  return _f3.arg().analyticalIntegral(code,rangeName) ;
202  }
203 
204  return 0 ;
205 }
206 
207 ////////////////////////////////////////////////////////////////////////////////
208 
209 Int_t RooBDecay::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t /*staticInitOK*/) const
210 {
211  if (matchArgs(directVars, generateVars, _t)) return 1;
212  return 0;
213 }
214 
215 ////////////////////////////////////////////////////////////////////////////////
216 
218 {
219  R__ASSERT(code==1);
220  Double_t gammamin = 1/_tau-TMath::Abs(_dgamma)/2;
221  while(1) {
222  Double_t t = -log(RooRandom::uniform())/gammamin;
223  if (_type == Flipped || (_type == DoubleSided && RooRandom::uniform() <0.5) ) t *= -1;
224  if ( t<_t.min() || t>_t.max() ) continue;
225 
226  Double_t dgt = _dgamma*t/2;
227  Double_t dmt = _dm*t;
228  Double_t ft = fabs(t);
229  Double_t f = exp(-ft/_tau)*(_f0*cosh(dgt)+_f1*sinh(dgt)+_f2*cos(dmt)+_f3*sin(dmt));
230  if(f < 0) {
231  cout << "RooBDecay::generateEvent(" << GetName() << ") ERROR: PDF value less than zero" << endl;
232  ::abort();
233  }
234  Double_t w = 1.001*exp(-ft*gammamin)*(TMath::Abs(_f0)+TMath::Abs(_f1)+sqrt(_f2*_f2+_f3*_f3));
235  if(w < f) {
236  cout << "RooBDecay::generateEvent(" << GetName() << ") ERROR: Envelope function less than p.d.f. " << endl;
237  cout << _f0 << endl;
238  cout << _f1 << endl;
239  cout << _f2 << endl;
240  cout << _f3 << endl;
241  ::abort();
242  }
243  if(w*RooRandom::uniform() > f) continue;
244  _t = t;
245  break;
246  }
247 }
RooArgSet * getVariables(Bool_t stripDisconnected=kTRUE) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
Definition: RooAbsArg.cxx:2073
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
RooArgSet * coefVars(Int_t coefIdx) const
Return set of parameters with are used exclusively by the coefficient functions.
Definition: RooBDecay.cxx:137
virtual ~RooBDecay()
Destructor.
Definition: RooBDecay.cxx:107
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...
Definition: RooBDecay.cxx:209
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
virtual Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported...
Definition: RooAbsReal.cxx:335
#define R__ASSERT(e)
Definition: TError.h:96
#define f(i)
Definition: RSha256.hxx:104
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
STL namespace.
Short_t Abs(Short_t d)
Definition: TMathBase.h:108
double cos(double)
RooRealProxy _f0
Definition: RooBDecay.h:61
double sqrt(double)
Most general description of B decay time distribution with effects of CP violation, mixing and life time differences.
Definition: RooBDecay.h:24
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t tau
Definition: TRolke.cxx:630
double sinh(double)
RooRealProxy _f2
Definition: RooBDecay.h:63
DecayType _type
Definition: RooBDecay.h:71
RooRealProxy _dm
Definition: RooBDecay.h:65
double sin(double)
RooRealProxy _f1
Definition: RooBDecay.h:62
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
double cosh(double)
void generateEvent(Int_t code)
Interface for generation of anan event using the algorithm corresponding to the specified code...
Definition: RooBDecay.cxx:217
Int_t _basisSinh
Definition: RooBDecay.h:67
RooRealProxy _tau
Definition: RooBDecay.h:59
#define ClassImp(name)
Definition: Rtypes.h:359
static Double_t uniform(TRandom *generator=randomGenerator())
Return a number uniformly distributed from (0,1)
Definition: RooRandom.cxx:84
Double_t min(const char *rname=0) const
Definition: RooRealProxy.h:56
double Double_t
Definition: RtypesCore.h:55
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
RooRealProxy _t
Definition: RooBDecay.h:58
RooRealProxy _f3
Definition: RooBDecay.h:64
int type
Definition: TGX11.cxx:120
Int_t declareBasis(const char *expression, const RooArgList &params)
Declare a basis function for use in this physics model.
Int_t getCoefAnalyticalIntegral(Int_t coef, RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Default implementation of function advertising integration capabilities.
Definition: RooBDecay.cxx:161
Double_t coefAnalyticalIntegral(Int_t coef, Int_t code, const char *rangeName=0) const
Default implementation of function implementing advertised integrals.
Definition: RooBDecay.cxx:185
RooRealProxy _dgamma
Definition: RooBDecay.h:60
Int_t _basisSin
Definition: RooBDecay.h:69
Double_t max(const char *rname=0) const
Definition: RooRealProxy.h:57
Int_t _basisCosh
Definition: RooBDecay.h:66
virtual Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral. ...
Definition: RooAbsReal.cxx:361
TF1 * f1
Definition: legend1.C:11
RooBDecay()
Definition: RooBDecay.h:32
Int_t _basisCos
Definition: RooBDecay.h:68
const RooAbsReal & arg() const
Definition: RooRealProxy.h:43
double exp(double)
virtual Double_t coefficient(Int_t basisIndex) const
Definition: RooBDecay.cxx:113
char name[80]
Definition: TGX11.cxx:109
double log(double)