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
20Single or double sided decay function that can be analytically convolved
21with any RooResolutionModel implementation. It declares the basis functions
22for 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 "RooRealVar.h"
35#include "RooRandom.h"
36
37#include "TError.h"
38
39using namespace std;
40
42
43////////////////////////////////////////////////////////////////////////////////
44/// Create a new RooDecay.
45/// \param[in] name Name of this object.
46/// \param[in] title Title (for *e.g.* plotting)
47/// \param[in] t Convolution variable (*e.g.* time).
48/// \param[in] tau Decay constant.
49/// \param[in] model Resolution model for the convolution.
50/// \param[in] type One of the decays types `SingleSided, Flipped, DoubleSided`
51RooDecay::RooDecay(const char *name, const char *title,
52 RooRealVar& t, RooAbsReal& tau,
53 const RooResolutionModel& model, DecayType type) :
54 RooAbsAnaConvPdf(name,title,model,t),
55 _t("t","time",this,t),
56 _tau("tau","decay time",this,tau),
57 _type(type)
58{
59 switch(type) {
60 case SingleSided:
61 _basisExp = declareBasis("exp(-@0/@1)",tau) ;
62 break ;
63 case Flipped:
64 _basisExp = declareBasis("exp(@0/@1)",tau) ;
65 break ;
66 case DoubleSided:
67 _basisExp = declareBasis("exp(-abs(@0)/@1)",tau) ;
68 break ;
69 }
70}
71
72////////////////////////////////////////////////////////////////////////////////
73/// Copy constructor
74
75RooDecay::RooDecay(const RooDecay& other, const char* name) :
77 _t("t",this,other._t),
78 _tau("tau",this,other._tau),
79 _type(other._type),
80 _basisExp(other._basisExp)
81{
82}
83
84////////////////////////////////////////////////////////////////////////////////
85/// Destructor
86
88{
89}
90
91////////////////////////////////////////////////////////////////////////////////
92
93double RooDecay::coefficient(Int_t /*basisIndex*/) const
94{
95 return 1 ;
96}
97
98////////////////////////////////////////////////////////////////////////////////
99
100Int_t RooDecay::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, bool /*staticInitOK*/) const
101{
102 if (matchArgs(directVars,generateVars,_t)) return 1 ;
103 return 0 ;
104}
105
106////////////////////////////////////////////////////////////////////////////////
107
109{
110 R__ASSERT(code==1) ;
111
112 // Generate delta-t dependent
113 while(1) {
114 double rand = RooRandom::uniform() ;
115 double tval(0) ;
116
117 switch(_type) {
118 case SingleSided:
119 tval = -_tau*log(rand);
120 break ;
121 case Flipped:
122 tval= +_tau*log(rand);
123 break ;
124 case DoubleSided:
125 tval = (rand<=0.5) ? -_tau*log(2*rand) : +_tau*log(2*(rand-0.5)) ;
126 break ;
127 }
128
129 if (tval<_t.max() && tval>_t.min()) {
130 _t = tval ;
131 break ;
132 }
133 }
134}
#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.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
bool matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:57
Single or double sided decay function that can be analytically convolved with any RooResolutionModel ...
Definition: RooDecay.h:22
double coefficient(Int_t basisIndex) const override
Definition: RooDecay.cxx:93
RooDecay()
Definition: RooDecay.h:28
@ DoubleSided
Definition: RooDecay.h:25
@ SingleSided
Definition: RooDecay.h:25
@ Flipped
Definition: RooDecay.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: RooDecay.cxx:100
Int_t _basisExp
Definition: RooDecay.h:45
RooRealProxy _t
Definition: RooDecay.h:42
RooRealProxy _tau
Definition: RooDecay.h:43
void generateEvent(Int_t code) override
Interface for generation of an event using the algorithm corresponding to the specified code.
Definition: RooDecay.cxx:108
~RooDecay() override
Destructor.
Definition: RooDecay.cxx:87
DecayType _type
Definition: RooDecay.h:44
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.
RVec< PromoteType< T > > log(const RVec< T > &v)
Definition: RVec.hxx:1748