]>
Commit | Line | Data |
---|---|---|
0ca57c2f | 1 | //-------------------------------------------------------------------------- |
2 | // | |
3 | // Environment: | |
4 | // This software is part of the EvtGen package developed jointly | |
5 | // for the BaBar and CLEO collaborations. If you use all or part | |
6 | // of it, please give an appropriate acknowledgement. | |
7 | // | |
8 | // Copyright Information: See EvtGen/COPYRIGHT | |
9 | // Copyright (C) 1998 Caltech, UCSB | |
10 | // | |
11 | // Module: EvtTVP.cc | |
12 | // | |
13 | // Description: Routine to calculate W -> (n pi) current | |
14 | // according to [Kuhn, Was, Acta.Phys.Polon B39 (2008) 147] | |
15 | // | |
16 | // Modification history: | |
17 | // AVL 6 July, 2012 Module created | |
18 | // | |
19 | //------------------------------------------------------------------------ | |
20 | // | |
21 | #include "EvtGenBase/EvtPatches.hh" | |
22 | #include <iostream> | |
23 | #include <iomanip> | |
24 | #include <fstream> | |
25 | #include <ctype.h> | |
26 | #include <stdlib.h> | |
27 | #include <string.h> | |
28 | #include "EvtGenBase/EvtParticle.hh" | |
29 | #include "EvtGenBase/EvtPDL.hh" | |
30 | #include "EvtGenBase/EvtGenKine.hh" | |
31 | #include "EvtGenModels/EvtTauHadnu.hh" | |
32 | #include "EvtGenBase/EvtDiracSpinor.hh" | |
33 | #include "EvtGenBase/EvtReport.hh" | |
34 | #include "EvtGenBase/EvtVector4C.hh" | |
35 | #include "EvtGenBase/EvtTensor4C.hh" | |
36 | #include "EvtGenBase/EvtIdSet.hh" | |
37 | #include "EvtGenBase/EvtParser.hh" | |
38 | ||
39 | #include "EvtGenModels/EvtWnPi.hh" | |
40 | ||
41 | using namespace std; | |
42 | ||
43 | // W+ -> pi_ current | |
44 | EvtVector4C EvtWnPi::WCurrent(EvtVector4R q1) { | |
45 | return q1; | |
46 | } | |
47 | ||
48 | // W+ -> pi+ pi0 current | |
49 | EvtVector4C EvtWnPi::WCurrent(EvtVector4R q1, EvtVector4R q2) { | |
50 | return BWr(q1+q2)*(q1-q2); | |
51 | } | |
52 | ||
53 | // W+ -> pi+ pi+ pi- current | |
54 | EvtVector4C EvtWnPi::WCurrent(EvtVector4R q1, EvtVector4R q2, EvtVector4R q3) { | |
55 | EvtVector4R Q=q1+q2+q3; | |
56 | double Q2=Q.mass2(); | |
57 | return BWa(Q)*( (q1-q3) - (Q*(Q*(q1-q3))/Q2)*BWr(q2+q3) + | |
58 | (q2-q3) - (Q*(Q*(q2-q3))/Q2)*BWr(q1+q3) ); | |
59 | } | |
60 | ||
61 | // W+ -> pi+ pi+ pi- pi- pi+ current with symmetrization | |
62 | EvtVector4C EvtWnPi::WCurrent(EvtVector4R q1, EvtVector4R q2, EvtVector4R q3, EvtVector4R q4, EvtVector4R q5) { | |
63 | // double Q2 = Qtot*Qtot; | |
64 | // return q1-Qtot*(q1*Qtot)/Q2; | |
65 | EvtVector4C V = JB(q1, q2, q3, q4, q5) + JB(q5, q2, q3, q4, q1) + JB(q1, q5, q3, q4, q2) + | |
66 | JB(q1,q2,q4,q3,q5)+JB(q5,q2,q4,q3,q1)+JB(q1,q5,q4,q3,q2); | |
67 | // cout<<"BC2: Qtot="<<Qtot<<", V="<<V<<endl; | |
68 | return V; | |
69 | } | |
70 | ||
71 | ||
72 | // a1 -> pi+ pi+ pi- BW | |
73 | EvtComplex EvtWnPi::BWa(EvtVector4R q) { | |
74 | double const _mA1=1.26, _GA1=0.4; | |
75 | EvtComplex I(0,1); | |
76 | double Q2 = q.mass2(); | |
77 | double GA1=_GA1*pi3G(Q2)/pi3G(_mA1*_mA1); | |
78 | EvtComplex denBA1(_mA1*_mA1 - Q2,-1.*_mA1*GA1); | |
79 | return _mA1*_mA1 / denBA1; | |
80 | } | |
81 | ||
82 | ||
83 | EvtComplex EvtWnPi::BWf(EvtVector4R q) { | |
84 | double const mf=0.8, Gf=0.6; | |
85 | EvtComplex I(0,1); | |
86 | double Q2 = q.mass2(); | |
87 | return mf*mf/(mf*mf-Q2-I*mf*Gf); | |
88 | } | |
89 | ||
90 | EvtComplex EvtWnPi::BWr(EvtVector4R q) { | |
91 | double _mRho = 0.775, _gammaRho=0.149, _mRhopr=1.364, _gammaRhopr=0.400, _beta=-0.108; | |
92 | double m1=EvtPDL::getMeanMass(EvtPDL::getId("pi+")), m2=EvtPDL::getMeanMass(EvtPDL::getId("pi+")); | |
93 | double mQ2=q.mass2(); | |
94 | ||
95 | // momenta in the rho->pipi decay | |
96 | double dRho= _mRho*_mRho - m1*m1 - m2*m2; | |
97 | double pPiRho = (1.0/_mRho)*sqrt((dRho*dRho)/4.0 - m1*m1*m2*m2); | |
98 | ||
99 | double dRhopr= _mRhopr*_mRhopr - m1*m1 - m2*m2; | |
100 | double pPiRhopr = (1.0/_mRhopr)*sqrt((dRhopr*dRhopr)/4.0 - m1*m1*m2*m2); | |
101 | ||
102 | double dQ= mQ2 - m1*m1 - m2*m2; | |
103 | double pPiQ = (1.0/sqrt(mQ2))*sqrt((dQ*dQ)/4.0 - m1*m1*m2*m2); | |
104 | ||
105 | ||
106 | double gammaRho = _gammaRho*_mRho/sqrt(mQ2)*pow((pPiQ/pPiRho),3); | |
107 | EvtComplex BRhoDem(_mRho*_mRho - mQ2,-1.0*_mRho*gammaRho); | |
108 | EvtComplex BRho= _mRho*_mRho / BRhoDem; | |
109 | ||
110 | double gammaRhopr = _gammaRhopr*_mRhopr/sqrt(mQ2)*pow((pPiQ/pPiRhopr),3); | |
111 | EvtComplex BRhoprDem(_mRhopr*_mRhopr - mQ2,-1.0*_mRho*gammaRhopr); | |
112 | EvtComplex BRhopr= _mRhopr*_mRhopr / BRhoprDem; | |
113 | ||
114 | return (BRho + _beta*BRhopr)/(1+_beta); | |
115 | } | |
116 | ||
117 | double EvtWnPi::pi3G(double m2) { | |
118 | double mPi = EvtPDL::getMeanMass(EvtPDL::getId("pi+")); | |
119 | double _mRho = 0.775; | |
120 | if ( m2 > (_mRho+mPi) ) { | |
121 | return m2*(1.623 + 10.38/m2 - 9.32/(m2*m2) + 0.65/(m2*m2*m2)); | |
122 | } | |
123 | else { | |
124 | double t1=m2-9.0*mPi*mPi; | |
125 | return 4.1*pow(t1,3.0)*(1.0 - 3.3*t1+5.8*t1*t1); | |
126 | }; | |
127 | } | |
128 | ||
129 | EvtVector4C EvtWnPi::JB( EvtVector4R p1, EvtVector4R p2, EvtVector4R p3, EvtVector4R p4, EvtVector4R p5) { | |
130 | EvtVector4R Qtot = p1+p2+p3+p4+p5, Qa=p1+p2+p3; | |
131 | EvtTensor4C T= (1/Qtot.mass2())*EvtGenFunctions::directProd(Qtot,Qtot) - EvtTensor4C::g(); | |
132 | EvtVector4R V13 = Qa*( p2*(p1-p3) )/Qa.mass2() - (p1-p3); | |
133 | EvtVector4R V23 = Qa*( p1*(p2-p3) )/Qa.mass2() - (p2-p3); | |
134 | return BWa(Qtot)*BWa(Qa)*BWf(p4+p5)*( | |
135 | T.cont1(V13)*BWr(p1+p3) + T.cont1(V23)*BWr(p2+p3) | |
136 | ); | |
137 | } |