Logo ROOT  
Reference Guide
RooDecay.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 RooDecay
18  \ingroup Roofit
19 
20 Single or double sided decay function that can be analytically convolved
21 with any RooResolutionModel implementation. It declares the basis functions
22 for the analytical convolution with a RooResolutionModel. See RooAbsAnaConvPdf.
23 \f[
24  \mathrm{basis} = \begin{cases}
25  \exp\left(-\frac{t}{\tau}\right) & \mathrm{SingleSided} \\
26  \exp\left( \frac{t}{\tau}\right) & \mathrm{Flipped} \\
27  \exp\left(-\frac{|t|}{\tau}\right) & \mathrm{DoubleSided}
28  \end{cases}
29 \f]
30 **/
31 
32 #include "RooDecay.h"
33 
34 #include "RooFit.h"
35 #include "RooRealVar.h"
36 #include "RooRandom.h"
37 
38 #include "TError.h"
39 
40 using namespace std;
41 
43 
44 ////////////////////////////////////////////////////////////////////////////////
45 /// Create a new RooDecay.
46 /// \param[in] name Name of this object.
47 /// \param[in] title Title (for *e.g.* plotting)
48 /// \param[in] t Convolution variable (*e.g.* time).
49 /// \param[in] tau Decay constant.
50 /// \param[in] model Resolution model for the convolution.
51 /// \param[in] type One of the decays types `SingleSided, Flipped, DoubleSided`
52 RooDecay::RooDecay(const char *name, const char *title,
53  RooRealVar& t, RooAbsReal& tau,
55  RooAbsAnaConvPdf(name,title,model,t),
56  _t("t","time",this,t),
57  _tau("tau","decay time",this,tau),
58  _type(type)
59 {
60  switch(type) {
61  case SingleSided:
62  _basisExp = declareBasis("exp(-@0/@1)",tau) ;
63  break ;
64  case Flipped:
65  _basisExp = declareBasis("exp(@0/@1)",tau) ;
66  break ;
67  case DoubleSided:
68  _basisExp = declareBasis("exp(-abs(@0)/@1)",tau) ;
69  break ;
70  }
71 }
72 
73 ////////////////////////////////////////////////////////////////////////////////
74 /// Copy constructor
75 
76 RooDecay::RooDecay(const RooDecay& other, const char* name) :
77  RooAbsAnaConvPdf(other,name),
78  _t("t",this,other._t),
79  _tau("tau",this,other._tau),
80  _type(other._type),
81  _basisExp(other._basisExp)
82 {
83 }
84 
85 ////////////////////////////////////////////////////////////////////////////////
86 /// Destructor
87 
89 {
90 }
91 
92 ////////////////////////////////////////////////////////////////////////////////
93 
94 Double_t RooDecay::coefficient(Int_t /*basisIndex*/) const
95 {
96  return 1 ;
97 }
98 
99 ////////////////////////////////////////////////////////////////////////////////
100 
101 Int_t RooDecay::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t /*staticInitOK*/) const
102 {
103  if (matchArgs(directVars,generateVars,_t)) return 1 ;
104  return 0 ;
105 }
106 
107 ////////////////////////////////////////////////////////////////////////////////
108 
110 {
111  R__ASSERT(code==1) ;
112 
113  // Generate delta-t dependent
114  while(1) {
115  Double_t rand = RooRandom::uniform() ;
116  Double_t tval(0) ;
117 
118  switch(_type) {
119  case SingleSided:
120  tval = -_tau*log(rand);
121  break ;
122  case Flipped:
123  tval= +_tau*log(rand);
124  break ;
125  case DoubleSided:
126  tval = (rand<=0.5) ? -_tau*log(2*rand) : +_tau*log(2*(rand-0.5)) ;
127  break ;
128  }
129 
130  if (tval<_t.max() && tval>_t.min()) {
131  _t = tval ;
132  break ;
133  }
134  }
135 }
RooFit.h
RooAbsAnaConvPdf::declareBasis
Int_t declareBasis(const char *expression, const RooArgList &params)
Declare a basis function for use in this physics model.
Definition: RooAbsAnaConvPdf.cxx:158
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
RooDecay::~RooDecay
virtual ~RooDecay()
Destructor.
Definition: RooDecay.cxx:88
RooDecay::SingleSided
@ SingleSided
Definition: RooDecay.h:25
log
double log(double)
RooDecay::getGenerator
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: RooDecay.cxx:101
RooDecay::coefficient
virtual Double_t coefficient(Int_t basisIndex) const
Definition: RooDecay.cxx:94
RooAbsReal
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
RooDecay::RooDecay
RooDecay()
Definition: RooDecay.h:28
RooDecay.h
RooDecay::_t
RooRealProxy _t
Definition: RooDecay.h:42
bool
RooResolutionModel
RooResolutionModel is the base class for PDFs that represent a resolution model that can be convolute...
Definition: RooResolutionModel.h:26
RooTemplateProxy::min
double min(const char *rname=0) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
Definition: RooTemplateProxy.h:275
RooDecay
Single or double sided decay function that can be analytically convolved with any RooResolutionModel ...
Definition: RooDecay.h:22
RooRandom.h
RooDecay::DecayType
DecayType
Definition: RooDecay.h:25
RooRealVar.h
RooDecay::_basisExp
Int_t _basisExp
Definition: RooDecay.h:45
Double_t
double Double_t
Definition: RtypesCore.h:59
R__ASSERT
#define R__ASSERT(e)
Definition: TError.h:118
RooDecay::_type
DecayType _type
Definition: RooDecay.h:44
name
char name[80]
Definition: TGX11.cxx:110
RooDecay::generateEvent
void generateEvent(Int_t code)
Interface for generation of an event using the algorithm corresponding to the specified code.
Definition: RooDecay.cxx:109
RooTemplateProxy::max
double max(const char *rname=0) const
Query upper limit of range. This requires the payload to be RooAbsRealLValue or derived.
Definition: RooTemplateProxy.h:277
make_cnn_model.model
model
Definition: make_cnn_model.py:6
RooAbsReal::matchArgs
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
Definition: RooAbsReal.cxx:3389
type
int type
Definition: TGX11.cxx:121
RooDecay::DoubleSided
@ DoubleSided
Definition: RooDecay.h:25
RooRealVar
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:39
RooDecay::Flipped
@ Flipped
Definition: RooDecay.h:25
RooDecay::_tau
RooRealProxy _tau
Definition: RooDecay.h:43
RooRandom::uniform
static Double_t uniform(TRandom *generator=randomGenerator())
Return a number uniformly distributed from (0,1)
Definition: RooRandom.cxx:83
RooAbsAnaConvPdf
RooAbsAnaConvPdf is the base class for PDFs that represent a physics model that can be analytically c...
Definition: RooAbsAnaConvPdf.h:34
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:29
int
TError.h