Updates EvtGen Code

[u/mrichter/AliRoot.git] / TEvtGen / EvtGenModels / EvtLb2Lll.cpp

diff --git a/TEvtGen/EvtGenModels/EvtLb2Lll.cxx b/TEvtGen/EvtGenModels/EvtLb2Lll.cpp
similarity index 96%
rename from TEvtGen/EvtGenModels/EvtLb2Lll.cxx
rename to TEvtGen/EvtGenModels/EvtLb2Lll.cpp
index 9565db82a5fd665db6c768ee507c318b69408078..bd4a596c8c3ea6d403be4dc725f422bd8fc36c91 100644 (file)
--- a/TEvtGen/EvtGenModels/EvtLb2Lll.cxx
+++ b/TEvtGen/EvtGenModels/EvtLb2Lll.cpp
@@ -8,11 +8,16 @@
 //
 // Modification history:
 //
+//  10/07/2012  MK   Fix calculation of N1, N2; based on hep-ph/021144
 //  09/02/2009  PR   Commented check for (anti-)Lambda0 names
 //  15/09/2004  PR   Module created according to PHSP model
 //  20/02/2005  PR   Added parameters, created matrix element (without polarization)
 //  04/03/2005  PR   LD contrib., corrected WC eff. according to Chen. Geng.
 //
+// Todo list:
+//
+//   - Properly handle antiparticles, needs change of u, ubar to v, vbar in
+//     hadronic current, or other way of putting that in
 //----------------------------------------------------------------------------------
 
 #ifdef WIN32
@@ -193,9 +198,9 @@ void EvtLb2Lll::init(){
   report(INFO,"EvtGen") << " EvtLb2Lll maximum probability was set to " << m_maxProbability << std::endl;
   m_poleSize=0;
 
-  // Initialize Wilson coeficients by Buras and Munz
+  // Initialize Wilson coefficients by Buras and Munz
   // TODO: should have common W.C. source for all decays in EvtGen
-  m_WC.CalculateAllCoeficients();
+  m_WC.CalculateAllCoefficients();
 
 }
 
