Logo ROOT   6.10/09
Reference Guide
RooBCPGenDecay.cxx
Go to the documentation of this file.
1 /*****************************************************************************
2  * Project: RooFit *
3  * Package: RooFitModels *
4  * @(#)root/roofit:$Id$
5  * Authors: *
6  * JGS, Jim Smith , University of Colorado, jgsmith@pizero.colorado.edu *
7  * History:
8  * 15-Aug-2002 JGS Created initial version
9  * 11-Sep-2002 JGS Mods to introduce mu (Mirna van Hoek, JGS, Nick Danielson)
10  * *
11  * Copyright (c) 2000-2005, Regents of the University of California, *
12  * University of Colorado *
13  * and Stanford University. All rights reserved. *
14  * *
15  * Redistribution and use in source and binary forms, *
16  * with or without modification, are permitted according to the terms *
17  * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
18  *****************************************************************************/
19 
20 /** \class RooBCPGenDecay
21  \ingroup Roofit
22 
23 Implement standard CP physics model with S and C (no mention of lambda)
24 Suitably stolen and modified from RooBCPEffDecay
25 **/
26 
27 #include "RooRealVar.h"
28 #include "RooRandom.h"
29 #include "RooBCPGenDecay.h"
30 #include "RooRealIntegral.h"
31 
32 using namespace std;
33 
35 
36 ////////////////////////////////////////////////////////////////////////////////
37 /// Constructor
38 
39 RooBCPGenDecay::RooBCPGenDecay(const char *name, const char *title,
40  RooRealVar& t, RooAbsCategory& tag,
42  RooAbsReal& avgMistag,
44  RooAbsReal& delMistag,
45  RooAbsReal& mu,
47  RooAbsAnaConvPdf(name,title,model,t),
48  _avgC("C","Coefficient of cos term",this,a),
49  _avgS("S","Coefficient of sin term",this,b),
50  _avgMistag("avgMistag","Average mistag rate",this,avgMistag),
51  _delMistag("delMistag","Delta mistag rate",this,delMistag),
52  _mu("mu","Tagg efficiency difference",this,mu),
53  _t("t","time",this,t),
54  _tau("tau","decay time",this,tau),
55  _dm("dm","mixing frequency",this,dm),
56  _tag("tag","CP state",this,tag),
57  _genB0Frac(0),
58  _type(type)
59 {
60  switch(type) {
61  case SingleSided:
62  _basisExp = declareBasis("exp(-@0/@1)",RooArgList(tau,dm)) ;
63  _basisSin = declareBasis("exp(-@0/@1)*sin(@0*@2)",RooArgList(tau,dm)) ;
64  _basisCos = declareBasis("exp(-@0/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
65  break ;
66  case Flipped:
67  _basisExp = declareBasis("exp(@0)/@1)",RooArgList(tau,dm)) ;
68  _basisSin = declareBasis("exp(@0/@1)*sin(@0*@2)",RooArgList(tau,dm)) ;
69  _basisCos = declareBasis("exp(@0/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
70  break ;
71  case DoubleSided:
72  _basisExp = declareBasis("exp(-abs(@0)/@1)",RooArgList(tau,dm)) ;
73  _basisSin = declareBasis("exp(-abs(@0)/@1)*sin(@0*@2)",RooArgList(tau,dm)) ;
74  _basisCos = declareBasis("exp(-abs(@0)/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
75  break ;
76  }
77 }
78 
79 ////////////////////////////////////////////////////////////////////////////////
80 /// Copy constructor
81 
83  RooAbsAnaConvPdf(other,name),
84  _avgC("C",this,other._avgC),
85  _avgS("S",this,other._avgS),
86  _avgMistag("avgMistag",this,other._avgMistag),
87  _delMistag("delMistag",this,other._delMistag),
88  _mu("mu",this,other._mu),
89  _t("t",this,other._t),
90  _tau("tau",this,other._tau),
91  _dm("dm",this,other._dm),
92  _tag("tag",this,other._tag),
93  _genB0Frac(other._genB0Frac),
94  _type(other._type),
95  _basisExp(other._basisExp),
96  _basisSin(other._basisSin),
97  _basisCos(other._basisCos)
98 {
99 }
100 
101 ////////////////////////////////////////////////////////////////////////////////
102 /// Destructor
103 
105 {
106 }
107 
108 ////////////////////////////////////////////////////////////////////////////////
109 /// B0 : _tag = +1
110 /// B0bar : _tag = -1
111 
113 {
114  if (basisIndex==_basisExp) {
115  //exp term: (1 -/+ dw + mu*_tag*w)
116  return (1 - _tag*_delMistag + _mu*_tag*(1. - 2.*_avgMistag)) ;
117  // = 1 + _tag*deltaDil/2 + _mu*avgDil
118  }
119 
120  if (basisIndex==_basisSin) {
121  //sin term: (+/- (1-2w) + _mu*(1 -/+ delw))*S
122  return (_tag*(1-2*_avgMistag) + _mu*(1. - _tag*_delMistag))*_avgS ;
123  // = (_tag*avgDil + _mu*(1 + tag*deltaDil/2)) * S
124  }
125 
126  if (basisIndex==_basisCos) {
127  //cos term: -(+/- (1-2w) + _mu*(1 -/+ delw))*C
128  return -1.*(_tag*(1-2*_avgMistag) + _mu*(1. - _tag*_delMistag))*_avgC ;
129  // = -(_tag*avgDil + _mu*(1 + _tag*deltaDil/2) )* C
130  }
131 
132  return 0 ;
133 }
134 
135 ////////////////////////////////////////////////////////////////////////////////
136 
137 Int_t RooBCPGenDecay::getCoefAnalyticalIntegral(Int_t /*code*/, RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
138 {
139  if (rangeName) return 0 ;
140  if (matchArgs(allVars,analVars,_tag)) return 1 ;
141  return 0 ;
142 }
143 
144 ////////////////////////////////////////////////////////////////////////////////
145 
146 Double_t RooBCPGenDecay::coefAnalyticalIntegral(Int_t basisIndex, Int_t code, const char* /*rangeName*/) const
147 {
148  switch(code) {
149  // No integration
150  case 0: return coefficient(basisIndex) ;
151 
152  // Integration over 'tag'
153  case 1:
154  if (basisIndex==_basisExp) {
155  return 2 ;
156  }
157 
158  if (basisIndex==_basisSin) {
159  return 2*_mu*_avgS ;
160  }
161  if (basisIndex==_basisCos) {
162  return -2*_mu*_avgC ;
163  }
164  break ;
165 
166  default:
167  assert(0) ;
168  }
169 
170  return 0 ;
171 }
172 
173 ////////////////////////////////////////////////////////////////////////////////
174 
175 Int_t RooBCPGenDecay::getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t staticInitOK) const
176 {
177  if (staticInitOK) {
178  if (matchArgs(directVars,generateVars,_t,_tag)) return 2 ;
179  }
180  if (matchArgs(directVars,generateVars,_t)) return 1 ;
181  return 0 ;
182 }
183 
184 ////////////////////////////////////////////////////////////////////////////////
185 
187 {
188  if (code==2) {
189  // Calculate the fraction of mixed events to generate
190  Double_t sumInt = RooRealIntegral("sumInt","sum integral",*this,RooArgSet(_t.arg(),_tag.arg())).getVal() ;
191  _tag = 1 ;
192  Double_t b0Int = RooRealIntegral("mixInt","mix integral",*this,RooArgSet(_t.arg())).getVal() ;
193  _genB0Frac = b0Int/sumInt ;
194  }
195 }
196 
197 ////////////////////////////////////////////////////////////////////////////////
198 /// Generate mix-state dependent
199 
201 {
202  if (code==2) {
203  Double_t rand = RooRandom::uniform() ;
204  _tag = (rand<=_genB0Frac) ? 1 : -1 ;
205  }
206 
207  // Generate delta-t dependent
208  while(1) {
209  Double_t rand = RooRandom::uniform() ;
210  Double_t tval(0) ;
211 
212  switch(_type) {
213  case SingleSided:
214  tval = -_tau*log(rand);
215  break ;
216  case Flipped:
217  tval= +_tau*log(rand);
218  break ;
219  case DoubleSided:
220  tval = (rand<=0.5) ? -_tau*log(2*rand) : +_tau*log(2*(rand-0.5)) ;
221  break ;
222  }
223 
224  // Accept event if T is in generated range
225  Double_t maxDil = 1.0 ;
226 // 2 in next line is conservative and inefficient - allows for delMistag=1!
227  Double_t maxAcceptProb = 2 + fabs(maxDil*_avgS) + fabs(maxDil*_avgC);
228  Double_t acceptProb = (1-_tag*_delMistag + _mu*_tag*(1. - 2.*_avgMistag))
229  + (_tag*(1-2*_avgMistag) + _mu*(1. - _tag*_delMistag))*_avgS*sin(_dm*tval)
230  - (_tag*(1-2*_avgMistag) + _mu*(1. - _tag*_delMistag))*_avgC*cos(_dm*tval);
231 
232  Bool_t accept = maxAcceptProb*RooRandom::uniform() < acceptProb ? kTRUE : kFALSE ;
233 
234  if (tval<_t.max() && tval>_t.min() && accept) {
235  _t = tval ;
236  break ;
237  }
238  }
239 
240 }
RooRealProxy _mu
Implement standard CP physics model with S and C (no mention of lambda) Suitably stolen and modified ...
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
RooRealProxy _avgC
RooRealProxy _dm
virtual ~RooBCPGenDecay()
Destructor.
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
TArc * a
Definition: textangle.C:12
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...
STL namespace.
double cos(double)
virtual Double_t coefficient(Int_t basisIndex) const
B0 : _tag = +1 B0bar : _tag = -1.
Double_t _genB0Frac
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
RooRealProxy _tau
const RooAbsCategory & arg() const
RooCategoryProxy _tag
friend class RooArgSet
Definition: RooAbsArg.h:469
double sin(double)
DecayType _type
RooRealProxy _delMistag
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
RooArgSet S(const RooAbsArg &v1)
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
static double C[]
void initGenerator(Int_t code)
Interface for one-time initialization to setup the generator for the specified code.
const Bool_t kFALSE
Definition: RtypesCore.h:92
#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
RooRealProxy _avgMistag
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
void generateEvent(Int_t code)
Generate mix-state dependent.
int type
Definition: TGX11.cxx:120
Double_t max(const char *rname=0) const
Definition: RooRealProxy.h:57
friend class RooRealIntegral
Definition: RooAbsPdf.h:294
RooAbsCategory is the common abstract base class for objects that represent a discrete value with a f...
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 Double_t b
Definition: TRolke.cxx:630
virtual Double_t coefAnalyticalIntegral(Int_t coef, Int_t code, const char *rangeName=0) const
Default implementation of function implementing advertised integrals.
RooRealProxy _t
const RooAbsReal & arg() const
Definition: RooRealProxy.h:43
const Bool_t kTRUE
Definition: RtypesCore.h:91
virtual Int_t getCoefAnalyticalIntegral(Int_t coef, RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Default implementation of function advertising integration capabilities.
RooRealProxy _avgS
double log(double)