Converting TEvtGen to native cmake
[u/mrichter/AliRoot.git] / TEvtGen / EvtGen / EvtGenModels / EvtBTo3piCP.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: EvtBTo3piCP.cc
12 //
13 // Description: Routine to decay B->pi+ pi- pi0
14 //              and has CP violation.
15 //
16 // Modification history:
17 //
18 //    RYD/VERSILLE     March 2, 1997        Module created
19 //
20 //------------------------------------------------------------------------
21 //
22 #include "EvtGenBase/EvtPatches.hh"
23 #include <stdlib.h>
24 #include "EvtGenBase/EvtParticle.hh"
25 #include "EvtGenBase/EvtCPUtil.hh"
26 #include "EvtGenBase/EvtPDL.hh"
27 #include "EvtGenBase/EvtReport.hh"
28 #include "EvtGenModels/EvtBTo3piCP.hh"
29 #include "EvtGenBase/EvtId.hh"
30 #include <string>
31 #include "EvtGenBase/EvtConst.hh"
32
33 #ifdef WIN32
34 extern "C" {
35   extern void __stdcall EVT3PIONS(double *,int *,double *,
36                                   double *,double *,double *,double *,
37                                   double *,double *,double *);
38 }
39 #else
40 extern "C" {
41   extern void evt3pions_(double *,int *,double *,
42                          double *,double *,double *,double *,
43                          double *,double *,double *);
44 }
45 #endif
46
47 EvtBTo3piCP::~EvtBTo3piCP() {}
48
49
50 std::string EvtBTo3piCP::getName(){
51
52   return "BTO3PI_CP";     
53
54 }
55
56
57 EvtDecayBase* EvtBTo3piCP::clone(){
58
59   return new EvtBTo3piCP;
60
61 }
62
63 void EvtBTo3piCP::init(){
64
65   // check that there are 2 arguments
66   checkNArg(2);
67   checkNDaug(3);
68
69   checkSpinParent(EvtSpinType::SCALAR);
70
71   checkSpinDaughter(0,EvtSpinType::SCALAR);
72   checkSpinDaughter(1,EvtSpinType::SCALAR);
73   checkSpinDaughter(2,EvtSpinType::SCALAR);
74 }
75
76
77
78 void EvtBTo3piCP::initProbMax(){
79
80   // perform common blocks initialization before
81   // first use
82   double alpha=getArg(1);
83   int iset;
84
85   iset=10000;
86
87   double p4piplus[4],p4piminus[4],p4gamm1[4],p4gamm2[4]; 
88
89   double realA,imgA,realbarA,imgbarA;
90
91 #ifdef WIN32
92   EVT3PIONS(&alpha,&iset,p4piplus,p4piminus,p4gamm1,p4gamm2,
93              &realA,&imgA,&realbarA,&imgbarA);
94 #else
95   evt3pions_(&alpha,&iset,p4piplus,p4piminus,p4gamm1,p4gamm2,
96              &realA,&imgA,&realbarA,&imgbarA);
97 #endif
98
99   setProbMax(1.5);
100
101 }
102
103 void EvtBTo3piCP::decay( EvtParticle *p){
104
105   //added by Lange Jan4,2000
106   static EvtId B0=EvtPDL::getId("B0");
107   static EvtId B0B=EvtPDL::getId("anti-B0");
108
109   double t;
110   EvtId other_b;
111
112   EvtCPUtil::getInstance()->OtherB(p,t,other_b,0.5);
113
114   EvtParticle *pip,*pim,*pi0;
115
116   p->makeDaughters(getNDaug(),getDaugs());
117
118   //  p->init_daug(SCALAR,&pip,SCALAR,&pim,SCALAR,&pi0);
119   pip=p->getDaug(0);
120   pim=p->getDaug(1);
121   pi0=p->getDaug(2);
122
123   EvtVector4R p4[3];
124
125   double dm=getArg(0);
126   double alpha=getArg(1);
127   int iset;
128
129   iset=0;
130
131   double p4piplus[4],p4piminus[4],p4gamm1[4],p4gamm2[4]; 
132
133   double realA,imgA,realbarA,imgbarA;
134
135 #ifdef WIN32
136   EVT3PIONS(&alpha,&iset,p4piplus,p4piminus,p4gamm1,p4gamm2,
137             &realA,&imgA,&realbarA,&imgbarA);
138 #else
139   evt3pions_(&alpha,&iset,p4piplus,p4piminus,p4gamm1,p4gamm2,
140              &realA,&imgA,&realbarA,&imgbarA);
141 #endif
142
143   p4[0].set(p4piplus[3],p4piplus[0],p4piplus[1],p4piplus[2]);
144   p4[1].set(p4piminus[3],p4piminus[0],p4piminus[1],p4piminus[2]);
145   p4[2].set(p4gamm1[3]+p4gamm2[3],p4gamm1[0]+p4gamm2[0],
146             p4gamm1[1]+p4gamm2[1],p4gamm1[2]+p4gamm2[2]);
147
148   if (pip->getId()==EvtPDL::getId("pi+")) {
149     pip->init( getDaug(0), p4[0] );
150     pim->init( getDaug(1), p4[1] );
151   }
152   else {
153     pip->init( getDaug(0), p4[1] );
154     pim->init( getDaug(1), p4[0] );  
155   }
156
157    pi0->init( getDaug(2), p4[2] );
158    
159    EvtComplex amp;
160
161    EvtComplex A(realA,imgA);
162    EvtComplex Abar(realbarA,imgbarA);
163
164    if (other_b==B0B){
165      amp=A*cos(dm*t/(2*EvtConst::c))+
166        EvtComplex(0.,1.)*Abar*sin(dm*t/(2*EvtConst::c));
167    }
168    if (other_b==B0){
169      amp=Abar*cos(dm*t/(2*EvtConst::c))+
170        EvtComplex(0.,1.)*A*sin(dm*t/(2*EvtConst::c));
171    }
172
173    vertex(amp);
174
175   return ;
176 }
177