Converting TEvtGen to native cmake
[u/mrichter/AliRoot.git] / TEvtGen / EvtGen / EvtGenBase / EvtIncoherentMixing.cpp
CommitLineData
0ca57c2f 1// $Id: EvtIncoherentMixing.cpp,v 1.13 2009-11-27 09:09:41 mwhitehe Exp $
2// Include files
3
4
5// local
6#include "EvtGenBase/EvtIncoherentMixing.hh"
7#include <stdlib.h>
8#include "EvtGenBase/EvtPDL.hh"
9#include "EvtGenBase/EvtId.hh"
10#include "EvtGenBase/EvtRandom.hh"
11
12//-----------------------------------------------------------------------------
13// Implementation file for class : EvtIncoherentMixing
14//
15// 2003-10-09 : Patrick Robbe
16//-----------------------------------------------------------------------------
17
18
19bool EvtIncoherentMixing::_doB0Mixing = false ;
20bool EvtIncoherentMixing::_doBsMixing = false ;
21bool EvtIncoherentMixing::_enableFlip = false ;
22double EvtIncoherentMixing::_dGammad = 0. ;
23double EvtIncoherentMixing::_deltamd = 0.502e12 ;
24// dGamma_s corresponds to DeltaGamma / Gamma = 10 %
25double EvtIncoherentMixing::_dGammas = 6.852e10 ;
26double EvtIncoherentMixing::_deltams = 20.e12 ;
27
28//=============================================================================
29// Standard constructor, initializes variables
30//=============================================================================
31EvtIncoherentMixing::EvtIncoherentMixing( ) {
32 _doB0Mixing = false ;
33 _doBsMixing = false ;
34 _dGammad = 0. ;
35 // dGammas corresponds to DeltaGamma / Gamma = 10 %
36 _dGammas = 6.852e10 ;
37 _deltamd = 0.502e12 ;
38 _deltams = 20.e12 ;
39 _enableFlip = false ;
40}
41//=============================================================================
42EvtIncoherentMixing::~EvtIncoherentMixing( )
43{
44}
45// ============================================================================
46void EvtIncoherentMixing::incoherentB0Mix( const EvtId id, double &t ,
47 int &mix )
48{
49 static EvtId B0 = EvtPDL::getId( "B0" ) ;
50 static EvtId B0B = EvtPDL::getId( "anti-B0" ) ;
51
52 if ( ( B0 != id ) && ( B0B != id ) ) {
53 report(ERROR,"EvtGen") << "Bad configuration in incoherentB0Mix"
54 << std::endl ;
55 ::abort() ;
56 }
57
58 double x = getdeltamd() * EvtPDL::getctau( B0 ) / EvtConst::c ;
59
60 double y = getdGammad() * ( EvtPDL::getctau( B0 ) / EvtConst::c ) / 2. ;
61
62 double fac = 1. ; // No CP violation
63
64 double mixprob = ( x*x + y*y ) / ( x*x + y*y + ( 1./fac ) *
65 ( 2. + x*x - y*y ) ) ;
66
67 int mixsign ;
68
69 // decide if state is mixed
70 mixsign = ( mixprob > EvtRandom::Flat( 0. , 1. ) ) ? -1 : 1 ;
71
72 double prob ;
73
74 do {
75 t = -log( EvtRandom::Flat() ) * EvtPDL::getctau( B0 ) / ( 1. - y ) ;
76 prob = ( 1. + exp( -2. * y * t / EvtPDL::getctau( B0 ) ) +
77 mixsign * 2. * exp( -y * t / EvtPDL::getctau( B0 ) ) *
78 cos( getdeltamd() * t / EvtConst::c ) ) / 2. ;
79 } while ( prob < 2. * EvtRandom::Flat() ) ;
80
81 mix = 0 ;
82 if ( mixsign == -1 ) mix = 1 ;
83
84 return ;
85}
86// ============================================================================
87void EvtIncoherentMixing::incoherentBsMix( const EvtId id, double &t ,
88 int &mix )
89{
90 static EvtId BS = EvtPDL::getId( "B_s0" ) ;
91 static EvtId BSB = EvtPDL::getId( "anti-B_s0" ) ;
92
93 if ( ( BS != id ) && ( BSB != id ) ) {
94 report(ERROR,"EvtGen") << "Bad configuration in incoherentBsMix"
95 << std::endl ;
96 ::abort() ;
97 }
98
99 double x = getdeltams() * EvtPDL::getctau( BS ) / EvtConst::c ;
100
101 double y = getdGammas() * ( EvtPDL::getctau( BS ) / EvtConst::c ) / 2. ;
102
103 double fac = 1. ; // No CP violation
104
105 double mixprob = ( x*x + y*y ) / ( x*x + y*y + ( 1./fac ) *
106 ( 2. + x*x - y*y ) ) ;
107
108 int mixsign ;
109
110 // decide if state is mixed
111 mixsign = ( mixprob > EvtRandom::Flat( 0. , 1. ) ) ? -1 : 1 ;
112
113 double prob ;
114
115 do {
116 t = -log( EvtRandom::Flat() ) * EvtPDL::getctau( BS ) / ( 1. - y ) ;
117 prob = ( 1. + exp( -2. * y * t / EvtPDL::getctau( BS ) ) +
118 mixsign * 2. * exp( -y * t / EvtPDL::getctau( BS ) ) *
119 cos( getdeltams() * t / EvtConst::c ) ) / 2. ;
120 } while ( prob < 2. * EvtRandom::Flat() ) ;
121
122 mix = 0 ;
123 if ( mixsign == -1 ) mix = 1 ;
124
125 return ;
126}
127
128// ========================================================================
129bool EvtIncoherentMixing::isBsMixed ( EvtParticle * p )
130{
131 if ( ! ( p->getParent() ) ) return false ;
132
133 static EvtId BS0=EvtPDL::getId("B_s0");
134 static EvtId BSB=EvtPDL::getId("anti-B_s0");
135
136 if ( ( p->getId() != BS0 ) && ( p->getId() != BSB ) ) return false ;
137
138 if ( ( p->getParent()->getId() == BS0 ) ||
139 ( p->getParent()->getId() == BSB ) ) return true ;
140
141 return false ;
142}
143
144// ========================================================================
145bool EvtIncoherentMixing::isB0Mixed ( EvtParticle * p )
146{
147 if ( ! ( p->getParent() ) ) return false ;
148
149 static EvtId B0 =EvtPDL::getId("B0");
150 static EvtId B0B=EvtPDL::getId("anti-B0");
151
152 if ( ( p->getId() != B0 ) && ( p->getId() != B0B ) ) return false ;
153
154 if ( ( p->getParent()->getId() == B0 ) ||
155 ( p->getParent()->getId() == B0B ) ) return true ;
156
157 return false ;
158}
159//============================================================================
160// Return the tag of the event (ie the anti-flavour of the produced
161// B meson). Flip the flavour of the event with probB probability
162//============================================================================
163void EvtIncoherentMixing::OtherB( EvtParticle * p ,
164 double & t ,
165 EvtId & otherb ,
166 double probB )
167{
168 //if(p->getId() == B0 || p->getId() == B0B)
169 //added by liming Zhang
170 enableFlip();
171 if ( ( isB0Mixed( p ) ) || ( isBsMixed( p ) ) ) {
172 p->getParent()->setLifetime() ;
173 t = p->getParent()->getLifetime() ;
174 }
175 else {
176 p->setLifetime() ;
177 t = p->getLifetime() ;
178 }
179
180 if ( flipIsEnabled() ) {
181 //std::cout << " liming << flipIsEnabled " << std::endl;
182 // Flip the flavour of the particle with probability probB
183 bool isFlipped = ( EvtRandom::Flat( 0. , 1. ) < probB ) ;
184
185 if ( isFlipped ) {
186 if ( ( isB0Mixed( p ) ) || ( isBsMixed( p ) ) ) {
187 p->getParent()
188 ->setId( EvtPDL::chargeConj( p->getParent()->getId() ) ) ;
189 p->setId( EvtPDL::chargeConj( p->getId() ) ) ;
190 }
191 else {
192 p->setId( EvtPDL::chargeConj( p->getId() ) ) ;
193 }
194 }
195 }
196
197 if ( ( isB0Mixed( p ) ) || ( isBsMixed( p ) ) ) {
198 // if B has mixed, tag flavour is charge conjugate of parent of B-meson
199 otherb = EvtPDL::chargeConj( p->getParent()->getId() ) ;
200 }
201 else {
202 // else it is opposite flavour than this B hadron
203 otherb = EvtPDL::chargeConj( p->getId() ) ;
204 }
205
206 return ;
207}
208//============================================================================
209// Return the tag of the event (ie the anti-flavour of the produced
210// B meson). No flip
211//============================================================================
212void EvtIncoherentMixing::OtherB( EvtParticle * p ,
213 double & t ,
214 EvtId & otherb )
215{
216 if ( ( isB0Mixed( p ) ) || ( isBsMixed( p ) ) ) {
217 p->getParent()->setLifetime() ;
218 t = p->getParent()->getLifetime() ;
219 }
220 else {
221 p->setLifetime() ;
222 t = p->getLifetime() ;
223 }
224
225 if ( ( isB0Mixed( p ) ) || ( isBsMixed( p ) ) ) {
226 // if B has mixed, tag flavour is charge conjugate of parent of B-meson
227 otherb = EvtPDL::chargeConj( p->getParent()->getId() ) ;
228 }
229 else {
230 // else it is opposite flavour than this B hadron
231 otherb = EvtPDL::chargeConj( p->getId() ) ;
232 }
233
234 return ;
235}
236
237
238// activate or desactivate the Bs mixing
239void EvtIncoherentMixing::setB0Mixing() { _doB0Mixing = true ; }
240void EvtIncoherentMixing::unsetB0Mixing() { _doB0Mixing = false ; }
241
242// activate or desactivate the B0 mixing
243void EvtIncoherentMixing::setBsMixing() { _doBsMixing = true ; }
244void EvtIncoherentMixing::unsetBsMixing() { _doBsMixing = false ; }
245
246// is mixing activated ?
247bool EvtIncoherentMixing::doB0Mixing() { return _doB0Mixing ; }
248bool EvtIncoherentMixing::doBsMixing() { return _doBsMixing ; }
249
250// set values for the mixing
251void EvtIncoherentMixing::setdGammad( double value ) { _dGammad = value ; }
252void EvtIncoherentMixing::setdeltamd( double value ) { _deltamd = value ; }
253void EvtIncoherentMixing::setdGammas( double value ) { _dGammas = value ; }
254void EvtIncoherentMixing::setdeltams( double value ) { _deltams = value ; }
255
256// get parameters for mixing
257double EvtIncoherentMixing::getdGammad() { return _dGammad ; }
258double EvtIncoherentMixing::getdeltamd() { return _deltamd ; }
259double EvtIncoherentMixing::getdGammas() { return _dGammas ; }
260double EvtIncoherentMixing::getdeltams() { return _deltams ; }
261
262bool EvtIncoherentMixing::flipIsEnabled() { return _enableFlip ; }
263void EvtIncoherentMixing::enableFlip() { _enableFlip = true ; }
264void EvtIncoherentMixing::disableFlip() { _enableFlip = false ; }