TEvtGen/EvtGenBase/EvtSemiLeptonicVectorAmp.cxx

   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: EvtSemiLeptonicVectorAmp.cc
  12 //
  13 // Description: Routine to implement semileptonic decays to vector
  14 //              mesons.
  15 //
  16 // Modification history:
  17 //
  18 //    DJL       April 17,1998       Module created
  19 //
  20 //------------------------------------------------------------------------
  21 //
  22 #include "EvtGenBase/EvtPatches.hh"
  23 #include "EvtGenBase/EvtPatches.hh"
  24 #include "EvtGenBase/EvtParticle.hh"
  25 #include "EvtGenBase/EvtGenKine.hh"
  26 #include "EvtGenBase/EvtPDL.hh"
  27 #include "EvtGenBase/EvtReport.hh"
  28 #include "EvtGenBase/EvtTensor4C.hh"
  29 #include "EvtGenBase/EvtVector4C.hh"
  30 #include "EvtGenBase/EvtDiracSpinor.hh"
  31 #include "EvtGenBase/EvtSemiLeptonicVectorAmp.hh"
  32 #include "EvtGenBase/EvtId.hh"
  33 #include "EvtGenBase/EvtAmp.hh"
  34 #include "EvtGenBase/EvtSemiLeptonicFF.hh"
  35 using std::endl;
  36
  37 void EvtSemiLeptonicVectorAmp::CalcAmp( EvtParticle *parent,
  38                                         EvtAmp& amp,
  39                                         EvtSemiLeptonicFF *FormFactors ) {
  40
  41   static EvtId EM=EvtPDL::getId("e-");
  42   static EvtId MUM=EvtPDL::getId("mu-");
  43   static EvtId TAUM=EvtPDL::getId("tau-");
  44   static EvtId EP=EvtPDL::getId("e+");
  45   static EvtId MUP=EvtPDL::getId("mu+");
  46   static EvtId TAUP=EvtPDL::getId("tau+");
  47
  48   static EvtId D0=EvtPDL::getId("D0");
  49   static EvtId D0B=EvtPDL::getId("anti-D0");
  50   static EvtId DP=EvtPDL::getId("D+");
  51   static EvtId DM=EvtPDL::getId("D-");
  52   static EvtId DSM=EvtPDL::getId("D_s-");
  53   static EvtId DSP=EvtPDL::getId("D_s+");
  54
  55   //Add the lepton and neutrino 4 momenta to find q2
  56
  57   EvtVector4R q = parent->getDaug(1)->getP4()
  58                     + parent->getDaug(2)->getP4();
  59   double q2 = (q.mass2());
  60
  61   double a1f,a2f,vf,a0f,a3f;
  62   double m_meson = parent->getDaug(0)->mass();
  63
  64   FormFactors->getvectorff(parent->getId(),
  65                            parent->getDaug(0)->getId(),
  66                            q2,
  67                            m_meson,
  68                            &a1f,
  69                            &a2f,
  70                            &vf,
  71                            &a0f);
  72
  73   double costhl_flag = 1.0;
  74
  75   if(parent->getId()==D0||parent->getId()==D0B||
  76      parent->getId()==DP||parent->getId()==DM) {
  77     costhl_flag = -1.0;
  78   }
  79   if(parent->getId()==DSP||parent->getId()==DSM) {
  80     costhl_flag = -1.0;
  81   }
  82   vf = vf * costhl_flag;
  83
  84   EvtVector4R p4b;
  85   p4b.set(parent->mass(),0.0,0.0,0.0);
  86
  87   EvtVector4R p4meson = parent->getDaug(0)->getP4();
  88
  89   EvtVector4C l1,l2;
  90
  91   EvtId l_num = parent->getDaug(1)->getId();
  92   double m_b = parent->mass();
  93
  94   a3f = ((m_b+m_meson)/(2.0*m_meson))*a1f -
  95         ((m_b-m_meson)/(2.0*m_meson))*a2f;
  96
  97   EvtTensor4C tds;
  98   if (l_num==EM||l_num==MUM||l_num==TAUM){
  99
 100     tds = a1f*(m_b+m_meson)*EvtTensor4C::g();
 101     tds.addDirProd((-a2f/(m_b+m_meson))*p4b,p4b+p4meson);
 102     tds+=EvtComplex(0.0,vf/(m_b+m_meson))
 103                   *dual(directProd(p4meson+p4b,p4b-p4meson));
 104     tds.addDirProd((a0f-a3f)*2.0*(m_meson/q2)*p4b,p4b-p4meson);
 105
 106     l1=EvtLeptonVACurrent(parent->getDaug(1)->spParent(0),
 107                           parent->getDaug(2)->spParentNeutrino());
 108     l2=EvtLeptonVACurrent(parent->getDaug(1)->spParent(1),
 109                           parent->getDaug(2)->spParentNeutrino());
 110   }
 111   else{
 112     if (l_num==EP||l_num==MUP||l_num==TAUP){
 113     tds = a1f*(m_b+m_meson)*EvtTensor4C::g();
 114     tds.addDirProd((-a2f/(m_b+m_meson))*p4b,p4b+p4meson);
 115     tds-=EvtComplex(0.0,vf/(m_b+m_meson))
 116                   *dual(directProd(p4meson+p4b,p4b-p4meson));
 117     tds.addDirProd((a0f-a3f)*2.0*(m_meson/q2)*p4b,p4b-p4meson);
 118
 119     l1=EvtLeptonVACurrent(parent->getDaug(2)->spParentNeutrino(),
 120                             parent->getDaug(1)->spParent(0));
 121     l2=EvtLeptonVACurrent(parent->getDaug(2)->spParentNeutrino(),
 122                             parent->getDaug(1)->spParent(1));
 123     }
 124     else{
 125       report(ERROR,"EvtGen") << "Wrong lepton number"<<endl;
 126     }
 127   }
 128
 129   EvtVector4C et0=tds.cont1( parent->getDaug(0)->epsParent(0).conj() );
 130   EvtVector4C et1=tds.cont1( parent->getDaug(0)->epsParent(1).conj() );
 131   EvtVector4C et2=tds.cont1( parent->getDaug(0)->epsParent(2).conj() );
 132
 133
 134   amp.vertex(0,0,l1.cont(et0));
 135   amp.vertex(0,1,l2.cont(et0));
 136
 137   amp.vertex(1,0,l1.cont(et1));
 138   amp.vertex(1,1,l2.cont(et1));
 139
 140   amp.vertex(2,0,l1.cont(et2));
 141   amp.vertex(2,1,l2.cont(et2));
 142
 143   return;
 144 }
 145