Logo ROOT   6.10/09
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
22 **/
23 
24 #include "RooFit.h"
25 
26 #include "Riostream.h"
27 #include "Riostream.h"
28 #include "RooDecay.h"
29 #include "RooRealVar.h"
30 #include "RooRandom.h"
31 
32 #include "TError.h"
33 
34 using namespace std;
35 
37 
38 ////////////////////////////////////////////////////////////////////////////////
39 /// Constructor
40 
41 RooDecay::RooDecay(const char *name, const char *title,
44  RooAbsAnaConvPdf(name,title,model,t),
45  _t("t","time",this,t),
46  _tau("tau","decay time",this,tau),
47  _type(type)
48 {
49  switch(type) {
50  case SingleSided:
51  _basisExp = declareBasis("exp(-@0/@1)",tau) ;
52  break ;
53  case Flipped:
54  _basisExp = declareBasis("exp(@0/@1)",tau) ;
55  break ;
56  case DoubleSided:
57  _basisExp = declareBasis("exp(-abs(@0)/@1)",tau) ;
58  break ;
59  }
60 }
61 
62 ////////////////////////////////////////////////////////////////////////////////
63 /// Copy constructor
64 
65 RooDecay::RooDecay(const RooDecay& other, const char* name) :
66  RooAbsAnaConvPdf(other,name),
67  _t("t",this,other._t),
68  _tau("tau",this,other._tau),
69  _type(other._type),
70  _basisExp(other._basisExp)
71 {
72 }
73 
74 ////////////////////////////////////////////////////////////////////////////////
75 /// Destructor
76 
78 {
79 }
80 
81 ////////////////////////////////////////////////////////////////////////////////
82 
83 Double_t RooDecay::coefficient(Int_t /*basisIndex*/) const
84 {
85  return 1 ;
86 }
87 
88 ////////////////////////////////////////////////////////////////////////////////
89 
90 Int_t RooDecay::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t /*staticInitOK*/) const
91 {
92  if (matchArgs(directVars,generateVars,_t)) return 1 ;
93  return 0 ;
94 }
95 
96 ////////////////////////////////////////////////////////////////////////////////
97 
99 {
100  R__ASSERT(code==1) ;
101 
102  // Generate delta-t dependent
103  while(1) {
104  Double_t rand = RooRandom::uniform() ;
105  Double_t tval(0) ;
106 
107  switch(_type) {
108  case SingleSided:
109  tval = -_tau*log(rand);
110  break ;
111  case Flipped:
112  tval= +_tau*log(rand);
113  break ;
114  case DoubleSided:
115  tval = (rand<=0.5) ? -_tau*log(2*rand) : +_tau*log(2*(rand-0.5)) ;
116  break ;
117  }
118 
119  if (tval<_t.max() && tval>_t.min()) {
120  _t = tval ;
121  break ;
122  }
123  }
124 }
RooRealProxy _tau
Definition: RooDecay.h:43
virtual ~RooDecay()
Destructor.
Definition: RooDecay.cxx:77
RooRealProxy _t
Definition: RooDecay.h:42
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
virtual Double_t coefficient(Int_t basisIndex) const
Definition: RooDecay.cxx:83
RooDecay()
Definition: RooDecay.h:28
#define R__ASSERT(e)
Definition: TError.h:96
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
STL namespace.
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:90
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
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
#define ClassImp(name)
Definition: Rtypes.h:336
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
int type
Definition: TGX11.cxx:120
void generateEvent(Int_t code)
Interface for generation of anan event using the algorithm corresponding to the specified code...
Definition: RooDecay.cxx:98
DecayType _type
Definition: RooDecay.h:44
Double_t max(const char *rname=0) const
Definition: RooRealProxy.h:57
double log(double)
Single or double sided decay function that can be analytically convolved with any RooResolutionModel ...
Definition: RooDecay.h:22