]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TEvtGen/EvtGenModels/EvtSVSCPiso.cpp
fine tuning of TOF tail (developing task)
[u/mrichter/AliRoot.git] / TEvtGen / EvtGenModels / EvtSVSCPiso.cpp
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: EvtSVSCPiso.cc
12//
13// Description: Routine to decay scalar -> vectors scalar
14// with CP violation and isospin amplitudes.
15// More specifically, it is indended to handle
16// decays like B->rho pi and B->a1 pi.
17//
18// Modification history:
19//
20// RYD/NK Febuary 16, 1998 Module created
21//
22//------------------------------------------------------------------------
23//
24#include "EvtGenBase/EvtPatches.hh"
25#include <stdlib.h>
26#include "EvtGenBase/EvtParticle.hh"
27#include "EvtGenBase/EvtRandom.hh"
28#include "EvtGenBase/EvtGenKine.hh"
29#include "EvtGenBase/EvtCPUtil.hh"
30#include "EvtGenBase/EvtPDL.hh"
31#include "EvtGenBase/EvtReport.hh"
32#include "EvtGenBase/EvtVector4C.hh"
33#include "EvtGenModels/EvtSVSCPiso.hh"
34#include "EvtGenBase/EvtId.hh"
35#include <string>
36#include "EvtGenBase/EvtConst.hh"
37
38EvtSVSCPiso::~EvtSVSCPiso() {}
39
40std::string EvtSVSCPiso::getName(){
41
42 return "SVS_CP_ISO";
43
44}
45
46
47EvtDecayBase* EvtSVSCPiso::clone(){
48
49 return new EvtSVSCPiso;
50
51}
52
53void EvtSVSCPiso::init(){
54
0ca57c2f 55 // check that there are 27 arguments
56 checkNArg(27);
da0e9ce3 57 checkNDaug(2);
58
59 checkSpinParent(EvtSpinType::SCALAR);
60
61 checkSpinDaughter(0,EvtSpinType::VECTOR);
62 checkSpinDaughter(1,EvtSpinType::SCALAR);
63
64}
65
66
67void EvtSVSCPiso::initProbMax(){
68
69//this might need some revision..
70
71if ((EvtPDL::chg3(getDaug(0)) > 0) && (EvtPDL::chg3(getDaug(1)) == 0)) {
72 setProbMax(2.0*(getArg(3)*getArg(3) + 4.0*getArg(23)*getArg(23)));
73}
74
75if ((EvtPDL::chg3(getDaug(0)) < 0) && (EvtPDL::chg3(getDaug(1)) == 0)) {
76 setProbMax(2.0*(getArg(5)*getArg(5) + 4.0*getArg(25)*getArg(25)));
77}
78
79if ((EvtPDL::chg3(getDaug(0)) == 0) && (EvtPDL::chg3(getDaug(1)) > 0)) {
80 setProbMax(2.0*(getArg(7)*getArg(7) + 4.0*getArg(23)*getArg(23)));
81}
82
83if ((EvtPDL::chg3(getDaug(0)) == 0) && (EvtPDL::chg3(getDaug(1)) < 0)) {
84 setProbMax(2.0*(getArg(9)*getArg(9) + 4.0*getArg(25)*getArg(25)));
85}
86
87if ((EvtPDL::chg3(getDaug(0)) > 0) && (EvtPDL::chg3(getDaug(1)) < 0)) {
88 setProbMax(2.0*(getArg(11)*getArg(11) + getArg(23)*getArg(23) +
89 getArg(19)*getArg(19) + getArg(13)*getArg(13) +
90 getArg(25)*getArg(25) + getArg(21)*getArg(21)));
91}
92
93if ((EvtPDL::chg3(getDaug(0)) < 0) && (EvtPDL::chg3(getDaug(1)) > 0)) {
94 setProbMax(2.0*(getArg(15)*getArg(15) + getArg(23)*getArg(23) +
95 getArg(19)*getArg(19) + getArg(17)*getArg(17) +
96 getArg(25)*getArg(25) + getArg(21)*getArg(21)));
97}
98
99if ((EvtPDL::chg3(getDaug(0)) == 0) && (EvtPDL::chg3(getDaug(1)) == 0)) {
100 setProbMax(2.0*(getArg(7)*getArg(7) + getArg(3)*getArg(3) + getArg(11)*getArg(11) +
101 getArg(15)*getArg(15) + 4.0*getArg(19)*getArg(19) + getArg(9)*getArg(9)+
102 getArg(5)*getArg(5) + getArg(13)*getArg(13) + getArg(17)*getArg(17) +
103 4.0*getArg(21)*getArg(21)));
104}
105
106}
107
108
109void EvtSVSCPiso::decay( EvtParticle *p){
110
111 //added by Lange Jan4,2000
112 static EvtId B0=EvtPDL::getId("B0");
113 static EvtId B0B=EvtPDL::getId("anti-B0");
114
115 double t;
116 EvtId other_b;
117 int charged(0);
118
119 int first_time=0;
120 int flip=0;
121 EvtId ds[2];
122
123
124//randomly generate the tag (B0 or B0B)
125
126 double tag = EvtRandom::Flat(0.0,1.0);
127 if (tag < 0.5) {
128
0ca57c2f 129 EvtCPUtil::getInstance()->OtherB(p,t,other_b,1.0);
130 other_b = B0;
da0e9ce3 131 }
132 else {
133
0ca57c2f 134 EvtCPUtil::getInstance()->OtherB(p,t,other_b,0.0);
135 other_b = B0B;
da0e9ce3 136 }
137
138 if (p->getNDaug()==0) first_time=1;
139
140 if (first_time){
0ca57c2f 141 if (EvtRandom::Flat(0.0,1.0)<getArg(2)) flip=1;
da0e9ce3 142 }
143 else{
144 if (getDaug(0)!=p->getDaug(0)->getId()) flip=1;
145 }
146
147 if (!flip) {
148 ds[0]=getDaug(0);
149 ds[1]=getDaug(1);
150 }
151 else{
152 ds[0]=EvtPDL::chargeConj(getDaug(0));
153 ds[1]=EvtPDL::chargeConj(getDaug(1));
154 }
155
156 p->initializePhaseSpace(getNDaug(),ds);
157
158 EvtParticle *v,*s;
159 v=p->getDaug(0);
160 s=p->getDaug(1);
161
162 EvtComplex amp;
163
164 EvtComplex A_f,Abar_f;
165 EvtComplex A_fbar,Abar_fbar;
166 EvtComplex Apm, Apm_bar, Amp, Amp_bar;
167
168 EvtComplex Tp0, Tp0_bar, T0p, T0p_bar,Tpm, Tpm_bar, Tmp, Tmp_bar;
169 EvtComplex P1, P1_bar, P0, P0_bar;
170
171 Tp0 = EvtComplex(getArg(3)*cos(getArg(4)),getArg(3)*sin(getArg(4)));
172 Tp0_bar = EvtComplex(getArg(5)*cos(getArg(6)),getArg(5)*sin(getArg(6)));
173 T0p = EvtComplex(getArg(7)*cos(getArg(8)),getArg(7)*sin(getArg(8)));
174 T0p_bar = EvtComplex(getArg(9)*cos(getArg(10)),getArg(9)*sin(getArg(10)));
175 Tpm = EvtComplex(getArg(11)*cos(getArg(12)),getArg(11)*sin(getArg(12)));
176 Tpm_bar = EvtComplex(getArg(13)*cos(getArg(14)),getArg(13)*sin(getArg(14)));
177 Tmp = EvtComplex(getArg(15)*cos(getArg(16)),getArg(15)*sin(getArg(16)));
178 Tmp_bar = EvtComplex(getArg(17)*cos(getArg(18)),getArg(17)*sin(getArg(18)));
179 P0 = EvtComplex(getArg(19)*cos(getArg(20)),getArg(19)*sin(getArg(20)));
180 P0_bar = EvtComplex(getArg(21)*cos(getArg(22)),getArg(21)*sin(getArg(22)));
181 P1 = EvtComplex(getArg(23)*cos(getArg(24)),getArg(23)*sin(getArg(24)));
182 P1_bar = EvtComplex(getArg(25)*cos(getArg(26)),getArg(25)*sin(getArg(26)));
183
184
185//***********************charged modes****************************
186
187 if ((EvtPDL::chg3(getDaug(0)) > 0 ) && (EvtPDL::chg3(getDaug(1)) == 0)) {
188
189//V+ S0, so T+0 + 2 P1
190
191 charged = 1;
192 A_f = Tp0 + 2.0*P1;
193 }
194
195 if ((EvtPDL::chg3(getDaug(0)) < 0 ) && (EvtPDL::chg3(getDaug(1)) == 0)) {
196
197//V- S0, so T+0_bar + 2P1_bar
198
199 charged = 1;
200 A_f = Tp0_bar + 2.0*P1_bar;
201 }
202
203 if ((EvtPDL::chg3(getDaug(0)) == 0 ) && (EvtPDL::chg3(getDaug(1)) > 0)) {
204
205//V0 S+, so T0+ - 2 P1
206
207 charged = 1;
208 A_f = T0p - 2.0*P1;
209 }
210
211 if ((EvtPDL::chg3(getDaug(0)) == 0 ) && (EvtPDL::chg3(getDaug(1)) < 0)) {
212
213//V0 S-, so T0+_bar - 2 P1_bar
214
215 charged = 1;
216 A_f = T0p_bar - 2.0*P1_bar;
217 }
218
219
220//***********************neutral modes***************************
221
222
223//V+ S-, so Af = T+- + P1 + P0
224Apm = Tpm + P1 + P0;
225Apm_bar = Tpm_bar + P1_bar + P0_bar;
226
227//V- S+, so Af = T-+ - P1 + P0
228Amp = Tmp - P1 + P0;
229Amp_bar = Tmp_bar - P1_bar + P0;
230
231
232 if ((EvtPDL::chg3(getDaug(0)) > 0 ) && (EvtPDL::chg3(getDaug(1)) < 0)) {
233
234//V+ S-
235 charged = 0;
236 A_f = Apm;
237 Abar_f = Apm_bar;
238 A_fbar = Amp;
239 Abar_fbar = Amp_bar;
240
241 }
242
243 if ((EvtPDL::chg3(getDaug(0)) < 0 ) && (EvtPDL::chg3(getDaug(1)) > 0)) {
244
245//V- S+
246 charged = 0;
247 A_f = Amp;
248 Abar_f = Amp_bar;
249 A_fbar = Apm;
250 Abar_fbar = Apm_bar;
251
252 }
253
254 if ((EvtPDL::chg3(getDaug(0)) == 0 ) && (EvtPDL::chg3(getDaug(1)) == 0)) {
255
256//V0 S0
257 charged = 0;
258 A_f = T0p + Tp0 - Tpm - Tmp - 2.0*P0 ;
259 Abar_f = T0p_bar + Tp0_bar - Tpm_bar - Tmp_bar - 2.0*P0_bar;
260 A_fbar = A_f;
261 Abar_fbar = Abar_f;
262
263 }
264
265if (charged==0) {
266
267 if (!flip) {
268 if (other_b==B0B){
269
270 amp=A_f*cos(getArg(1)*t/(2*EvtConst::c))+
271 EvtComplex(cos(-2.0*getArg(0)),sin(-2.0*getArg(0)))*
272 EvtComplex(0.0,1.0)*Abar_f*sin(getArg(1)*t/(2*EvtConst::c));
273 }
274 if (other_b==B0){
275
276 amp=A_f*EvtComplex(cos(2.0*getArg(0)),sin(2.0*getArg(0)))*
277 EvtComplex(0.0,1.0)*sin(getArg(1)*t/(2*EvtConst::c))+
278 Abar_f*cos(getArg(1)*t/(2*EvtConst::c));
279 }
280 }
281 else{
282 if (other_b==B0B){
283
284 amp=A_fbar*cos(getArg(1)*t/(2*EvtConst::c))+
285 EvtComplex(cos(-2.0*getArg(0)),sin(-2.0*getArg(0)))*
286 EvtComplex(0.0,1.0)*Abar_fbar*sin(getArg(1)*t/(2*EvtConst::c));
287 }
288 if (other_b==B0){
289
290 amp=A_fbar*EvtComplex(cos(2.0*getArg(0)),sin(2.0*getArg(0)))*
291 EvtComplex(0.0,1.0)*sin(getArg(1)*t/(2*EvtConst::c))+
292 Abar_fbar*cos(getArg(1)*t/(2*EvtConst::c));
293 }
294 }
295
296}
297else amp = A_f;
298
299 EvtVector4R p4_parent;
300
301 p4_parent=v->getP4()+s->getP4();
302
303 double norm=1.0/v->getP4().d3mag();
304
305 vertex(0,amp*norm*p4_parent*(v->epsParent(0)));
306 vertex(1,amp*norm*p4_parent*(v->epsParent(1)));
307 vertex(2,amp*norm*p4_parent*(v->epsParent(2)));
308
309 return ;
310}
311