Fix for definitions for CINT
[u/mrichter/AliRoot.git] / TEvtGen / EvtGenModels / EvtPi0Dalitz.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: EvtPi0Dalitz.cc
12 //
13 // Description: pi0 -> e+ e- gamma 
14 //
15 // Modification history:
16 //
17 //    DJL/RYD    June 30, 1998          Module created
18 //
19 //------------------------------------------------------------------------
20 //
21 #include "EvtGenBase/EvtPatches.hh"
22 #include <stdlib.h>
23 #include <fstream>
24 #include <stdio.h>
25 #include <string>
26 #include "EvtGenBase/EvtGenKine.hh"
27 #include "EvtGenBase/EvtParticle.hh"
28 #include "EvtGenBase/EvtPDL.hh"
29 #include "EvtGenBase/EvtReport.hh"
30 #include "EvtGenModels/EvtPi0Dalitz.hh"
31 #include "EvtGenBase/EvtVector4C.hh"
32 #include "EvtGenBase/EvtDiracSpinor.hh"
33 #include "EvtGenBase/EvtVector4C.hh"
34 #include "EvtGenBase/EvtTensor4C.hh"
35 using std::fstream;
36
37 EvtPi0Dalitz::~EvtPi0Dalitz() {}
38
39 std::string EvtPi0Dalitz::getName(){
40
41     return "PI0_DALITZ";
42
43 }
44
45 EvtDecayBase* EvtPi0Dalitz::clone(){
46
47   return new EvtPi0Dalitz;
48
49 }
50
51
52 void EvtPi0Dalitz::initProbMax(){
53
54   setProbMax(3.5);
55
56 }
57
58
59 void EvtPi0Dalitz::init(){
60
61   // check that there are 0 arguments
62   checkNArg(0);
63   checkNDaug(3);
64
65     
66   checkSpinParent(EvtSpinType::SCALAR);
67
68   checkSpinDaughter(0,EvtSpinType::DIRAC);
69   checkSpinDaughter(1,EvtSpinType::DIRAC);
70   checkSpinDaughter(2,EvtSpinType::PHOTON);
71
72 }
73
74
75 void EvtPi0Dalitz::decay( EvtParticle *p){
76
77   EvtParticle *ep, *em, *gamma;
78   setWeight(p->initializePhaseSpace(getNDaug(),getDaugs(),false,
79                                     0.00000002,0,1));
80
81   ep=p->getDaug(0);
82   em=p->getDaug(1);
83   gamma=p->getDaug(2);
84
85   // the next four lines generates events with a weight such that
86   // the efficiency for selecting them is good. The parameter below of
87   // 0.1 is the size of the peak at low q^2 (in arbitrary units).
88   // The value of 0.1 is appropriate for muons. 
89   // when you use this remember to remove the cut on q^2!
90    
91
92   //ep em invariant mass^2
93   double m2=(ep->getP4()+em->getP4()).mass2();
94   EvtVector4R q=ep->getP4()+em->getP4();
95   //Just use the prob summed over spins...
96
97   EvtTensor4C w,v;
98
99   v=2.0*(gamma->getP4()*q)*EvtGenFunctions::directProd(q,gamma->getP4()) 
100     - (gamma->getP4()*q)*(gamma->getP4()*q)*EvtTensor4C::g()
101     -m2*EvtGenFunctions::directProd(gamma->getP4(),gamma->getP4());
102  
103   w=4.0*( EvtGenFunctions::directProd(ep->getP4(),em->getP4()) + 
104           EvtGenFunctions::directProd(em->getP4(),ep->getP4())
105            -EvtTensor4C::g()*(ep->getP4()*em->getP4()-ep->getP4().mass2()));
106
107   double prob=(real(cont(v,w)))/(m2*m2);
108   prob *=(1.0/( (0.768*0.768-m2)*(0.768*0.768-m2)
109            +0.768*0.768*0.151*0.151));
110
111   //  report(INFO,"EvtGen") << "prob is "<<prob<<endl;
112   setProb(prob);
113
114   return;
115 }
116
117