Logo ROOT  
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
22Most general description of B decay time distribution with effects
23of CP violation, mixing and life time differences. This function can
24be analytically convolved with any RooResolutionModel implementation
25**/
26
27#include "Riostream.h"
28#include "TMath.h"
29#include "RooBDecay.h"
30#include "RooRealVar.h"
31#include "RooRandom.h"
32
33#include "TError.h"
34
35using namespace std;
36
38
39////////////////////////////////////////////////////////////////////////////////
40
41RooBDecay::RooBDecay(const char *name, const char* title,
42 RooRealVar& t, RooAbsReal& tau, RooAbsReal& dgamma,
44 RooAbsReal& dm, const RooResolutionModel& model, DecayType type) :
45 RooAbsAnaConvPdf(name, title, model, t),
46 _t("t", "time", this, t),
47 _tau("tau", "Average Decay Time", this, tau),
48 _dgamma("dgamma", "Delta Gamma", this, dgamma),
49 _f0("f0", "Cosh Coefficient", this, f0),
50 _f1("f1", "Sinh Coefficient", this, f1),
51 _f2("f2", "Cos Coefficient", this, f2),
52 _f3("f3", "Sin Coefficient", this, f3),
53 _dm("dm", "Delta Mass", this, dm),
54 _type(type)
55
56{
57 //Constructor
58 switch(type)
59 {
60 case SingleSided:
61 _basisCosh = declareBasis("exp(-@0/@1)*cosh(@0*@2/2)", RooArgList(tau,dgamma));
62 _basisSinh = declareBasis("exp(-@0/@1)*sinh(@0*@2/2)", RooArgList(tau,dgamma));
63 _basisCos = declareBasis("exp(-@0/@1)*cos(@0*@2)",RooArgList(tau, dm));
64 _basisSin = declareBasis("exp(-@0/@1)*sin(@0*@2)",RooArgList(tau, dm));
65 break;
66 case Flipped:
67 _basisCosh = declareBasis("exp(@0/@1)*cosh(@0*@2/2)", RooArgList(tau,dgamma));
68 _basisSinh = declareBasis("exp(@0/@1)*sinh(@0*@2/2)", RooArgList(tau,dgamma));
69 _basisCos = declareBasis("exp(@0/@1)*cos(@0*@2)",RooArgList(tau, dm));
70 _basisSin = declareBasis("exp(@0/@1)*sin(@0*@2)",RooArgList(tau, dm));
71 break;
72 case DoubleSided:
73 _basisCosh = declareBasis("exp(-abs(@0)/@1)*cosh(@0*@2/2)", RooArgList(tau,dgamma));
74 _basisSinh = declareBasis("exp(-abs(@0)/@1)*sinh(@0*@2/2)", RooArgList(tau,dgamma));
75 _basisCos = declareBasis("exp(-abs(@0)/@1)*cos(@0*@2)",RooArgList(tau, dm));
76 _basisSin = declareBasis("exp(-abs(@0)/@1)*sin(@0*@2)",RooArgList(tau, dm));
77 break;
78 }
79}
80
81////////////////////////////////////////////////////////////////////////////////
82///Copy constructor
83
84RooBDecay::RooBDecay(const RooBDecay& other, const char* name) :
85 RooAbsAnaConvPdf(other, name),
86 _t("t", this, other._t),
87 _tau("tau", this, other._tau),
88 _dgamma("dgamma", this, other._dgamma),
89 _f0("f0", this, other._f0),
90 _f1("f1", this, other._f1),
91 _f2("f2", this, other._f2),
92 _f3("f3", this, other._f3),
93 _dm("dm", this, other._dm),
94 _basisCosh(other._basisCosh),
95 _basisSinh(other._basisSinh),
96 _basisCos(other._basisCos),
97 _basisSin(other._basisSin),
98 _type(other._type)
99{
100}
101
102////////////////////////////////////////////////////////////////////////////////
103///Destructor
104
106{
107}
108
109////////////////////////////////////////////////////////////////////////////////
110
111double RooBDecay::coefficient(Int_t basisIndex) const
112{
113 if(basisIndex == _basisCosh)
114 {
115 return _f0;
116 }
117 if(basisIndex == _basisSinh)
118 {
119 return _f1;
120 }
121 if(basisIndex == _basisCos)
122 {
123 return _f2;
124 }
125 if(basisIndex == _basisSin)
126 {
127 return _f3;
128 }
129
130 return 0 ;
131}
132
133////////////////////////////////////////////////////////////////////////////////
134
136{
137 if(basisIndex == _basisCosh)
138 {
139 return _f0.arg().getVariables();
140 }
141 if(basisIndex == _basisSinh)
142 {
143 return _f1.arg().getVariables();
144 }
145 if(basisIndex == _basisCos)
146 {
147 return _f2.arg().getVariables();
148 }
149 if(basisIndex == _basisSin)
150 {
151 return _f3.arg().getVariables();
152 }
153
154 return 0 ;
155}
156
157////////////////////////////////////////////////////////////////////////////////
158
159Int_t RooBDecay::getCoefAnalyticalIntegral(Int_t coef, RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
160{
161 if(coef == _basisCosh)
162 {
163 return _f0.arg().getAnalyticalIntegral(allVars,analVars,rangeName) ;
164 }
165 if(coef == _basisSinh)
166 {
167 return _f1.arg().getAnalyticalIntegral(allVars,analVars,rangeName) ;
168 }
169 if(coef == _basisCos)
170 {
171 return _f2.arg().getAnalyticalIntegral(allVars,analVars,rangeName) ;
172 }
173 if(coef == _basisSin)
174 {
175 return _f3.arg().getAnalyticalIntegral(allVars,analVars,rangeName) ;
176 }
177
178 return 0 ;
179}
180
181////////////////////////////////////////////////////////////////////////////////
182
183double RooBDecay::coefAnalyticalIntegral(Int_t coef, Int_t code, const char* rangeName) const
184{
185 if(coef == _basisCosh)
186 {
187 return _f0.arg().analyticalIntegral(code,rangeName) ;
188 }
189 if(coef == _basisSinh)
190 {
191 return _f1.arg().analyticalIntegral(code,rangeName) ;
192 }
193 if(coef == _basisCos)
194 {
195 return _f2.arg().analyticalIntegral(code,rangeName) ;
196 }
197 if(coef == _basisSin)
198 {
199 return _f3.arg().analyticalIntegral(code,rangeName) ;
200 }
201
202 return 0 ;
203}
204
205////////////////////////////////////////////////////////////////////////////////
206
207Int_t RooBDecay::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, bool /*staticInitOK*/) const
208{
209 if (matchArgs(directVars, generateVars, _t)) return 1;
210 return 0;
211}
212
213////////////////////////////////////////////////////////////////////////////////
214
216{
217 R__ASSERT(code==1);
218 double gammamin = 1/_tau-TMath::Abs(_dgamma)/2;
219 while(1) {
220 double t = -log(RooRandom::uniform())/gammamin;
221 if (_type == Flipped || (_type == DoubleSided && RooRandom::uniform() <0.5) ) t *= -1;
222 if ( t<_t.min() || t>_t.max() ) continue;
223
224 double dgt = _dgamma*t/2;
225 double dmt = _dm*t;
226 double ft = fabs(t);
227 double f = exp(-ft/_tau)*(_f0*cosh(dgt)+_f1*sinh(dgt)+_f2*cos(dmt)+_f3*sin(dmt));
228 if(f < 0) {
229 cout << "RooBDecay::generateEvent(" << GetName() << ") ERROR: PDF value less than zero" << endl;
230 ::abort();
231 }
232 double w = 1.001*exp(-ft*gammamin)*(TMath::Abs(_f0)+TMath::Abs(_f1)+sqrt(_f2*_f2+_f3*_f3));
233 if(w < f) {
234 cout << "RooBDecay::generateEvent(" << GetName() << ") ERROR: Envelope function less than p.d.f. " << endl;
235 cout << _f0 << endl;
236 cout << _f1 << endl;
237 cout << _f2 << endl;
238 cout << _f3 << endl;
239 ::abort();
240 }
241 if(w*RooRandom::uniform() > f) continue;
242 _t = t;
243 break;
244 }
245}
#define f(i)
Definition: RSha256.hxx:104
#define ClassImp(name)
Definition: Rtypes.h:375
#define R__ASSERT(e)
Definition: TError.h:118
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.
RooArgSet * getVariables(bool stripDisconnected=true) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
Definition: RooAbsArg.cxx:2057
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
virtual double analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
Definition: RooAbsReal.cxx:434
bool 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:408
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:57
Most general description of B decay time distribution with effects of CP violation,...
Definition: RooBDecay.h:25
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...
Definition: RooBDecay.cxx:207
RooArgSet * coefVars(Int_t coefIdx) const override
Return set of parameters with are used exclusively by the coefficient functions.
Definition: RooBDecay.cxx:135
RooRealProxy _t
Definition: RooBDecay.h:58
RooRealProxy _f2
Definition: RooBDecay.h:63
Int_t _basisCosh
Definition: RooBDecay.h:66
DecayType _type
Definition: RooBDecay.h:71
Int_t _basisSin
Definition: RooBDecay.h:69
RooRealProxy _tau
Definition: RooBDecay.h:59
void generateEvent(Int_t code) override
Interface for generation of an event using the algorithm corresponding to the specified code.
Definition: RooBDecay.cxx:215
RooBDecay()
Definition: RooBDecay.h:32
RooRealProxy _f0
Definition: RooBDecay.h:61
double coefAnalyticalIntegral(Int_t coef, Int_t code, const char *rangeName=0) const override
Default implementation of function implementing advertised integrals.
Definition: RooBDecay.cxx:183
@ DoubleSided
Definition: RooBDecay.h:29
@ SingleSided
Definition: RooBDecay.h:29
double coefficient(Int_t basisIndex) const override
Definition: RooBDecay.cxx:111
RooRealProxy _dgamma
Definition: RooBDecay.h:60
RooRealProxy _f3
Definition: RooBDecay.h:64
Int_t _basisSinh
Definition: RooBDecay.h:67
RooRealProxy _dm
Definition: RooBDecay.h:65
~RooBDecay() override
Destructor.
Definition: RooBDecay.cxx:105
Int_t getCoefAnalyticalIntegral(Int_t coef, RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const override
Default implementation of function advertising integration capabilities.
Definition: RooBDecay.cxx:159
Int_t _basisCos
Definition: RooBDecay.h:68
RooRealProxy _f1
Definition: RooBDecay.h:62
static double uniform(TRandom *generator=randomGenerator())
Return a number uniformly distributed from (0,1)
Definition: RooRandom.cxx:81
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 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.
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
RVec< PromoteType< T > > cosh(const RVec< T > &v)
Definition: RVec.hxx:1767
RVec< PromoteType< T > > cos(const RVec< T > &v)
Definition: RVec.hxx:1759
RVec< PromoteType< T > > log(const RVec< T > &v)
Definition: RVec.hxx:1748
RVec< PromoteType< T > > sinh(const RVec< T > &v)
Definition: RVec.hxx:1766
RVec< PromoteType< T > > exp(const RVec< T > &v)
Definition: RVec.hxx:1744
RVec< PromoteType< T > > sin(const RVec< T > &v)
Definition: RVec.hxx:1758
TF1 * f1
Definition: legend1.C:11
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition: TMathBase.h:123