@@ -298,10 +303,10 @@ void EvtLb2Lll::calcAmp(EvtAmp *amp,EvtParticle *parent){
   int i,j,spins[4];
   char ch;
 
-  double r,M_L,M_Lb,M_s,M_c,M_b,q2,alpha,alpha_s,eta,M_W,M_t;
+  double r,M_L,M_Lb,M_s,M_c,M_b,q2,alpha,M_W,M_t;
   double M_psi[2]={0,0},Gamma_psi[2]={0,0},k_psi[2]={0,0};
   double F0_1,F0_2,a_F1,a_F2,b_F1,b_F2,F1,F2;
-  double f_1,f_2,f_3,g_1,g_2,g_3,f_1T,f_2T,f_3T,g_1T,g_2T,g_3T,f_TV,f_TS,g_TV,g_TS,f_T,g_T;
+  double f_1,f_2,f_3,g_1,g_2,g_3,f_1T,f_2T,f_3T,g_1T,g_2T,g_3T,f_TV,f_TS,g_TV(0.0),g_TS,f_T,g_T;
   EvtComplex A1,A2,A3,B1,B2,B3,D1,D2,D3,E1,E2,E3,N1,N2,H1,H2;
   EvtComplex C_SL,C_BR,C_LLtot,C_LRtot,C_LL,C_LR,C_RL,C_RR,C_LRLR,C_RLLR,C_LRRL,C_RLRL,C_T,C_TE;
   EvtComplex Yld,C_7eff,C_9eff;
@@ -320,8 +325,6 @@ void EvtLb2Lll::calcAmp(EvtAmp *amp,EvtParticle *parent){
   M_c     = 1.35;
   M_b     = 4.8;
   alpha   = 1./137.036;
-  alpha_s = 0.217484;
-  eta     = 0.556289;
   M_W     = 80.425;
   M_t     = 174.3;
   M_psi[0] = 3.096916;
@@ -481,12 +484,13 @@ void EvtLb2Lll::calcAmp(EvtAmp *amp,EvtParticle *parent){
   E2 = 0.5*(C_RR-C_RL)*(f_2-g_2) + 0.5*(C_LRtot-C_LLtot)*(f_2+g_2);
   E3 = 0.5*(C_RR-C_RL)*(f_3-g_3) + 0.5*(C_LRtot-C_LLtot)*(f_3+g_3);
 
-  N1 = (f_1*(M_Lb+M_L)+f_3*q2)/M_b*(C_LRLR+C_RLLR+C_LRRL+C_RLRL);
-  N2 = (f_1*(M_Lb+M_L)+f_3*q2)/M_b*(C_LRLR+C_RLLR-C_LRRL-C_RLRL);
+  N1 = (f_1*(M_Lb-M_L)+f_3*q2)/M_b*(C_LRLR+C_RLLR+C_LRRL+C_RLRL);  // Should be mLb - mL
+  N2 = (f_1*(M_Lb-M_L)+f_3*q2)/M_b*(C_LRLR+C_RLLR-C_LRRL-C_RLRL);
 
   H1 = (g_1*(M_Lb+M_L)-g_3*q2)/M_b*(C_LRLR-C_RLLR+C_LRRL-C_RLRL);
   H2 = (g_1*(M_Lb+M_L)-g_3*q2)/M_b*(C_LRLR-C_RLLR-C_LRRL+C_RLRL);
 
+
   for(i=0;i<4;i++){
     lbar_Gmu_l   [i/2][i%2] = EvtLeptonVCurrent(parent->getDaug(1)->spParent(i/2),parent->getDaug(2)->spParent(i%2));
     lbar_GmuG5_l [i/2][i%2] = EvtLeptonACurrent(parent->getDaug(1)->spParent(i/2),parent->getDaug(2)->spParent(i%2));
@@ -549,15 +553,18 @@ void EvtLb2Lll::calcAmp(EvtAmp *amp,EvtParticle *parent){
     //std::cout << "Matrix = " << Matrix[j/2][j%2][i/2][i%2] << std::endl;
     Matrix[j/2][j%2][i/2][i%2] += cont(lbar_Smunu_l[i/2][i%2]  , 4*C_T*f_T*hbar_Smunu_h[j/2][j%2]);
     //std::cout << "Matrix = " << Matrix[j/2][j%2][i/2][i%2] << std::endl;
-    Matrix[j/2][j%2][i/2][i%2] += cont(lbar_Smunu_l[i/2][i%2]  , -4*C_T*f_TV*i1*(directProd(q_mu,hbar_Gmu_h[j/2][j%2])-directProd(hbar_Gmu_h[j/2][j%2],q_mu)));
+    Matrix[j/2][j%2][i/2][i%2] += 
+      cont(lbar_Smunu_l[i/2][i%2]  , 
+           -4*C_T*f_TV*i1*(EvtGenFunctions::directProd(q_mu,hbar_Gmu_h[j/2][j%2])-
+                           EvtGenFunctions::directProd(hbar_Gmu_h[j/2][j%2],q_mu)));
     //std::cout << "Matrix = " << Matrix[j/2][j%2][i/2][i%2] << std::endl;
-    Matrix[j/2][j%2][i/2][i%2] += cont(lbar_Smunu_l[i/2][i%2]  , -4*C_T*f_TS*i1*(directProd(P_mu,q_mu)-directProd(q_mu,P_mu))*hbar_h[j/2][j%2]);
+    Matrix[j/2][j%2][i/2][i%2] += cont(lbar_Smunu_l[i/2][i%2]  , -4*C_T*f_TS*i1*(EvtGenFunctions::directProd(P_mu,q_mu)-EvtGenFunctions::directProd(q_mu,P_mu))*hbar_h[j/2][j%2]);
     //std::cout << "Matrix = " << Matrix[j/2][j%2][i/2][i%2] << std::endl;
     Matrix[j/2][j%2][i/2][i%2] += cont(lbar_ESmunu_l[i/2][i%2] , 4*C_TE*f_T*i1*hbar_Smunu_h[j/2][j%2]);
     //std::cout << "Matrix = " << Matrix[j/2][j%2][i/2][i%2] << std::endl;
-    Matrix[j/2][j%2][i/2][i%2] += cont(lbar_ESmunu_l[i/2][i%2] , 4*C_TE*f_TV*(directProd(q_mu,hbar_Gmu_h[j/2][j%2])-directProd(hbar_Gmu_h[j/2][j%2],q_mu)));
+    Matrix[j/2][j%2][i/2][i%2] += cont(lbar_ESmunu_l[i/2][i%2] , 4*C_TE*f_TV*(EvtGenFunctions::directProd(q_mu,hbar_Gmu_h[j/2][j%2])-EvtGenFunctions::directProd(hbar_Gmu_h[j/2][j%2],q_mu)));
     //std::cout << "Matrix = " << Matrix[j/2][j%2][i/2][i%2] << std::endl;
-    Matrix[j/2][j%2][i/2][i%2] += cont(lbar_ESmunu_l[i/2][i%2] , 4*C_TE*f_TS*(directProd(P_mu,q_mu)-directProd(q_mu,P_mu))*hbar_h[j/2][j%2]);
+    Matrix[j/2][j%2][i/2][i%2] += cont(lbar_ESmunu_l[i/2][i%2] , 4*C_TE*f_TS*(EvtGenFunctions::directProd(P_mu,q_mu)-EvtGenFunctions::directProd(q_mu,P_mu))*hbar_h[j/2][j%2]);
     //std::cout << "Matrix = " << Matrix[j/2][j%2][i/2][i%2] << std::endl;
     //Matrix[j/2][j%2][i/2][i%2] *= G_F*alpha/4/sqrt(2)/EvtConst::pi*V_tb*conj(V_ts);
     //std::cout << "Matrix = " << Matrix[j/2][j%2][i/2][i%2] << std::endl;