]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TEvtGen/EvtGenModels/EvtBHadronic.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Removing the flat makefiles
[u/mrichter/AliRoot.git] / TEvtGen / EvtGenModels / EvtBHadronic.cxx
CommitLineData
da0e9ce3 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: EvtBHadronic.cc
12//
13// Description: Model for hadronic B decays. Uses naive factorization.
14//
15// Modification history:
16//
17// RYD June 14, 1997 Module created
18//
19//------------------------------------------------------------------------
20//
21#include "EvtGenBase/EvtPatches.hh"
22#include <stdlib.h>
23#include "EvtGenBase/EvtParticle.hh"
24#include "EvtGenModels/EvtISGW2FF.hh"
25#include "EvtGenBase/EvtGenKine.hh"
26#include "EvtGenBase/EvtPDL.hh"
27#include "EvtGenModels/EvtISGW2FF.hh"
28#include "EvtGenBase/EvtTensor4C.hh"
29#include "EvtGenBase/EvtVector4C.hh"
30#include "EvtGenBase/EvtVector4R.hh"
31#include "EvtGenModels/EvtBHadronic.hh"
32#include "EvtGenBase/EvtReport.hh"
33#include <string>
34using std::endl;
35
36EvtBHadronic::~EvtBHadronic() {}
37
38std::string EvtBHadronic::getName(){
39
40 return "BHADRONIC";
41
42}
43
44
45EvtDecayBase* EvtBHadronic::clone(){
46
47 return new EvtBHadronic;
48
49}
50
51
52void EvtBHadronic::init(){
53
54 // check that there are 2 argument
55 checkNArg(2);
56
57
58}
59
60void EvtBHadronic::decay( EvtParticle *p){
61
62 //added by Lange Jan4,2000
63 static EvtId B0=EvtPDL::getId("B0");
64 static EvtId D0=EvtPDL::getId("D0");
65 static EvtId DST0=EvtPDL::getId("D*0");
66 static EvtId D3P10=EvtPDL::getId("D'_10");
67 static EvtId D3P20=EvtPDL::getId("D_2*0");
68 static EvtId D3P00=EvtPDL::getId("D_0*0");
69 static EvtId D1P10=EvtPDL::getId("D_10");
70
71 static EvtISGW2FF ffmodel;
72
73 p->initializePhaseSpace(getNDaug(),getDaugs());
74
75 EvtVector4R p4[MAX_DAUG];
76 double m;
77
78 m = p->mass();
79
80
81 int i,j;
82
83 for ( i=0; i<getNDaug(); i++) {
84 p4[i]=p->getDaug(i)->getP4();
85 }
86
87 int bcurrent,wcurrent;
88 int nbcurrent=0;
89 int nwcurrent=0;
90
91 bcurrent=(int)getArg(0);
92 wcurrent=(int)getArg(1);
93
94 EvtVector4C jb[5],jw[5];
95 EvtTensor4C g,tds;
96
97 EvtVector4R p4b;
98 p4b.set(p->mass(),0.0,0.0,0.0);
99
100 EvtVector4R q;
101 double q2;
102
103 EvtComplex ep_meson_bb[5];
104 double f,gf,ap,am;
105 double fp,fm;
106 double hf,kf,bp,bm;
107 EvtVector4R pp,pm;
108 EvtVector4C ep_meson_b[5];
109
110 switch (bcurrent) {
111
112 // D0
113 case 1:
114 q=p4b-p4[0];
115 q2=q*q;
116 nbcurrent=1;
117 ffmodel.getscalarff(B0,D0,EvtPDL::getMeanMass(D0),
118 EvtPDL::getMeanMass(getDaugs()[1]),&fp,&fm);
119 jb[0]=EvtVector4C(fp*(p4b+p4[0])-fm*q);
120 break;
121 // D*
122 case 2:
123 q=p4b-p4[0];
124 q2=q*q;
125 nbcurrent=3;
126 ffmodel.getvectorff(B0,DST0,EvtPDL::getMeanMass(DST0),q2,&f,&gf,&ap,&am);
127
128 g.setdiag(1.0,-1.0,-1.0,-1.0);
129 tds = -f*g
130 -ap*(directProd(p4b,p4b)+directProd(p4b,p4[0]))
131 -gf*EvtComplex(0.0,1.0)*dual(directProd(p4[0]+p4b,p4b-p4[0]))
132 -am*((directProd(p4b,p4b)-directProd(p4b,p4[0])));
133 jb[0]=tds.cont1(p->getDaug(0)->epsParent(0).conj());
134 jb[1]=tds.cont1(p->getDaug(0)->epsParent(1).conj());
135 jb[2]=tds.cont1(p->getDaug(0)->epsParent(2).conj());
136 break;
137 // D1
138 case 3:
139 q=p4b-p4[0];
140 q2=q*q;
141 nbcurrent=3;
142 ffmodel.getvectorff(B0,D3P10,EvtPDL::getMeanMass(D3P10),q2,&f,&gf,&ap,&am);
143
144 g.setdiag(1.0,-1.0,-1.0,-1.0);
145 tds = -f*g
146 -ap*(directProd(p4b,p4b)+directProd(p4b,p4[0]))
147 -gf*EvtComplex(0.0,1.0)*dual(directProd(p4[0]+p4b,p4b-p4[0]))
148 -am*((directProd(p4b,p4b)-directProd(p4b,p4[0])));
149 jb[0]=tds.cont1(p->getDaug(0)->epsParent(0).conj());
150 jb[1]=tds.cont1(p->getDaug(0)->epsParent(1).conj());
151 jb[2]=tds.cont1(p->getDaug(0)->epsParent(2).conj());
152
153 break;
154 // D2*
155 case 4:
156 q=p4b-p4[0];
157 q2=q*q;
158 nbcurrent=5;
159 ffmodel.gettensorff(B0,D3P20,EvtPDL::getMeanMass(D3P20),q2,&hf,&kf,&bp,&bm);
160 g.setdiag(1.0,-1.0,-1.0,-1.0);
161
162 ep_meson_b[0] = ((p->getDaug(0)->epsTensorParent(0)).cont2(p4b)).conj();
163 ep_meson_b[1] = ((p->getDaug(0)->epsTensorParent(1)).cont2(p4b)).conj();
164 ep_meson_b[2] = ((p->getDaug(0)->epsTensorParent(2)).cont2(p4b)).conj();
165 ep_meson_b[3] = ((p->getDaug(0)->epsTensorParent(3)).cont2(p4b)).conj();
166 ep_meson_b[4] = ((p->getDaug(0)->epsTensorParent(4)).cont2(p4b)).conj();
167
168 pp=p4b+p4[0];
169 pm=p4b-p4[0];
170
171 ep_meson_bb[0]=ep_meson_b[0]*(p4b);
172 ep_meson_bb[1]=ep_meson_b[1]*(p4b);
173 ep_meson_bb[2]=ep_meson_b[2]*(p4b);
174 ep_meson_bb[3]=ep_meson_b[3]*(p4b);
175 ep_meson_bb[4]=ep_meson_b[4]*(p4b);
176
177 jb[0]=EvtComplex(0.0,(1.0)*hf)*dual(directProd(pp,pm)).cont2(ep_meson_b[0])
178 -kf*ep_meson_b[0]
179 -bp*ep_meson_bb[0]*pp-bm*ep_meson_bb[0]*pm;
180
181 jb[1]=EvtComplex(0.0,(1.0)*hf)*dual(directProd(pp,pm)).cont2(ep_meson_b[1])
182 -kf*ep_meson_b[1]
183 -bp*ep_meson_bb[1]*pp-bm*ep_meson_bb[1]*pm;
184
185 jb[2]=EvtComplex(0.0,(1.0)*hf)*dual(directProd(pp,pm)).cont2(ep_meson_b[2])
186 -kf*ep_meson_b[2]
187 -bp*ep_meson_bb[2]*pp-bm*ep_meson_bb[2]*pm;
188
189 jb[3]=EvtComplex(0.0,(1.0)*hf)*dual(directProd(pp,pm)).cont2(ep_meson_b[3])
190 -kf*ep_meson_b[3]
191 -bp*ep_meson_bb[3]*pp-bm*ep_meson_bb[3]*pm;
192
193 jb[4]=EvtComplex(0.0,(1.0)*hf)*dual(directProd(pp,pm)).cont2(ep_meson_b[4])
194 -kf*ep_meson_b[4]
195 -bp*ep_meson_bb[4]*pp-bm*ep_meson_bb[4]*pm;
196 break;
197 // D_0*
198 case 5:
199 q=p4b-p4[0];
200 q2=q*q;
201 double f,gf,ap,am;
202 nbcurrent=3;
203 ffmodel.getvectorff(B0,D1P10,EvtPDL::getMeanMass(D1P10),q2,&f,&gf,&ap,&am);
204 g.setdiag(1.0,-1.0,-1.0,-1.0);
205 tds = -f*g
206 -ap*(directProd(p4b,p4b)+directProd(p4b,p4[0]))
207 +gf*EvtComplex(0.0,1.0)*dual(directProd(p4[0]+p4b,p4b-p4[0]))
208 -am*((directProd(p4b,p4b)-directProd(p4b,p4[0])));
209 jb[0]=tds.cont1(p->getDaug(0)->epsParent(0).conj());
210 jb[1]=tds.cont1(p->getDaug(0)->epsParent(1).conj());
211 jb[2]=tds.cont1(p->getDaug(0)->epsParent(2).conj());
212 break;
213 // D_1
214 case 6:
215 q=p4b-p4[0];
216 q2=q*q;
217 nbcurrent=1;
218 ffmodel.getscalarff(B0,D3P00,EvtPDL::getMeanMass(D3P00),q2,&fp,&fm);
219 jb[0]=fp*(p4b+p4[0])+fm*q;
220 break;
221 default:
222 report(ERROR,"EvtGen")<<"In EvtBHadronic, unknown hadronic current."<<endl;
223
224 }
225
226 double norm;
227
228 switch (wcurrent) {
229
230
231 case 1: case 3: case 4:
232 nwcurrent=1;
233 jw[0]=p4[getNDaug()-1];
234 break;
235
236 case 2: case 5: case 6:
237 nwcurrent=3;
238 norm=1.0/sqrt(p4[1].get(0)*p4[1].get(0)/p4[1].mass2()-1.0);
239 jw[0]=norm*p->getDaug(getNDaug()-1)->epsParent(0);
240 jw[1]=norm*p->getDaug(getNDaug()-1)->epsParent(1);
241 jw[2]=norm*p->getDaug(getNDaug()-1)->epsParent(2);
242 break;
243
244
245 default:
246 report(ERROR,"EvtGen")<<"In EvtBHadronic, unknown W current."<<endl;
247
248 }
249
250 if (nbcurrent==1&&nwcurrent==1){
251 vertex(jb[0]*jw[0]);
252 return;
253 }
254
255 if (nbcurrent==1){
256 for(j=0;j<nwcurrent;j++){
257 vertex(j,jb[0]*jw[j]);
258 }
259 return;
260 }
261
262 if (nwcurrent==1){
263 for(j=0;j<nbcurrent;j++){
264 vertex(j,jb[j]*jw[0]);
265 }
266 return;
267 }
268
269 for(j=0;j<nbcurrent;j++){
270 for(i=0;i<nwcurrent;i++){
271 vertex(j,i,jb[j]*jw[i]);
272 }
273 }
274 return;
275
276}
277
278
279
280
281
282
283