Converting TEvtGen to native cmake
[u/mrichter/AliRoot.git] / TEvtGen / EvtGen / EvtGenModels / EvtKstarnunu.cpp
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: EvtKstarnunu.cc
12 //
13 // Description: B ==> K* nu nubar
14 //
15 // Modification history:
16 //
17 //    RYD     September 29, 1997        Module created
18 //
19 //------------------------------------------------------------------------
20 //
21 // Package:     EVT
22 // Module:      EvtKstarnunu
23 // 
24 //              
25 // Implementation:
26 //     
27 //
28 // Author:      RYD 
29 // Created:     Sept. 29 1997
30 //
31 //
32 #include "EvtGenBase/EvtPatches.hh"
33 #include <stdlib.h>
34 #include <iostream>
35 #include <string>
36 #include "EvtGenBase/EvtParticle.hh"
37 #include "EvtGenBase/EvtPDL.hh"
38 #include "EvtGenBase/EvtGenKine.hh"
39 #include "EvtGenBase/EvtDiracSpinor.hh"
40 #include "EvtGenBase/EvtTensor4C.hh"
41 #include "EvtGenModels/EvtKstarnunu.hh"
42 #include "EvtGenBase/EvtReport.hh"
43 #include "EvtGenBase/EvtVector4C.hh"
44
45 EvtKstarnunu::~EvtKstarnunu() {}
46
47 std::string EvtKstarnunu::getName(){
48
49   return "KSTARNUNU";     
50
51 }
52
53
54 EvtDecayBase* EvtKstarnunu::clone(){
55
56   return new EvtKstarnunu;
57
58 }
59
60 void EvtKstarnunu::init(){
61
62   // check that there are 0 arguments
63   checkNArg(0);
64   checkNDaug(3);
65
66   //We expect the parent to be a scalar 
67   //and the daughters to be K neutrino netrino
68
69   checkSpinParent(EvtSpinType::SCALAR);
70
71   checkSpinDaughter(0,EvtSpinType::VECTOR);
72   checkSpinDaughter(1,EvtSpinType::NEUTRINO);
73   checkSpinDaughter(2,EvtSpinType::NEUTRINO);
74
75 }
76
77
78 void EvtKstarnunu::decay(EvtParticle *p){
79   
80   static EvtId NUE=EvtPDL::getId("nu_e");
81   static EvtId NUM=EvtPDL::getId("nu_mu");
82   static EvtId NUT=EvtPDL::getId("nu_tau");
83   static EvtId NUEB=EvtPDL::getId("anti-nu_e");
84   static EvtId NUMB=EvtPDL::getId("anti-nu_mu");
85   static EvtId NUTB=EvtPDL::getId("anti-nu_tau");
86
87   p->initializePhaseSpace(getNDaug(),getDaugs());
88   
89   double m_b = p->mass();
90
91   EvtParticle *meson, *neutrino1, *neutrino2;
92   meson = p->getDaug(0);
93   neutrino1 = p->getDaug(1);
94   neutrino2 = p->getDaug(2);
95   EvtVector4R momnu1 = neutrino1->getP4();
96   EvtVector4R momnu2 = neutrino2->getP4();
97   EvtVector4R momkstar = meson->getP4();
98   
99   double v0_0, a1_0, a2_0;
100   double m2v0, a1_b, a2_b;
101   v0_0 = 0.47;
102   a1_0 = 0.37;
103   a2_0 = 0.40;
104   m2v0 = 5.*5.;
105   a1_b = -0.023;
106   a2_b = 0.034;
107
108   EvtVector4R q = momnu1+momnu2;
109   double q2 = q.mass2();   
110
111   double v0, a1, a2;
112   v0 = v0_0/(1-q2/m2v0);
113   a1 = a1_0*(1+a1_b*q2);
114   a2 = a2_0*(1+a2_b*q2);
115
116   EvtVector4R p4b; p4b.set(m_b,0.,0.,0.);  // Do calcs in mother rest frame
117   
118   double m_k = meson->mass();
119
120   EvtTensor4C tds=(-2*v0/(m_b+m_k))*dual(EvtGenFunctions::directProd(p4b,momkstar))
121     - EvtComplex(0.0,1.0)*
122     ( (m_b+m_k)*a1*EvtTensor4C::g()
123       - (a2/(m_b+m_k))*EvtGenFunctions::directProd(p4b-momkstar,p4b+momkstar)); 
124   
125   EvtVector4C l;
126
127   if (getDaug(1)==NUE||getDaug(1)==NUM||getDaug(1)==NUT) {
128         l=EvtLeptonVACurrent(neutrino1->spParentNeutrino(),
129                              neutrino2->spParentNeutrino());
130   }
131   if (getDaug(1)==NUEB||getDaug(1)==NUMB||getDaug(1)==NUTB) {
132         l=EvtLeptonVACurrent(neutrino2->spParentNeutrino(),
133                              neutrino1->spParentNeutrino());
134   }
135
136   EvtVector4C et0,et1,et2; 
137   et0 = tds.cont1( meson->epsParent(0).conj() );
138   et1 = tds.cont1( meson->epsParent(1).conj() );
139   et2 = tds.cont1( meson->epsParent(2).conj() );
140
141   vertex(0,l*et0);
142   vertex(1,l*et1);
143   vertex(2,l*et2);
144
145   return;
146 }
147