Converting TEvtGen to native cmake

[u/mrichter/AliRoot.git] / TEvtGen / EvtGen / EvtGenModels / EvtbTosllScalarAmp.cpp

diff --git a/TEvtGen/EvtGen/EvtGenModels/EvtbTosllScalarAmp.cpp b/TEvtGen/EvtGen/EvtGenModels/EvtbTosllScalarAmp.cpp
new file mode 100644 (file)
index 0000000..52294f9
--- /dev/null
+++ b/TEvtGen/EvtGen/EvtGenModels/EvtbTosllScalarAmp.cpp
@@ -0,0 +1,183 @@
+//--------------------------------------------------------------------------
+//
+// Environment:
+//      This software is part of the EvtGen package developed jointly
+//      for the BaBar and CLEO collaborations.  If you use all or part
+//      of it, please give an appropriate acknowledgement.
+//
+// Copyright Information: See EvtGen/COPYRIGHT
+//      Copyright (C) 2000      Caltech, UCSB
+//
+// Module: EvtbTosllScalarAmp.cc
+//
+// Description: Routine to implement bTosll decays to pseudo-scalar
+//              mesons. 
+//
+// Modification history:
+//
+//    Ryd       January 5,2000       Module created
+//
+//------------------------------------------------------------------------
+//
+#include "EvtGenBase/EvtPatches.hh"
+#include "EvtGenBase/EvtConst.hh"
+#include "EvtGenBase/EvtPatches.hh"
+#include "EvtGenBase/EvtParticle.hh"
+#include "EvtGenBase/EvtGenKine.hh"
+#include "EvtGenBase/EvtPDL.hh"
+#include "EvtGenBase/EvtReport.hh"
+#include "EvtGenBase/EvtVector4C.hh"
+#include "EvtGenBase/EvtTensor4C.hh"
+#include "EvtGenBase/EvtDiracSpinor.hh"
+#include "EvtGenModels/EvtbTosllScalarAmp.hh"
+#include "EvtGenBase/EvtId.hh"
+#include "EvtGenBase/EvtIdSet.hh"
+#include "EvtGenBase/EvtAmp.hh"
+#include "EvtGenModels/EvtbTosllAmp.hh"
+#include "EvtGenModels/EvtbTosllFF.hh"
+
+void EvtbTosllScalarAmp::CalcAmp( EvtParticle *parent,
+                                 EvtAmp& amp,
+                                 EvtbTosllFF *formFactors ) {
+
+  //Add the lepton and neutrino 4 momenta to find q2
+
+  EvtVector4R q = parent->getDaug(1)->getP4() 
+                    + parent->getDaug(2)->getP4();
+  double q2 = (q.mass2());
+
+  double fp(0.),f0(0.),ft(0.);
+  double mesonmass = parent->getDaug(0)->mass();
+  double parentmass = parent->mass();
+
+  formFactors->getScalarFF(parent->getId(),
+                           parent->getDaug(0)->getId(),
+                           q2,
+                           mesonmass,
+                           fp,f0,ft);
+
+  EvtId daught = parent->getDaug(0)->getId();
+  bool btod = false;
+  bool nnlo = true;
+  if (daught == EvtPDL::getId(std::string("pi+")) ||
+         daught == EvtPDL::getId(std::string("pi-")) ||
+          daught == EvtPDL::getId(std::string("pi0")) ||
+          daught == EvtPDL::getId(std::string("eta")) ||
+          daught == EvtPDL::getId(std::string("eta'"))
+         )  btod = true;
+
+
+  EvtVector4R p4b;
+  p4b.set(parent->mass(),0.0,0.0,0.0);
+ 
+  EvtVector4C l11,l12;
+  EvtVector4C l21,l22;
+
+  EvtVector4C a11,a12;
+  EvtVector4C a21,a22;
+
+  EvtId l_num = parent->getDaug(1)->getId();
+
+  EvtVector4C T1,T2;
+
+  EvtVector4R phat=p4b/parentmass;
+  EvtVector4R qhat=q/parentmass;
+
+  EvtComplex c7eff = EvtbTosllAmp::GetC7Eff(q2,nnlo);
+  EvtComplex c9eff = EvtbTosllAmp::GetC9Eff(q2,nnlo,btod);
+  EvtComplex c10eff = EvtbTosllAmp::GetC10Eff(q2,nnlo);
+
+  //double mbhat=1;
+  double mbhat=4.4/(parentmass);      
+  //double mkhat = 0.15;
+  double mkhat=mesonmass/(parentmass);
+  double shat=q2/(parentmass*parentmass);
+
+  double fm=(f0-fp)*(1-mkhat*mkhat)/shat;
+
+  EvtComplex aprime;
+  aprime = c9eff*fp+2.0*mbhat*c7eff*ft/(1+mkhat);
+  EvtComplex bprime;
+  bprime = c9eff*fm-2*mbhat*c7eff*ft*(1-mkhat)/shat;
+
+  EvtComplex cprime;
+  cprime = c10eff*fp;
+  EvtComplex dprime; 
+  dprime = c10eff*fm;
+
+  static EvtIdSet leptons("e-","mu-","tau-");
+  static EvtIdSet antileptons("e+","mu+","tau+");
+
+  if (leptons.contains(l_num)){
+
+    T1=aprime*phat+bprime*qhat;
+    T2=cprime*phat+dprime*qhat;
+
+    l11=EvtLeptonVCurrent(parent->getDaug(1)->spParent(0),
+                          parent->getDaug(2)->spParent(0));
+    l21=EvtLeptonVCurrent(parent->getDaug(1)->spParent(1),
+                          parent->getDaug(2)->spParent(0));
+    l12=EvtLeptonVCurrent(parent->getDaug(1)->spParent(0),
+                          parent->getDaug(2)->spParent(1));
+    l22=EvtLeptonVCurrent(parent->getDaug(1)->spParent(1),
+                          parent->getDaug(2)->spParent(1));
+    a11=EvtLeptonACurrent(parent->getDaug(1)->spParent(0),
+                          parent->getDaug(2)->spParent(0));
+    a21=EvtLeptonACurrent(parent->getDaug(1)->spParent(1),
+                          parent->getDaug(2)->spParent(0));
+    a12=EvtLeptonACurrent(parent->getDaug(1)->spParent(0),
+                          parent->getDaug(2)->spParent(1));
+    a22=EvtLeptonACurrent(parent->getDaug(1)->spParent(1),
+                          parent->getDaug(2)->spParent(1));
+  }
+  else{
+    if (antileptons.contains(l_num)){
+
+    T1=aprime*phat+bprime*qhat;
+    T2=cprime*phat+dprime*qhat;
+
+    l11=EvtLeptonVCurrent(parent->getDaug(1)->spParent(1),
+                          parent->getDaug(2)->spParent(1));
+    l21=EvtLeptonVCurrent(parent->getDaug(1)->spParent(0),
+                          parent->getDaug(2)->spParent(1));
+    l12=EvtLeptonVCurrent(parent->getDaug(1)->spParent(1),
+                          parent->getDaug(2)->spParent(0));
+    l22=EvtLeptonVCurrent(parent->getDaug(1)->spParent(0),
+                          parent->getDaug(2)->spParent(0));
+
+    a11=EvtLeptonACurrent(parent->getDaug(1)->spParent(1),
+                          parent->getDaug(2)->spParent(1));
+    a21=EvtLeptonACurrent(parent->getDaug(1)->spParent(0),
+                          parent->getDaug(2)->spParent(1));
+    a12=EvtLeptonACurrent(parent->getDaug(1)->spParent(1),
+                          parent->getDaug(2)->spParent(0));
+    a22=EvtLeptonACurrent(parent->getDaug(1)->spParent(0),
+                          parent->getDaug(2)->spParent(0));
+
+    }
+    else{
+      report(ERROR,"EvtGen") << "Wrong lepton number\n";
+    }
+  }
+
+  amp.vertex(0,0,l11*T1+a11*T2);
+  amp.vertex(0,1,l12*T1+a12*T2);
+  amp.vertex(1,0,l21*T1+a21*T2);
+  amp.vertex(1,1,l22*T1+a22*T2);
+
+}
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+