00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032 #include "EvtGenBase/EvtPatches.hh"
00033 #include <stdlib.h>
00034 #include <iostream>
00035 #include <string>
00036 #include "EvtGenBase/EvtParticle.hh"
00037 #include "EvtGenBase/EvtPDL.hh"
00038 #include "EvtGenBase/EvtGenKine.hh"
00039 #include "EvtGenBase/EvtDiracSpinor.hh"
00040 #include "EvtGenBase/EvtTensor4C.hh"
00041 #include "EvtGenModels/EvtKstarnunu.hh"
00042 #include "EvtGenBase/EvtReport.hh"
00043 #include "EvtGenBase/EvtVector4C.hh"
00044
00045 EvtKstarnunu::~EvtKstarnunu() {}
00046
00047 void EvtKstarnunu::getName(std::string& model_name){
00048
00049 model_name="KSTARNUNU";
00050
00051 }
00052
00053
00054 EvtDecayBase* EvtKstarnunu::clone(){
00055
00056 return new EvtKstarnunu;
00057
00058 }
00059
00060 void EvtKstarnunu::init(){
00061
00062
00063 checkNArg(0);
00064 checkNDaug(3);
00065
00066
00067
00068
00069 checkSpinParent(EvtSpinType::SCALAR);
00070
00071 checkSpinDaughter(0,EvtSpinType::VECTOR);
00072 checkSpinDaughter(1,EvtSpinType::NEUTRINO);
00073 checkSpinDaughter(2,EvtSpinType::NEUTRINO);
00074
00075 }
00076
00077
00078 void EvtKstarnunu::decay(EvtParticle *p){
00079
00080 static EvtId NUE=EvtPDL::getId("nu_e");
00081 static EvtId NUM=EvtPDL::getId("nu_mu");
00082 static EvtId NUT=EvtPDL::getId("nu_tau");
00083 static EvtId NUEB=EvtPDL::getId("anti-nu_e");
00084 static EvtId NUMB=EvtPDL::getId("anti-nu_mu");
00085 static EvtId NUTB=EvtPDL::getId("anti-nu_tau");
00086
00087 p->initializePhaseSpace(getNDaug(),getDaugs());
00088
00089 double m_b = p->mass();
00090
00091 EvtParticle *meson, *neutrino1, *neutrino2;
00092 meson = p->getDaug(0);
00093 neutrino1 = p->getDaug(1);
00094 neutrino2 = p->getDaug(2);
00095 EvtVector4R momnu1 = neutrino1->getP4();
00096 EvtVector4R momnu2 = neutrino2->getP4();
00097 EvtVector4R momkstar = meson->getP4();
00098
00099 double v0_0, a0_0, a1_0, a2_0;
00100 double m2v0, m2a0, a1_b, a2_b;
00101 v0_0 = 0.47;
00102 a0_0 = 0.30;
00103 a1_0 = 0.37;
00104 a2_0 = 0.40;
00105 m2v0 = 5.*5.;
00106 m2a0 = 4.8*4.8;
00107 a1_b = -0.023;
00108 a2_b = 0.034;
00109
00110 EvtVector4R q = momnu1+momnu2;
00111 double q2 = q.mass2();
00112
00113 double v0, a1, a2;
00114 v0 = v0_0/(1-q2/m2v0);
00115 a1 = a1_0*(1+a1_b*q2);
00116 a2 = a2_0*(1+a2_b*q2);
00117
00118 EvtVector4R p4b; p4b.set(m_b,0.,0.,0.);
00119
00120 double m_k = meson->mass();
00121
00122 EvtTensor4C tds=(-2*v0/(m_b+m_k))*dual(directProd(p4b,momkstar))
00123 - EvtComplex(0.0,1.0)*
00124 ( (m_b+m_k)*a1*EvtTensor4C::g()
00125 - (a2/(m_b+m_k))*directProd(p4b-momkstar,p4b+momkstar));
00126
00127 EvtVector4C l;
00128
00129 if (getDaug(1)==NUE||getDaug(1)==NUM||getDaug(1)==NUT) {
00130 l=EvtLeptonVACurrent(neutrino1->spParentNeutrino(),
00131 neutrino2->spParentNeutrino());
00132 }
00133 if (getDaug(1)==NUEB||getDaug(1)==NUMB||getDaug(1)==NUTB) {
00134 l=EvtLeptonVACurrent(neutrino2->spParentNeutrino(),
00135 neutrino1->spParentNeutrino());
00136 }
00137
00138 EvtVector4C et0,et1,et2;
00139 et0 = tds.cont1( meson->epsParent(0).conj() );
00140 et1 = tds.cont1( meson->epsParent(1).conj() );
00141 et2 = tds.cont1( meson->epsParent(2).conj() );
00142
00143 vertex(0,l*et0);
00144 vertex(1,l*et1);
00145 vertex(2,l*et2);
00146
00147 return;
00148 }
00149