TEvtGen/EvtGenBase/EvtSemiLeptonicTensorAmp.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: EvtSemiLeptonicTensorAmp.cc
  12 //
  13 // Description: Routine to implement semileptonic decays to pseudo-scalar
  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/EvtSemiLeptonicTensorAmp.hh"
  32 #include "EvtGenBase/EvtId.hh"
  33 #include "EvtGenBase/EvtAmp.hh"
  34 #include "EvtGenBase/EvtSemiLeptonicFF.hh"
  35
  36 void EvtSemiLeptonicTensorAmp::CalcAmp( EvtParticle *parent,
  37                                         EvtAmp& amp,
  38                                         EvtSemiLeptonicFF *FormFactors ) {
  39   static EvtId EM=EvtPDL::getId("e-");
  40   static EvtId MUM=EvtPDL::getId("mu-");
  41   static EvtId TAUM=EvtPDL::getId("tau-");
  42   static EvtId EP=EvtPDL::getId("e+");
  43   static EvtId MUP=EvtPDL::getId("mu+");
  44   static EvtId TAUP=EvtPDL::getId("tau+");
  45
  46   static EvtId D0=EvtPDL::getId("D0");
  47   static EvtId D0B=EvtPDL::getId("anti-D0");
  48   static EvtId DP=EvtPDL::getId("D+");
  49   static EvtId DM=EvtPDL::getId("D-");
  50   static EvtId DSM=EvtPDL::getId("D_s-");
  51   static EvtId DSP=EvtPDL::getId("D_s+");
  52
  53   //Add the lepton and neutrino 4 momenta to find q2
  54
  55   EvtVector4R q = parent->getDaug(1)->getP4()
  56                     + parent->getDaug(2)->getP4();
  57   double q2 = (q.mass2());
  58
  59   double hf,kf,bpf,bmf;
  60
  61   FormFactors->gettensorff(parent->getId(),
  62                            parent->getDaug(0)->getId(),
  63                            q2,
  64                            parent->getDaug(0)->mass(),
  65                            &hf,
  66                            &kf,
  67                            &bpf,
  68                            &bmf);
  69
  70
  71   double costhl_flag = 1.0;
  72
  73   if(parent->getId()==D0||parent->getId()==D0B||
  74      parent->getId()==DP||parent->getId()==DM) {
  75     costhl_flag = -1.0;
  76   }
  77   if(parent->getId()==DSP||parent->getId()==DSM) {
  78     costhl_flag = -1.0;
  79   }
  80   hf = hf * costhl_flag;
  81
  82   EvtVector4R p4b;
  83   p4b.set(parent->mass(),0.0,0.0,0.0);
  84
  85   EvtVector4R p4meson = parent->getDaug(0)->getP4();
  86
  87   EvtVector4C l1,l2;
  88
  89   EvtId l_num = parent->getDaug(1)->getId();
  90
  91   EvtVector4C ep_meson_b[5];
  92
  93   ep_meson_b[0] = ((parent->getDaug(0)->epsTensorParent(0)).cont2(p4b)).conj();
  94   ep_meson_b[1] = ((parent->getDaug(0)->epsTensorParent(1)).cont2(p4b)).conj();
  95   ep_meson_b[2] = ((parent->getDaug(0)->epsTensorParent(2)).cont2(p4b)).conj();
  96   ep_meson_b[3] = ((parent->getDaug(0)->epsTensorParent(3)).cont2(p4b)).conj();
  97   ep_meson_b[4] = ((parent->getDaug(0)->epsTensorParent(4)).cont2(p4b)).conj();
  98
  99   EvtVector4R pp,pm;
 100
 101   pp=p4b+p4meson;
 102   pm=p4b-p4meson;
 103
 104   //lange - October 31,2002 - try to lessen the mass dependence of probmax
 105   double q2max = p4b.mass2() + p4meson.mass2() - 2.0*p4b.mass()*p4meson.mass();
 106   double q2maxin=1.0/q2max;
 107
 108   EvtComplex ep_meson_bb[5];
 109
 110   ep_meson_bb[0]=ep_meson_b[0]*(p4b);
 111   ep_meson_bb[1]=ep_meson_b[1]*(p4b);
 112   ep_meson_bb[2]=ep_meson_b[2]*(p4b);
 113   ep_meson_bb[3]=ep_meson_b[3]*(p4b);
 114   ep_meson_bb[4]=ep_meson_b[4]*(p4b);
 115
 116
 117   EvtVector4C tds0,tds1,tds2,tds3,tds4;
 118
 119   EvtTensor4C tds;
 120   if (l_num==EM||l_num==MUM||l_num==TAUM){
 121     EvtTensor4C tdual=EvtComplex(0.0,hf)*dual(directProd(pp,pm));
 122     tds0=tdual.cont2(ep_meson_b[0])
 123       -kf*ep_meson_b[0]
 124       -bpf*ep_meson_bb[0]*pp-bmf*ep_meson_bb[0]*pm;
 125     tds0*=q2maxin;
 126
 127     tds1=tdual.cont2(ep_meson_b[1])
 128       -kf*ep_meson_b[1]
 129       -bpf*ep_meson_bb[1]*pp-bmf*ep_meson_bb[1]*pm;
 130     tds1*=q2maxin;
 131
 132     tds2=tdual.cont2(ep_meson_b[2])
 133       -kf*ep_meson_b[2]
 134       -bpf*ep_meson_bb[2]*pp-bmf*ep_meson_bb[2]*pm;
 135     tds2*=q2maxin;
 136
 137     tds3=tdual.cont2(ep_meson_b[3])
 138       -kf*ep_meson_b[3]
 139       -bpf*ep_meson_bb[3]*pp-bmf*ep_meson_bb[3]*pm;
 140     tds3*=q2maxin;
 141
 142     tds4=tdual.cont2(ep_meson_b[4])
 143       -kf*ep_meson_b[4]
 144       -bpf*ep_meson_bb[4]*pp-bmf*ep_meson_bb[4]*pm;
 145     tds4*=q2maxin;
 146
 147
 148     l1=EvtLeptonVACurrent(parent->getDaug(1)->spParent(0),
 149                           parent->getDaug(2)->spParentNeutrino());
 150     l2=EvtLeptonVACurrent(parent->getDaug(1)->spParent(1),
 151                           parent->getDaug(2)->spParentNeutrino());
 152   }
 153   else{
 154     if (l_num==EP||l_num==MUP||l_num==TAUP){
 155      EvtTensor4C tdual=EvtComplex(0.0,-hf)*dual(directProd(pp,pm));
 156       tds0=tdual.cont2(ep_meson_b[0])
 157         -kf*ep_meson_b[0]
 158         -bpf*ep_meson_bb[0]*pp-bmf*ep_meson_bb[0]*pm;
 159       tds0*=q2maxin;
 160
 161       tds1=tdual.cont2(ep_meson_b[1])
 162         -kf*ep_meson_b[1]
 163         -bpf*ep_meson_bb[1]*pp-bmf*ep_meson_bb[1]*pm;
 164       tds1*=q2maxin;
 165
 166       tds2=tdual.cont2(ep_meson_b[2])
 167         -kf*ep_meson_b[2]
 168         -bpf*ep_meson_bb[2]*pp-bmf*ep_meson_bb[2]*pm;
 169       tds2*=q2maxin;
 170
 171       tds3=tdual.cont2(ep_meson_b[3])
 172         -kf*ep_meson_b[3]
 173         -bpf*ep_meson_bb[3]*pp-bmf*ep_meson_bb[3]*pm;
 174       tds3*=q2maxin;
 175
 176       tds4=tdual.cont2(ep_meson_b[4])
 177         -kf*ep_meson_b[4]
 178         -bpf*ep_meson_bb[4]*pp-bmf*ep_meson_bb[4]*pm;
 179       tds4*=q2maxin;
 180
 181       l1=EvtLeptonVACurrent(parent->getDaug(2)->spParentNeutrino(),
 182                             parent->getDaug(1)->spParent(0));
 183       l2=EvtLeptonVACurrent(parent->getDaug(2)->spParentNeutrino(),
 184                             parent->getDaug(1)->spParent(1));
 185     }
 186     else{
 187       report(ERROR,"EvtGen") << "dfnb89agngri wrong lepton number\n";
 188     }
 189   }
 190
 191   amp.vertex(0,0,l1*tds0);
 192   amp.vertex(0,1,l2*tds0);
 193
 194   amp.vertex(1,0,l1*tds1);
 195   amp.vertex(1,1,l2*tds1);
 196
 197   amp.vertex(2,0,l1*tds2);
 198   amp.vertex(2,1,l2*tds2);
 199
 200   amp.vertex(3,0,l1*tds3);
 201   amp.vertex(3,1,l2*tds3);
 202
 203   amp.vertex(4,0,l1*tds4);
 204   amp.vertex(4,1,l2*tds4);
 205
 206   return;
 207
 208 }
 209