Logo ROOT  
Reference Guide
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
Loading...
Searching...
No Matches
RooBCPEffDecay.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 RooBCPEffDecay
18 \ingroup Roofit
19
20PDF describing decay time distribution of B meson including effects of standard model CP violation.
21This function can be analytically convolved with any RooResolutionModel implementation.
22*/
23
24
25#include "Riostream.h"
26#include "RooRealVar.h"
27#include "RooRandom.h"
28#include "RooBCPEffDecay.h"
29#include "RooRealIntegral.h"
30
31
32/// \brief Constructor for RooBCPEffDecay.
33///
34/// Creates an instance of RooBCPEffDecay with the specified parameters.
35///
36/// \param[in] name The name of the PDF.
37/// \param[in] title The title of the PDF.
38/// \param[in] t The time variable.
39/// \param[in] tag The CP state category.
40/// \param[in] tau The decay time parameter.
41/// \param[in] dm The mixing frequency parameter.
42/// \param[in] avgMistag The average mistag rate parameter.
43/// \param[in] CPeigenval The CP eigen value parameter.
44/// \param[in] absLambda The absolute value of the complex lambda parameter.
45/// \param[in] argLambda The argument of the complex lambda parameter.
46/// \param[in] effRatio The B0/B0bar efficiency ratio.
47/// \param[in] delMistag Delta mistag rate parameter.
48/// \param[in] model The resolution model.
49/// \param[in] type The decay type.
50
51RooBCPEffDecay::RooBCPEffDecay(const char *name, const char *title,
53 RooAbsReal& tau, RooAbsReal& dm,
57 const RooResolutionModel& model, DecayType type) :
58 RooAbsAnaConvPdf(name,title,model,t),
59 _absLambda("absLambda","Absolute value of lambda",this,absLambda),
60 _argLambda("argLambda","Arg(Lambda)",this,argLambda),
61 _effRatio("effRatio","B0/B0bar efficiency ratio",this,effRatio),
62 _CPeigenval("CPeigenval","CP eigen value",this,CPeigenval),
63 _avgMistag("avgMistag","Average mistag rate",this,avgMistag),
64 _delMistag("delMistag","Delta mistag rate",this,delMistag),
65 _t("t","time",this,t),
66 _tau("tau","decay time",this,tau),
67 _dm("dm","mixing frequency",this,dm),
68 _tag("tag","CP state",this,tag),
69 _genB0Frac(0),
70 _type(type)
71{
72 switch(type) {
73 case SingleSided:
74 _basisExp = declareBasis("exp(-@0/@1)",RooArgList(tau,dm)) ;
75 _basisSin = declareBasis("exp(-@0/@1)*sin(@0*@2)",RooArgList(tau,dm)) ;
76 _basisCos = declareBasis("exp(-@0/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
77 break ;
78 case Flipped:
79 _basisExp = declareBasis("exp(@0)/@1)",RooArgList(tau,dm)) ;
80 _basisSin = declareBasis("exp(@0/@1)*sin(@0*@2)",RooArgList(tau,dm)) ;
81 _basisCos = declareBasis("exp(@0/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
82 break ;
83 case DoubleSided:
84 _basisExp = declareBasis("exp(-abs(@0)/@1)",RooArgList(tau,dm)) ;
85 _basisSin = declareBasis("exp(-abs(@0)/@1)*sin(@0*@2)",RooArgList(tau,dm)) ;
86 _basisCos = declareBasis("exp(-abs(@0)/@1)*cos(@0*@2)",RooArgList(tau,dm)) ;
87 break ;
88 }
89}
90
91////////////////////////////////////////////////////////////////////////////////
92/// Copy constructor.
93
96 _absLambda("absLambda",this,other._absLambda),
97 _argLambda("argLambda",this,other._argLambda),
98 _effRatio("effRatio",this,other._effRatio),
99 _CPeigenval("CPeigenval",this,other._CPeigenval),
100 _avgMistag("avgMistag",this,other._avgMistag),
101 _delMistag("delMistag",this,other._delMistag),
102 _t("t",this,other._t),
103 _tau("tau",this,other._tau),
104 _dm("dm",this,other._dm),
105 _tag("tag",this,other._tag),
106 _genB0Frac(other._genB0Frac),
107 _type(other._type),
108 _basisExp(other._basisExp),
109 _basisSin(other._basisSin),
110 _basisCos(other._basisCos)
111{
112}
113
114////////////////////////////////////////////////////////////////////////////////
115/// B0 : _tag = +1
116///
117/// B0bar : _tag = -1
118/// \param[in] basisIndex
119
121{
122 if (basisIndex==_basisExp) {
123 //exp term: (1 -/+ dw)(1+a^2)/2
124 return (1 - _tag*_delMistag)*(1+_absLambda*_absLambda)/2 ;
125 // = 1 + _tag*deltaDil/2
126 }
127
128 if (basisIndex==_basisSin) {
129 //sin term: +/- (1-2w)*ImLambda(= -etaCP * |l| * sin2b)
131 // = _tag*avgDil * ...
132 }
133
134 if (basisIndex==_basisCos) {
135 //cos term: +/- (1-2w)*(1-a^2)/2
136 return -1*_tag*(1-2*_avgMistag)*(1-_absLambda*_absLambda)/2 ;
137 // = -_tag*avgDil * ...
138 }
139
140 return 0 ;
141}
142
143////////////////////////////////////////////////////////////////////////////////
144
145Int_t RooBCPEffDecay::getCoefAnalyticalIntegral(Int_t /*code*/, RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
146{
147 if (rangeName) return 0 ;
148
149 if (matchArgs(allVars,analVars,_tag)) return 1 ;
150 return 0 ;
151}
152
153////////////////////////////////////////////////////////////////////////////////
154
155double RooBCPEffDecay::coefAnalyticalIntegral(Int_t basisIndex, Int_t code, const char* /*rangeName*/) const
156{
157 switch(code) {
158 // No integration
159 case 0: return coefficient(basisIndex) ;
160
161 // Integration over 'tag'
162 case 1:
163 if (basisIndex==_basisExp) {
164 return (1+_absLambda*_absLambda) ;
165 }
166
168 return 0 ;
169 }
170 break ;
171
172 default:
173 assert(0) ;
174 }
175
176 return 0 ;
177}
178
179////////////////////////////////////////////////////////////////////////////////
180
182{
183 if (staticInitOK) {
184 if (matchArgs(directVars,generateVars,_t,_tag)) return 2 ;
185 }
186 if (matchArgs(directVars,generateVars,_t)) return 1 ;
187 return 0 ;
188}
189
190////////////////////////////////////////////////////////////////////////////////
191
193{
194 if (code==2) {
195 // Calculate the fraction of mixed events to generate
196 double sumInt = RooRealIntegral("sumInt","sum integral",*this,RooArgSet(_t.arg(),_tag.arg())).getVal() ;
197 _tag = 1 ;
198 double b0Int = RooRealIntegral("mixInt","mix integral",*this,RooArgSet(_t.arg())).getVal() ;
200 }
201}
202
203////////////////////////////////////////////////////////////////////////////////
204/// Generates mix-state dependent.
205/// \param[in] code
206
208{
209 if (code==2) {
210 double rand = RooRandom::uniform() ;
211 _tag = (rand<=_genB0Frac) ? 1 : -1 ;
212 }
213
214 // Generate delta-t dependent
215 while(true) {
216 double rand = RooRandom::uniform() ;
217 double tval(0) ;
218
219 switch(_type) {
220 case SingleSided:
221 tval = -_tau*log(rand);
222 break ;
223 case Flipped:
224 tval= +_tau*log(rand);
225 break ;
226 case DoubleSided:
227 tval = (rand<=0.5) ? -_tau*log(2*rand) : +_tau*log(2*(rand-0.5)) ;
228 break ;
229 }
230
231 // Accept event if T is in generated range
232 double maxDil = 1.0 ;
233 double al2 = _absLambda*_absLambda ;
234 double maxAcceptProb = (1+al2) + std::abs(maxDil*_CPeigenval*_absLambda*_argLambda) + std::abs(maxDil*(1-al2)/2);
235 double acceptProb = (1+al2)/2*(1-_tag*_delMistag)
237 - (_tag*(1-2*_avgMistag))*(1-al2)/2*cos(_dm*tval);
238
240
241 if (tval<_t.max() && tval>_t.min() && accept) {
242 _t = tval ;
243 break ;
244 }
245 }
246
247}
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
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
Base class for PDFs that represent a physics model that can be analytically convolved with a resoluti...
Int_t declareBasis(const char *expression, const RooArgList &params)
Declare a basis function for use in this physics model.
friend class RooRealIntegral
Definition RooAbsArg.h:572
A space to attach TBranches.
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
bool matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
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:24
PDF describing decay time distribution of B meson including effects of standard model CP violation.
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...
double coefAnalyticalIntegral(Int_t coef, Int_t code, const char *rangeName=nullptr) const override
Default implementation of function implementing advertised integrals.
double coefficient(Int_t basisIndex) const override
B0 : _tag = +1.
Int_t getCoefAnalyticalIntegral(Int_t coef, RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=nullptr) const override
Default implementation of function advertising integration capabilities.
RooRealProxy _t
RooRealProxy _avgMistag
RooRealProxy _delMistag
void generateEvent(Int_t code) override
Generates mix-state dependent.
RooRealProxy _dm
RooRealProxy _absLambda
RooRealProxy _argLambda
RooCategoryProxy _tag
RooRealProxy _CPeigenval
RooRealProxy _tau
void initGenerator(Int_t code) override
Interface for one-time initialization to setup the generator for the specified code.
static double uniform(TRandom *generator=randomGenerator())
Return a number uniformly distributed from (0,1)
Definition RooRandom.cxx:77
Variable that can be changed from the outside.
Definition RooRealVar.h:37
RooResolutionModel is the base class for PDFs that represent a resolution model that can be convolute...
double max(const char *rname=nullptr) 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.
double min(const char *rname=nullptr) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.