BOSS 7.1.1
BESIII Offline Software System
Loading...
Searching...
No Matches
EvtD0ToKSLKK.cc
Go to the documentation of this file.
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: EvtD0ToKSLKK.cc
12//
13// Description: Routine to handle D0 -> KS0 K+ K- and KL0 K+ K-
14//
15// Modification history:
16//
17// Liaoyuan Dong Oct 24 15:17:18 2023 Module created
18//
19//------------------------------------------------------------------------
20//
22#include <fstream>
23#include <stdlib.h>
26#include "EvtGenBase/EvtPDL.hh"
31#include <string>
34using std::endl;
35
37
38void EvtD0ToKSLKK::getName(std::string& model_name){
39 model_name="D0ToKSLKK";
40}
41
45
47// checkNArg(1);
48 Narg = getNArg();
49 if(Narg == 0){
50 Uspin = 0;
51 }else{
52 Uspin = getArg(0);
53 }
55 Daug0Id = EvtPDL::getStdHep(getDaug(0));
56 checkNDaug(3);
61
62 phi[0] = 0;
63 rho[0] = 1;
64 phi[1] = 3.367900972;
65 rho[1] = 0.5567227875;
66 phi[2] = 3.568151995;
67 rho[2] = 1.122904371;
68 phi[3] = 3.353326851;
69 rho[3] = 1.414616255;
70 phi[4] = -0.7616302559;
71 rho[4] = 1.623197367;
72
73 modetype[0]= 23;
74 modetype[1]= 23;
75 modetype[2]= 13;
76 modetype[3]= 13;
77 modetype[4]= 13;
78
79 if(Uspin==1&&Daug0Id==130){
80 modetype[0]= 232;
81 modetype[1]= 232;
82 modetype[2]= 13;
83 modetype[3]= 13;
84 modetype[4]= 13;
85 }
86
87/*
88 for (int i=0; i<5; i++) {
89 cout << i << " rho= " << rho[i] << " phi= " << phi[i] << endl;
90 }
91*/
92 width[0] = 0.272;
93 width[1] = 0.004249;
94 width[2] = 0.272;
95 width[3] = 0.258;
96 width[4] = 0.4;
97
98 mass[0] = 0.919;
99 mass[1] = 1.019461;
100 mass[2] = 0.919;
101 mass[3] = 1.439;
102 mass[4] = 1.439;
103
104 mDM = 1.86484;
105 mK0 = 0.497614;
106 mKa = 0.49368;
107 mPi = 0.13957;
108 mEta = 0.547862;
109 mKa2 = 0.24371994;//0.49368^2;
110 mPi2 = 0.01947978;//0.13957^2;
111 mEta2 = 0.30015277;//0.547862^2;
112 mass_EtaP = 0.95778;
113 mass_Kaon = 0.49368;
114 mass_KS = 0.4976;
115
116 math_pi = 3.1415926;
117 mass_Pion2 = 0.0194797849;
118 mass_2Pion = 0.27914;
119 math_2pi = 6.2831852;
120 rD2 = 25.0; // 5*5
121 rRes2 = 9.0; // 3*3
122 g2 = 0.23; //K*0(1430)
123
124
125 GS1 = 0.636619783;
126 GS2 = 0.01860182466;
127 GS3 = 0.1591549458; // 1/(2*math_2pi)
128 GS4 = 0.00620060822; // mass_Pion2/math_pi
129
130 rho_omega = 0.00294;
131 phi_omega = -0.02;
132
133 int GG[4][4] = { {1,0,0,0}, {0,-1,0,0}, {0,0,-1,0}, {0,0,0,-1} };
134 for (int i=0; i<4; i++) {
135 for (int j=0; j<4; j++) {
136 G[i][j] = GG[i][j];
137 }
138 }
139}
140
142 setProbMax(1700.0);
143 if(Daug0Id==310) setProbMax(1700.0);
144 else if(Daug0Id==130&&Uspin==1) setProbMax(1785.0);
145 else setProbMax(1700.0);
146}
147
149/*
150 double maxprob = 0.0;
151 for(int ir=0;ir<=60000000;ir++){
152
153 p->initializePhaseSpace(getNDaug(),getDaugs());
154 EvtVector4R D1 = p->getDaug(0)->getP4();
155 EvtVector4R D2 = p->getDaug(1)->getP4();
156 EvtVector4R D3 = p->getDaug(2)->getP4();
157
158 //charge = EvtPDL::getStdHep(p->getId());
159
160 double P1[4], P2[4], P3[4];
161 P1[0] = D1.get(0); P1[1] = D1.get(1); P1[2] = D1.get(2); P1[3] = D1.get(3);
162 P2[0] = D2.get(0); P2[1] = D2.get(1); P2[2] = D2.get(2); P2[3] = D2.get(3);
163 P3[0] = D3.get(0); P3[1] = D3.get(1); P3[2] = D3.get(2); P3[3] = D3.get(3);
164
165 if(Daug0Id==310) SorL = true; else SorL = false;
166 double value;
167 int spin[5]={0,1,0,0,1};
168 if(SorL){
169 int g0[5]={5,1,3,1,5};
170 double r0[5] = {3,3,3,3,3};
171 double r1[5] = {5,5,5,5,5};
172 int nstates=5;
173 calEva(P1, P2, P3, mass, width, rho, phi, g0, spin, modetype, r0, r1, value, 0, nstates,charge,SorL);
174 }else if((!SorL)&&Uspin==1){
175 int g0[5]={5,1,3,1,5};
176 double r0[5] = {-1.566394443,-1.33043736,3,3,3};
177 double r1[5] = {0.1844175671,-1.397710917,5,5,5};
178 int nstates=5;
179 calEva(P1, P2, P3, mass, width, rho, phi, g0, spin, modetype, r0, r1, value, 0, nstates,charge,SorL);
180 }else{
181 int g0[5]={5,1,3,1,5};
182 double r0[5] = {3,3,3,3,3};
183 double r1[5] = {5,5,5,5,5};
184 int nstates=5;
185 calEva(P1, P2, P3, mass, width, rho, phi, g0, spin, modetype, r0, r1, value, 0, nstates,charge,SorL);
186 }
187
188 if (value<0) continue;
189 if(value>maxprob) {
190 maxprob=value;
191 cout << "ir= " << ir << endl;
192 cout << "double P1[4] = {" << P1[0] <<","<< P1[1] <<","<< P1[2] <<","<< P1[3] <<"};"<< endl;
193 cout << "double P2[4] = {" << P2[0] <<","<< P2[1] <<","<< P2[2] <<","<< P2[3] <<"};"<< endl;
194 cout << "double P3[4] = {" << P3[0] <<","<< P3[1] <<","<< P3[2] <<","<< P3[3] <<"};"<< endl;
195 cout << "MAX====> " << maxprob << endl;
196 }
197 }
198 printf("MAXprob = %.10f\n",maxprob);
199*/
200
201
203 EvtVector4R D1 = p->getDaug(0)->getP4();
204 EvtVector4R D2 = p->getDaug(1)->getP4();
205 EvtVector4R D3 = p->getDaug(2)->getP4();
206
207 //charge = EvtPDL::getStdHep(p->getId());
208
209 double P1[4], P2[4], P3[4];
210 P1[0] = D1.get(0); P1[1] = D1.get(1); P1[2] = D1.get(2); P1[3] = D1.get(3);
211 P2[0] = D2.get(0); P2[1] = D2.get(1); P2[2] = D2.get(2); P2[3] = D2.get(3);
212 P3[0] = D3.get(0); P3[1] = D3.get(1); P3[2] = D3.get(2); P3[3] = D3.get(3);
213
214 // P1[0] = 0.662259; P1[1] = -0.277017; P1[2] = 0.211748; P1[3] = -0.263429;
215 // P2[0] = 0.664422; P2[1] = 0.283792; P2[2] = -0.243945; P2[3] = 0.240194;
216 // P3[0] = 0.561591; P3[1] = 0.188743; P3[2] = -0.0572267;P3[3] = -0.181021;
217
218 // P1[0] = 0.663036; P1[1] = -0.188126; P1[2] = 0.265033; P1[3] = -0.293883;
219 // P2[0] = 0.694463; P2[1] = 0.201378; P2[2] = -0.399611; P2[3] = -0.195755;
220 // P3[0] = 0.530119; P3[1] = 0.122859; P3[2] = -0.123733; P3[3] = -0.083098;
221
222
223 if(Daug0Id==310) SorL = true; else SorL = false;
224 double value;
225 int spin[5]={0,1,0,0,1};
226 if(SorL){
227 int g0[5]={5,1,3,1,5};
228 double r0[5] = {3,3,3,3,3};
229 double r1[5] = {5,5,5,5,5};
230 int nstates=5;
231 calEva(P1, P2, P3, mass, width, rho, phi, g0, spin, modetype, r0, r1, value, 0, nstates,charge,SorL);
232 }else if((!SorL)&&Uspin==1){
233 int g0[5]={5,1,3,1,5};
234 double r0[5] = {-1.566394443,-1.33043736,3,3,3};
235 double r1[5] = {0.1844175671,-1.397710917,5,5,5};
236 int nstates=5;
237 calEva(P1, P2, P3, mass, width, rho, phi, g0, spin, modetype, r0, r1, value, 0, nstates,charge,SorL);
238 }else{
239 int g0[5]={5,1,3,1,5};
240 double r0[5] = {3,3,3,3,3};
241 double r1[5] = {5,5,5,5,5};
242 int nstates=5;
243 calEva(P1, P2, P3, mass, width, rho, phi, g0, spin, modetype, r0, r1, value, 0, nstates,charge,SorL);
244 }
245
246 setProb(value);
247
248 return ;
249}
250
251
252void EvtD0ToKSLKK::Com_Multi(double a1[2], double a2[2], double res[2])
253{
254 res[0] = a1[0]*a2[0]-a1[1]*a2[1];
255 res[1] = a1[1]*a2[0]+a1[0]*a2[1];
256}
257void EvtD0ToKSLKK::Com_Divide(double a1[2], double a2[2], double res[2])
258{
259 double tmp = a2[0]*a2[0]+a2[1]*a2[1];
260 res[0] = (a1[0]*a2[0]+a1[1]*a2[1])/tmp;
261 res[1] = (a1[1]*a2[0]-a1[0]*a2[1])/tmp;
262}
263//------------base---------------------------------
264double EvtD0ToKSLKK::SCADot(double a1[4], double a2[4])
265{
266 double _cal = a1[0]*a2[0]-a1[1]*a2[1]-a1[2]*a2[2]-a1[3]*a2[3];
267 return _cal;
268}
269double EvtD0ToKSLKK::barrier(int l, double sa, double sb, double sc, double r, double mass)
270{
271 double q = (sa+sb-sc)*(sa+sb-sc)/(4*sa)-sb;
272 //if(q < 0) q = 1e-16;
273 if(q < 0) q = -q;
274 double z;
275 z = q*r*r;
276 double sa0;
277 sa0 = mass*mass;
278 double q0 = (sa0+sb-sc)*(sa0+sb-sc)/(4*sa0)-sb;
279 //if(q0 < 0) q0 = 1e-16;
280 if(q0 < 0) q0 = -q0;
281 double z0 = q0*r*r;
282 double F = 0.0;
283 if(l == 0) F = 1;
284 if(l == 1) F = sqrt((1+z0)/(1+z));
285 if(l == 2) F = sqrt((9+3*z0+z0*z0)/(9+3*z+z*z));
286 return F;
287}
288//------------------spin-------------------------------------------
289void EvtD0ToKSLKK::calt1(double daug1[4], double daug2[4], double t1[4])
290{
291 double p, pq, tmp;
292 double pa[4], qa[4];
293 for(int i=0; i<4; i++) {
294 pa[i] = daug1[i] + daug2[i];
295 qa[i] = daug1[i] - daug2[i];
296 }
297 p = SCADot(pa,pa);
298 pq = SCADot(pa,qa);
299 tmp = pq/p;
300 for(int i=0; i<4; i++) {
301 t1[i] = qa[i] - tmp*pa[i];
302 }
303}
304void EvtD0ToKSLKK::calt2(double daug1[4], double daug2[4], double t2[4][4])
305{
306 double p, r;
307 double pa[4], t1[4];
308 calt1(daug1,daug2,t1);
309 r = SCADot(t1,t1)/3.0;
310 for(int i=0; i<4; i++) {
311 pa[i] = daug1[i] + daug2[i];
312 }
313 p = SCADot(pa,pa);
314 for(int i=0; i<4; i++) {
315 for(int j=0; j<4; j++) {
316 t2[i][j] = t1[i]*t1[j] - r*(G[i][j]-pa[i]*pa[j]/p);
317 }
318 }
319}
320//-------------------prop--------------------------------------------
321void EvtD0ToKSLKK::propagatorCBW(double mass, double width, double sx, double prop[2])
322{
323 double a[2], b[2];
324 a[0] = 1;
325 a[1] = 0;
326 b[0] = mass*mass-sx;
327 b[1] = -mass*width;
328 Com_Divide(a,b,prop);
329}
330double EvtD0ToKSLKK::wid(double mass2, double mass, double sa, double sb, double sc, double r2, int l)
331{
332 double widm = 0.;
333 double m = sqrt(sa);
334 double tmp = sb-sc;
335 double tmp1 = sa+tmp;
336 double q = 0.25*tmp1*tmp1/sa-sb;
337 //if(q<0) q = 1e-16;
338 if(q<0) q = -q;
339 double tmp2 = mass2+tmp;
340 double q0 = 0.25*tmp2*tmp2/mass2-sb;
341 //if(q0<0) q0 = 1e-16;
342 if(q0<0) q0 = -q0;
343 double z = q*r2;
344 double z0 = q0*r2;
345 double t = q/q0;
346 if(l == 0) {widm = sqrt(t)*mass/m;}
347 else if(l == 1) {widm = t*sqrt(t)*mass/m*(1+z0)/(1+z);}
348 else if(l == 2) {widm = t*t*sqrt(t)*mass/m*(9+3*z0+z0*z0)/(9+3*z+z*z);}
349 return widm;
350}
351double EvtD0ToKSLKK::widl1(double mass2, double mass, double sa, double sb, double sc, double r2)
352{
353 double widm = 0.;
354 double m = sqrt(sa);
355 double tmp = sb-sc;
356 double tmp1 = sa+tmp;
357 double q = 0.25*tmp1*tmp1/sa-sb;
358 //if(q<0) q = 1e-16;
359 if(q<0) q = -q;
360 double tmp2 = mass2+tmp;
361 double q0 = 0.25*tmp2*tmp2/mass2-sb;
362 //if(q0<0) q0 = 1e-16;
363 if(q0<0) q0 = -q0;
364 double z = q*r2;
365 double z0 = q0*r2;
366 double F = (1+z0)/(1+z);
367 double t = q/q0;
368 widm = t*sqrt(t)*mass/m*F;
369 return widm;
370}
371void EvtD0ToKSLKK::propagatorRBW(double mass, double width, double sa, double sb, double sc, double r2, int l, double prop[2])
372{
373 double a[2], b[2];
374 double mass2 = mass*mass;
375
376 a[0] = 1;
377 a[1] = 0;
378 b[0] = mass2-sa;
379 b[1] = -mass*width*wid(mass2,mass,sa,sb,sc,r2,l);
380 Com_Divide(a,b,prop);
381}
382
383void EvtD0ToKSLKK::propagatorFlatte(double mass, double width, double sa, double prop[2]){
384
385 double q2_Pi, q2_Ka;
386 double rhoPi[2], rhoKa[2];
387
388 q2_Pi = 0.25*sa-mPi*mPi;
389 q2_Ka = 0.25*sa-mKa*mKa;
390
391 if (q2_Pi > 0) {
392 rhoPi[0] = 2.0*sqrt(q2_Pi/sa);
393 rhoPi[1] = 0.0;
394 }
395 if (q2_Pi <= 0) {
396 rhoPi[0] = 0.0;
397 rhoPi[1] = 2.0*sqrt(-q2_Pi/sa);
398 }
399
400 if (q2_Ka > 0) {
401 rhoKa[0] = 2.0*sqrt(q2_Ka/sa);
402 rhoKa[1] = 0.0;
403 }
404 if (q2_Ka <= 0) {
405 rhoKa[0] = 0.0;
406 rhoKa[1] = 2.0*sqrt(-q2_Ka/sa);
407 }
408
409 //from paper PLB 598(2004) 149 : The simga pole in J/psi -> omega pi+ pi-
410 //Double_t g1 = 0.138;
411 //Double_t g2 = 0.6141;// g1*4.45;
412 //M = 0.970
413
414 //PLB 607(2005) 243 : Resonances in J/psi -> phi pi+ pi- and phi K+ K-
415 //M = 0.965
416 //Double_t g1 = 0.165;
417 //Double_t g2 = 0.69465;// g1*4.21;
418
419 //from paper PRD 86 052006 (2012) :LHCb barB^0_s -> J/psi pi+pi-
420 //Double_t g1 = 0.199;
421 //Double_t g2 = 0.5970;// g1*3.0;
422
423 double a[2], b[2];
424 a[0] = 1;
425 a[1] = 0;
426 b[0] = mass*mass - sa + 0.165*rhoPi[1] + 0.69465*rhoKa[1];
427 b[1] = - (0.165*rhoPi[0] + 0.69465*rhoKa[0]);
428 Com_Divide(a,b,prop);
429}
430
431
432
433//------------GS---used by rho----------------------------
434void EvtD0ToKSLKK::propagatorGS(double mass, double width, double sa, double sb, double sc, double r2, double prop[2])
435{
436 double a[2], b[2];
437 double mass2 = mass*mass;
438 double tmp = sb-sc;
439 double tmp1 = sa+tmp;
440 double q2 = 0.25*tmp1*tmp1/sa-sb;
441 //if(q2<0) q2 = 1e-16;
442 if(q2<0) q2 = -q2;
443
444 double tmp2 = mass2+tmp;
445 double q02 = 0.25*tmp2*tmp2/mass2-sb;
446 //if(q02<0) q02 = 1e-16;
447 if(q02<0) q02 = -q02;
448
449 double q = sqrt(q2);
450 double q0 = sqrt(q02);
451 double m = sqrt(sa);
452 double q03 = q0*q02;
453 double tmp3 = log(mass+2*q0)+1.2760418309; // log(mass_2Pion) = 1.2760418309;
454
455 double h = GS1*q/m*(log(m+2*q)+1.2760418309);
456 double h0 = GS1*q0/mass*tmp3;
457 double dh = h0*(0.125/q02-0.5/mass2)+GS3/mass2;
458 double d = GS2/q02*tmp3+GS3*mass/q0-GS4*mass/q03;
459 double f = mass2/q03*(q2*(h-h0)+(mass2-sa)*q02*dh);
460
461 a[0] = 1.0+d*width/mass;
462 a[1] = 0.0;
463 b[0] = mass2-sa+width*f;
464 b[1] = -mass*width*widl1(mass2,mass,sa,sb,sc,r2);
465 Com_Divide(a,b,prop);
466}
467void EvtD0ToKSLKK::KPiSLASS(double sa, double sb, double sc, double prop[2]) {
468 const double m1430 = 1.441;
469 const double sa0 = 2.076481; // m1430*m1430;
470 const double w1430 = 0.193;
471 const double Lass1 = 0.25/sa0;
472 double tmp = sb-sc;
473 double tmp1 = sa0+tmp;
474 double q0 = Lass1*tmp1*tmp1-sb;
475 //if(q0<0) q0 = 1e-16;
476 if(q0<0) q0 = -q0;
477 double tmp2 = sa+tmp;
478 double qs = 0.25*tmp2*tmp2/sa-sb;
479 double q = sqrt(qs);
480 double width = w1430*q*m1430/sqrt(sa*q0);
481 double temp_R = atan(m1430*width/(sa0-sa));
482 if(temp_R<0) temp_R += math_pi;
483 double deltaR = -109.7*math_pi/180.0 + temp_R;
484 double temp_F = atan(0.226*q/(2.0-3.8194*qs)); // 2.0*0.113 = 0.226; -33.8*0.113 = -3.8194
485 if(temp_F<0) temp_F += math_pi;
486 double deltaF = 0.1*math_pi/180.0 + temp_F;
487 double deltaS = deltaR + 2.0*deltaF;
488 double t1 = 0.96*sin(deltaF);
489 double t2 = sin(deltaR);
490 double CF[2], CS[2];
491 CF[0] = cos(deltaF);
492 CF[1] = sin(deltaF);
493 CS[0] = cos(deltaS);
494 CS[1] = sin(deltaS);
495 prop[0] = t1*CF[0] + t2*CS[0];
496 prop[1] = t1*CF[1] + t2*CS[1];
497}
498
499void EvtD0ToKSLKK::Flatte_rhoab(double sa, double sb, double sc, double rho[2]){
500 double q = (sa+sb-sc)*(sa+sb-sc)/(4*sa)-sb;
501 if(q>0) {
502 rho[0]=2* sqrt(q/sa);
503 rho[1]=0;
504 }
505 else if(q<0){
506 rho[0]=0;
507 rho[1]=2*sqrt(-q/sa);
508 }
509}
510
511void EvtD0ToKSLKK::propagatorKstr1430(double mass, double sx, double *sb, double *sc, double prop[2]) //K*1430 Flatte
512{
513 double unit[2]={1.0};
514 double ci[2]={0,1};
515 double rho1[2];
516 Flatte_rhoab(sx,sb[0],sc[0],rho1);
517 double rho2[2];
518 Flatte_rhoab(sx,sb[1],sc[1],rho2);
519 double gKPi_Kstr1430=0.2990, gEtaPK_Kstr1430=0.0529;
520 double tmp1[2]={gKPi_Kstr1430,0};
521 double tmp11[2];
522 double tmp2[2]={gEtaPK_Kstr1430,0};
523 double tmp22[2];
524 Com_Multi(tmp1,rho1,tmp11);
525 Com_Multi(tmp2,rho2,tmp22);
526 double tmp3[2]={tmp11[0]+tmp22[0],tmp11[1]+tmp22[1]};
527 double tmp31[2];
528 Com_Multi(tmp3, ci,tmp31);
529 double tmp4[2]={mass*mass-sx-tmp31[0], -1.0*tmp31[1]};
530 Com_Divide( unit,tmp4, prop);
531}
532
533void EvtD0ToKSLKK::propagatora0980p(double mass, double sx, double *sb, double *sc, double prop[2]) //a0980p Flatte
534{
535 double unit[2]={1.0};
536 double ci[2]={0,1};
537 double rho1[2];
538 Flatte_rhoab(sx,sb[0],sc[0],rho1);
539 double rho2[2];
540 Flatte_rhoab(sx,sb[1],sc[1],rho2);
541 double gKK_a0980p=0.341*0.892, gEtapi_a0980p=0.341;
542 double tmp1[2]={gKK_a0980p,0};
543 double tmp11[2];
544 double tmp2[2]={gEtapi_a0980p,0};
545 double tmp22[2];
546 Com_Multi(tmp1,rho1,tmp11);
547 Com_Multi(tmp2,rho2,tmp22);
548 double tmp3[2]={tmp11[0]+tmp22[0],tmp11[1]+tmp22[1]};
549 double tmp31[2];
550 Com_Multi(tmp3, ci,tmp31);
551 double tmp4[2]={mass*mass-sx-tmp31[0], -1.0*tmp31[1]};
552 Com_Divide( unit,tmp4, prop);
553}
554
555void EvtD0ToKSLKK::propagatora0980pfloated(double mass, double sx, double *sb, double *sc, double gKK, double prop[2]) //a0980p Flatte
556{
557 double unit[2]={1.0};
558 double ci[2]={0,1};
559 double rho1[2];
560 Flatte_rhoab(sx,sb[0],sc[0],rho1);
561 double rho2[2];
562 Flatte_rhoab(sx,sb[1],sc[1],rho2);
563 double gKK_a0980p=0.341*gKK, gEtapi_a0980p=0.341;
564 double tmp1[2]={gKK_a0980p,0};
565 double tmp11[2];
566 double tmp2[2]={gEtapi_a0980p,0};
567 double tmp22[2];
568 Com_Multi(tmp1,rho1,tmp11);
569 Com_Multi(tmp2,rho2,tmp22);
570 double tmp3[2]={tmp11[0]+tmp22[0],tmp11[1]+tmp22[1]};
571 double tmp31[2];
572 Com_Multi(tmp3, ci,tmp31);
573 double tmp4[2]={mass*mass-sx-tmp31[0], -1.0*tmp31[1]};
574 Com_Divide( unit,tmp4, prop);
575}
576
577
578void EvtD0ToKSLKK::propagatora0980wm(double mass, double width, double sx, double sb, double sc, double prop[2]) //a0980wm Flatte
579{
580 double unit[2]={1.0};
581 double ci[2]={0,1};
582 double rho1[2];
583 double tmp1[2]={mass,0};
584 double tmp2[2]={width,0};
585 double tmp11[2],tmp22[2];
586 Flatte_rhoab(sx,sb,sc,rho1);
587 Com_Multi(tmp1,rho1,tmp11);
588 Com_Multi(tmp11,tmp2,tmp22);
589 double tmp3[2];
590 Com_Multi(tmp22, ci,tmp3);
591 double tmp4[2]={mass*mass-sx-tmp3[0], -1.0*tmp3[1]};
592 Com_Divide( unit,tmp4, prop);
593}
594
595void EvtD0ToKSLKK::propagatora09800(double mass, double sx, double *sb, double *sc, double prop[2]) //a09800 Flatte
596{
597 double unit[2]={1.0};
598 double ci[2]={0,1};
599 double rho1[2];
600 Flatte_rhoab(sx,sb[0],sc[0],rho1);
601 double rho2[2];
602 Flatte_rhoab(sx,sb[1],sc[1],rho2);
603 double rho3[2];
604 Flatte_rhoab(sx,sb[2],sc[2],rho3);
605 double gKK_a09800=0.341*0.892, gEtapi_a09800=0.341, gK0K0_a09800=0.341*0.892;
606 double tmp1[2]={gKK_a09800,0};
607 double tmp11[2];
608 double tmp2[2]={gEtapi_a09800,0};
609 double tmp22[2];
610 double tmp3[2]={gK0K0_a09800,0};
611 double tmp33[2];
612 Com_Multi(tmp1,rho1,tmp11);
613 Com_Multi(tmp2,rho2,tmp22);
614 Com_Multi(tmp3,rho3,tmp33);
615 double tmp4[2]={tmp11[0]+tmp22[0]+tmp33[0],tmp11[1]+tmp22[1]+tmp33[1]};
616 double tmp41[2];
617 Com_Multi(tmp4, ci,tmp41);
618 double tmp5[2]={mass*mass-sx-tmp41[0], -1.0*tmp41[1]};
619 Com_Divide( unit,tmp5, prop);
620}
621
622void EvtD0ToKSLKK::propagatora09800floated(double mass, double sx, double *sb, double *sc, double gKK, double prop[2]) //a09800 Flatte
623{
624 double unit[2]={1.0};
625 double ci[2]={0,1};
626 double rho1[2];
627 Flatte_rhoab(sx,sb[0],sc[0],rho1);
628 double rho2[2];
629 Flatte_rhoab(sx,sb[1],sc[1],rho2);
630 double rho3[2];
631 Flatte_rhoab(sx,sb[2],sc[2],rho3);
632 double gKK_a09800=0.341*gKK, gEtapi_a09800=0.341, gK0K0_a09800=0.341*gKK;
633 double tmp1[2]={gKK_a09800,0};
634 double tmp11[2];
635 double tmp2[2]={gEtapi_a09800,0};
636 double tmp22[2];
637 double tmp3[2]={gK0K0_a09800,0};
638 double tmp33[2];
639 Com_Multi(tmp1,rho1,tmp11);
640 Com_Multi(tmp2,rho2,tmp22);
641 Com_Multi(tmp3,rho3,tmp33);
642 double tmp4[2]={tmp11[0]+tmp22[0]+tmp33[0],tmp11[1]+tmp22[1]+tmp33[1]};
643 double tmp41[2];
644 Com_Multi(tmp4, ci,tmp41);
645 double tmp5[2]={mass*mass-sx-tmp41[0], -1.0*tmp41[1]};
646 Com_Divide( unit,tmp5, prop);
647}
648
649
650
651
652void EvtD0ToKSLKK::propagatora098002channel(double mass, double sx, double *sb, double *sc, double prop[2]) //a09800 Flatte
653{
654 double unit[2]={1.0};
655 double ci[2]={0,1};
656 double rho1[2];
657 Flatte_rhoab(sx,sb[0],sc[0],rho1);
658 double rho2[2];
659 Flatte_rhoab(sx,sb[1],sc[1],rho2);
660 double gKK_a09800=0.341*0.892, gEtapi_a09800=0.341;
661 double tmp1[2]={gKK_a09800,0};
662 double tmp11[2];
663 double tmp2[2]={gEtapi_a09800,0};
664 double tmp22[2];
665 Com_Multi(tmp1,rho1,tmp11);
666 Com_Multi(tmp2,rho2,tmp22);
667 double tmp4[2]={tmp11[0]+tmp22[0],tmp11[1]+tmp22[1]};
668 double tmp41[2];
669 Com_Multi(tmp4, ci,tmp41);
670 double tmp5[2]={mass*mass-sx-tmp41[0], -1.0*tmp41[1]};
671 Com_Divide( unit,tmp5, prop);
672}
673
674
675void EvtD0ToKSLKK::rhoab(double sa, double sb, double sc, double res[2]) {
676 double tmp = sa+sb-sc;
677 double q = 0.25*tmp*tmp/sa-sb;
678 if(q>=0) {
679 res[0]=2.0*sqrt(q/sa);
680 res[1]=0.0;
681 } else {
682 res[0]=0.0;
683 res[1]=2.0*sqrt(-q/sa);
684 }
685}
686void EvtD0ToKSLKK::rho4Pi(double sa, double res[2]) {
687 double temp = 1.0-0.3116765584/sa; // 0.3116765584=0.13957*0.13957*16
688 if(temp>=0) {
689 res[0]=sqrt(temp)/(1.0+exp(9.8-3.5*sa));
690 res[1]=0.0;
691 } else {
692 res[0]=0.0;
693 res[1]=sqrt(-temp)/(1.0+exp(9.8-3.5*sa));//?????????????????????
694 }
695}
696
697void EvtD0ToKSLKK::propagatorsigma500(double sa, double sb, double sc, double prop[2]) {//for pipi_s wave
698 double f = 0.5843+1.6663*sa;
699 const double M = 0.9264;
700 const double mass2 = 0.85821696; // M*M
701 const double mpi2d2 = 0.00973989245;
702 double g1 = f*(sa-mpi2d2)/(mass2-mpi2d2)*exp((mass2-sa)/1.082);
703 double rho1s[2], rho1M[2], rho2s[2], rho2M[2], rho1[2], rho2[2];
704 rhoab(sa,sb,sc,rho1s);
705 rhoab(mass2,sb,sc,rho1M);
706 rho4Pi(sa,rho2s);
707 rho4Pi(mass2,rho2M);
708 Com_Divide(rho1s,rho1M,rho1);
709 Com_Divide(rho2s,rho2M,rho2);
710 double a[2], b[2];
711 a[0] = 1.0;
712 a[1] = 0.0;
713 b[0] = mass2-sa+M*(g1*rho1[1]+0.0024*rho2[1]);
714 b[1] = -M*(g1*rho1[0]+0.0024*rho2[0]);
715 Com_Divide(a,b,prop);
716}
717
718void EvtD0ToKSLKK::getprop(double sa, double sb, double sc, double mass, double width, double prop[2]){//rho propagator
719 double prop1[2], prop2[2];
720 propagatorGS(mass,width,sa,sb,sc,9.0,prop1);
721 propagatorRBW(0.783,0.008,sa,sb,sc,3.0,1,prop2);
722 double coef_omega[2];
723 coef_omega[0] = rho_omega*cos(phi_omega),
724 coef_omega[1] = rho_omega*sin(phi_omega);
725 double one[2]; one[0] = 1; one[1] = 0;
726 double temp[2];
727 Com_Multi(coef_omega,prop2,temp);
728 temp[0] = one[0] + 0.783*0.783*temp[0];
729 temp[1] = one[1] + 0.783*0.783*temp[1];
730 Com_Multi(prop1,temp,prop);
731}
732double EvtD0ToKSLKK::DDalitz(double P1[4], double P2[4], double P3[4], int Ang, double mass){
733 double pR[4], pD[4];
734 double temp_PDF, v_re;
735 temp_PDF = 0.0;
736 v_re = 0.0;
737 double B[2], s1, s2, s3, sR, sD;
738 for(int i=0; i<4; i++){
739 pR[i] = P1[i] + P2[i];
740 pD[i] = pR[i] + P3[i];
741 }
742 s1 = SCADot(P1,P1);
743 s2 = SCADot(P2,P2);
744 s3 = SCADot(P3,P3);
745 sR = SCADot(pR,pR);
746 sD = SCADot(pD,pD);
747 int G[4][4];
748 for(int i=0; i!=4; i++){
749 for(int j=0; j!=4; j++){
750 if(i==j){
751 if(i==0) G[i][j] = 1;
752 else G[i][j] = -1;
753 }
754 else G[i][j] = 0;
755 }
756 }
757 if(Ang == 0){
758 B[0] = 1;
759 B[1] = 1;
760 temp_PDF = 1;
761 }
762 if(Ang == 1){
763 B[0] = barrier(1,sR,s1,s2,3.0,mass);
764 B[1] = barrier(1,sD,sR,s3,5.0,mDM);
765 //B[0] = Barrier(1,sR,s1,s2,9.0);
766 //B[1] = Barrier(1,sD,sR,s3,25.0);
767 double t1[4], T1[4];
768 calt1(P1,P2,t1);
769 calt1(pR,P3,T1);
770 temp_PDF = 0;
771 for(int i=0; i<4; i++){
772 temp_PDF += t1[i]*T1[i]*G[i][i];
773 }
774 }
775 if(Ang == 2){
776 B[0] = barrier(2,sR,s1,s2,3.0,mass);
777 B[1] = barrier(2,sD,sR,s3,5.0,mDM);
778 //B[0] = Barrier(2,sR,s1,s2,9.0);
779 //B[1] = Barrier(2,sD,sR,s3,25.0);
780 double t2[4][4], T2[4][4];
781 calt2(P1,P2,t2);
782 calt2(pR,P3,T2);
783 temp_PDF = 0;
784 for(int i=0; i<4; i++){
785 for(int j=0; j<4; j++){
786 temp_PDF += t2[i][j]*T2[j][i]*G[i][i]*G[j][j];
787 }
788 }
789 }
790 v_re = temp_PDF*B[0]*B[1];
791 return v_re;
792}
793
794void EvtD0ToKSLKK::calPDF(double *Ks0, double *K1, double *K2, double* mass1, double* width1,
795 double* amp, double* phase,
796 int* g0, int* spin, int* modetype,
797 double* r0, double* r1, int first, int last, double PDF[2])
798{
799 double P12[4], P23[4], P13[4];
800 double cof[2], amp_PDF[2];
801 double s12, s13, s23;
802 for(int i=0; i<4; i++){
803 P12[i] = K1[i] + Ks0[i];
804 P13[i] = K2[i] + Ks0[i];
805 P23[i] = K1[i] + K2[i];
806
807//printf("abcde: %lf, %lf, %lf, %lf, %lf\n",P12[i], P13[i], K1[i], K2[i], Ks0[i]);
808 }
809 s12 = SCADot(P12,P12);
810 s13 = SCADot(P13,P13);
811 s23 = SCADot(P23,P23);
812 double s1,s2,s3;
813 s1 = SCADot(Ks0,Ks0);
814 s2 = SCADot(K1,K1);
815 s3 = SCADot(K2,K2);
816 double pro[2], temp_PDF, amp_tmp[2];
817 double Amp_KPiS[2];
818 //double mass1sq;
819 amp_PDF[0] = 0;
820 amp_PDF[1] = 0;
821 PDF[0] = 0;
822 PDF[1] = 0;
823 amp_tmp[0] = 0;
824 amp_tmp[1] = 0;
825
826//printf("ssssss: %lf, %lf, %lf, %lf, %lf, %lf\n",s1,s2,s3,s12,s13,s23);
827
828
829 for(int i=first; i<last; i++) {
830 amp_tmp[0] = 0;
831 amp_tmp[1] = 0;
832 //mass1sq = mass1[i]*mass1[i];
833 cof[0] = amp[i]*cos(phase[i]);
834 cof[1] = amp[i]*sin(phase[i]);
835 temp_PDF = 0;
836
837
838 if(modetype[i] == 12){
839 temp_PDF = DDalitz(Ks0, K1, K2, spin[i], mass1[i]);
840 if(g0[i]==1) propagatorRBW(mass1[i],width1[i],s12,mass_KS*mass_KS,mKa2,rRes2,spin[i],pro);
841 if(g0[i]==2) { // a0980p Flatte
842 double s11[2]={mass_KS*mass_KS, mPi*mPi};
843 double s22[2]={mKa*mKa, mEta*mEta};
844 propagatora0980p(mass1[i],s12,s11,s22,pro);
845 }
846 if(g0[i]==3) propagatora0980wm(mass1[i],width1[i],s12,mass_KS*mass_KS,mKa*mKa,pro);
847 if(g0[i]==0){
848 pro[0] = 1;
849 pro[1] = 0;
850 }
851 amp_tmp[0] = temp_PDF*pro[0];
852 amp_tmp[1] = temp_PDF*pro[1];
853 }
854
855
856
857 if(modetype[i] == 13){
858 temp_PDF = DDalitz(Ks0, K2, K1, spin[i], mass1[i]);
859 if(g0[i]==1) propagatorRBW(mass1[i],width1[i],s13,mass_KS*mass_KS,mKa2,rRes2,spin[i],pro);
860 if(g0[i]==2) { // a0980p Flatte
861 double s11[2]={mass_KS*mass_KS, mPi*mPi};
862 double s22[2]={mKa*mKa, mEta*mEta};
863 propagatora0980p(mass1[i],s13,s11,s22,pro);
864 }
865 if(g0[i]==3) propagatora0980wm(mass1[i],width1[i],s13,mass_KS*mass_KS,mKa*mKa,pro);
866 if(g0[i]==4) { // a0980p Flatte
867 double s11[2]={mass_KS*mass_KS, mPi*mPi};
868 double s22[2]={mKa*mKa, mEta*mEta};
869 propagatora0980pfloated(mass1[i],s13,s11,s22,r0[i],pro);
870 }
871 if(g0[i]==5) propagatorGS(mass1[i],width1[i],s13,mass_KS*mass_KS,mKa2,rRes2,pro);
872 if(g0[i]==0){
873 pro[0] = 1;
874 pro[1] = 0;
875 }
876//printf("%lf, %lf\n",pro[0],pro[1]);
877 amp_tmp[0] = temp_PDF*pro[0];
878 amp_tmp[1] = temp_PDF*pro[1];
879
880 }
881
882 if(modetype[i] == 23){
883 temp_PDF = DDalitz(K1, K2, Ks0, spin[i], mass1[i]);
884 if(g0[i]==0){
885 pro[0] = 1;
886 pro[1] = 0;
887 }
888 if(g0[i]==1) propagatorRBW(mass1[i],width1[i],s23,mKa2,mKa2,rRes2,spin[i],pro);
889 if(g0[i]==2) { // a09800 Flatte
890 double s11[3]={mKa*mKa, mPi*mPi, mass_KS*mass_KS};
891 double s22[3]={mKa*mKa, mEta*mEta, mass_KS*mass_KS};
892 propagatora09800(mass1[i],s23,s11,s22,pro);
893 }
894 if(g0[i]==3) { // a09800 Flatte 2 channel
895 double s11[3]={mKa*mKa, mPi*mPi};
896 double s22[3]={mKa*mKa, mEta*mEta};
897 propagatora098002channel(mass1[i],s23,s11,s22,pro);
898 }
899 if(g0[i]==4){ propagatorFlatte(mass1[i],width1[i],s23,pro); }// Only for f0(980)
900 if(g0[i]==5) propagatora0980wm(mass1[i],width1[i],s23,mKa*mKa,mKa*mKa,pro);
901 if(g0[i]==6) { // a09800 Flatte
902 double s11[3]={mKa*mKa, mPi*mPi, mass_KS*mass_KS};
903 double s22[3]={mKa*mKa, mEta*mEta, mass_KS*mass_KS};
904 propagatora09800floated(mass1[i],s23,s11,s22,r0[i],pro);
905 }
906 amp_tmp[0] = temp_PDF*pro[0];
907 amp_tmp[1] = temp_PDF*pro[1];
908 }
909
910 //------------------------------------KLKK-------------------------------------------------------
911 if(modetype[i] == 232){
912 temp_PDF = DDalitz(K1, K2, Ks0, spin[i], mass1[i]);
913 if(g0[i]==0){
914 pro[0] = 1;
915 pro[1] = 0;
916 }
917 if(g0[i]==1) propagatorRBW(mass1[i],width1[i],s23,mKa2,mKa2,rRes2,spin[i],pro);
918 if(g0[i]==2) { // a09800 Flatte
919 double s11[3]={mKa*mKa, mPi*mPi, mass_KS*mass_KS};
920 double s22[3]={mKa*mKa, mEta*mEta, mass_KS*mass_KS};
921 propagatora09800(mass1[i],s23,s11,s22,pro);
922 }
923 if(g0[i]==3) { // a09800 Flatte 2 channel
924 double s11[3]={mKa*mKa, mPi*mPi};
925 double s22[3]={mKa*mKa, mEta*mEta};
926 propagatora098002channel(mass1[i],s23,s11,s22,pro);
927 }
928 if(g0[i]==4){ propagatorFlatte(mass1[i],width1[i],s23,pro); }// Only for f0(980)
929 if(g0[i]==5) propagatora0980wm(mass1[i],width1[i],s23,mKa*mKa,mKa*mKa,pro);
930 if(g0[i]==6) { // a09800 Flatte
931 double s11[3]={mKa*mKa, mPi*mPi, mass_KS*mass_KS};
932 double s22[3]={mKa*mKa, mEta*mEta, mass_KS*mass_KS};
933 propagatora09800floated(mass1[i],s23,s11,s22,r0[i],pro);
934 }
935 if(g0[i]==7) propagatorGS(mass1[i],width1[i],s23,mKa*mKa,mKa2,rRes2,pro);
936
937 double uspinbr[2], coff[2];
938 coff[0] = r0[i]*cos(r1[i]);
939 coff[1] = r0[i]*sin(r1[i]);
940 double unit[2]; unit[0] = 1.0; unit[1] = 0.0;
941
942 uspinbr[0] = unit[0] - 2*((0.22650*0.22650)/(1.-(0.22650*0.22650)))*coff[0];
943 uspinbr[1] = unit[1] - 2*((0.22650*0.22650)/(1.-(0.22650*0.22650)))*coff[1];
944
945 double amp_tmpp[2];
946 amp_tmpp[0] = temp_PDF*pro[0];
947 amp_tmpp[1] = temp_PDF*pro[1];
948
949 Com_Multi(uspinbr,amp_tmpp,amp_tmp);
950 }
951
952 if(modetype[i] == 100){
953 amp_tmp[0] = 1.0;
954 amp_tmp[1] = 0.0;
955 }
956
957
958 //if(modetype[i] == 132){
959 // KPiSLASS(s13,s1,s3,Amp_KPiS);
960 // amp_tmp[0] = Amp_KPiS[0];
961 // amp_tmp[1] = Amp_KPiS[1];
962 //}
963
964 Com_Multi(amp_tmp,cof,amp_PDF);
965 PDF[0] += amp_PDF[0];
966 PDF[1] += amp_PDF[1];
967//printf("%.10lf, %.10lf\n",PDF[0],PDF[1]);
968 }
969}
970
971
972void EvtD0ToKSLKK::calEva(double* Ks0, double* K1, double* K2, double* mass1, double* width1, double* amp, double* phase, int* g0, int* spin, int* modetype, double* r0, double* r1, double &Result, int first, int last, int charge, bool SorL)
973{
974 double PDFD0[2],PDFD0bar[2];
975 if(charge>0){
976 calPDF(Ks0,K1,K2,mass1,width1,amp,phase,g0,spin,modetype,r0,r1,first,last,PDFD0);
977 Ks0[0]= Ks0[0]; Ks0[1]= (-1.0)*Ks0[1]; Ks0[2]= (-1.0)*Ks0[2]; Ks0[3]= (-1.0)*Ks0[3];
978 K1[0] = K1[0]; K1[1] = (-1.0)*K1[1]; K1[2] = (-1.0)*K1[2]; K1[3] = (-1.0)*K1[3];
979 K2[0] = K2[0]; K2[1] = (-1.0)*K2[1]; K2[2] = (-1.0)*K2[2]; K2[3] = (-1.0)*K2[3];
980
981 //calPDF(Ks0,K2,K1,mass1,width1,amp,phase,g0,spin,modetype,r0,r1,first,last,PDFD0bar);
982 calPDF(Ks0,K1,K2,mass1,width1,amp,phase,g0,spin,modetype,r0,r1,first,last,PDFD0bar);
983 }
984 if(charge<0){
985 calPDF(Ks0,K1,K2,mass1,width1,amp,phase,g0,spin,modetype,r0,r1,first,last,PDFD0bar);
986 Ks0[0]= Ks0[0]; Ks0[1]= (-1.0)*Ks0[1]; Ks0[2]= (-1.0)*Ks0[2]; Ks0[3]= (-1.0)*Ks0[3];
987 K1[0] = K1[0]; K1[1] = (-1.0)*K1[1]; K1[2] = (-1.0)*K1[2]; K1[3] = (-1.0)*K1[3];
988 K2[0] = K2[0]; K2[1] = (-1.0)*K2[1]; K2[2] = (-1.0)*K2[2]; K2[3] = (-1.0)*K2[3];
989 //calPDF(Ks0,K2,K1,mass1,width1,amp,phase,g0,spin,modetype,r0,r1,first,last,PDFD0);
990 calPDF(Ks0,K1,K2,mass1,width1,amp,phase,g0,spin,modetype,r0,r1,first,last,PDFD0);
991 }
992/*
993 //-----------------------Quantum Correlation Correction---------------------------------
994 double r_tag = 0.0586;
995 double R_tag = 1;
996 double delta_tag = 192.1/180.0*3.1415926;
997 double qcf[2];
998 qcf[0] = r_tag*R_tag*cos(-1.0*delta_tag);
999 qcf[1] = r_tag*R_tag*sin(-1.0*delta_tag);
1000 double ampD0_part1[2], qcfampD0bar[2];
1001 Com_Multi(qcf,PDFD0bar,qcfampD0bar);
1002 if(SorL){
1003 ampD0_part1[0] = PDFD0[0] - qcfampD0bar[0];
1004 ampD0_part1[1] = PDFD0[1] - qcfampD0bar[1];
1005 } else{
1006 ampD0_part1[0] = PDFD0[0] + qcfampD0bar[0];
1007 ampD0_part1[1] = PDFD0[1] + qcfampD0bar[1];
1008 }
1009 double value = ampD0_part1[0]*ampD0_part1[0]+ampD0_part1[1]*ampD0_part1[1] + r_tag*r_tag*(1-R_tag*R_tag)*(PDFD0bar[0]*PDFD0bar[0]+PDFD0bar[1]*PDFD0bar[1]);
1010*/
1011 double value = PDFD0[0]*PDFD0[0] + PDFD0[1]*PDFD0[1];
1012 if(value <=0) value = 1e-20;
1013 Result = value;
1014
1015// printf("value: %lf\n",value);
1016
1017}
1018
double sin(const BesAngle a)
Definition BesAngle.h:210
double cos(const BesAngle a)
Definition BesAngle.h:213
double mass
TFile f("ana_bhabha660a_dqa_mcPat_zy_old.root")
TF1 * g1
EvtComplex exp(const EvtComplex &c)
double mPi
*******INTEGER m_nBinMax INTEGER m_NdiMax !No of bins in histogram for cell exploration division $ !Last vertex $ !Last active cell $ !Last cell in buffer $ !No of sampling when dividing cell $ !No of function total $ !Flag for random ceel for $ !Flag for type of for WtMax $ !Flag which decides whether vertices are included in the sampling $ entire domain is hyp !Maximum effective eevents per saves r n generator level $ !Flag for chat level in !Latex Output unit
Definition FoamA.h:90
****INTEGER imax DOUBLE PRECISION m_pi *DOUBLE PRECISION m_amfin DOUBLE PRECISION m_Chfin DOUBLE PRECISION m_Xenph DOUBLE PRECISION m_sinw2 DOUBLE PRECISION m_GFermi DOUBLE PRECISION m_MfinMin DOUBLE PRECISION m_ta2 INTEGER m_out INTEGER m_KeyFSR INTEGER m_KeyQCD *COMMON c_Semalib $ !copy of input $ !CMS energy $ !beam mass $ !final mass $ !beam charge $ !final charge $ !smallest final mass $ !Z mass $ !Z width $ !EW mixing angle $ !Gmu Fermi $ alphaQED at q
Definition KKsem.h:33
TCrossPart * CS
Definition Mcgpj.cxx:51
TTree * t
Definition binning.cxx:23
EvtDecayBase * clone()
virtual ~EvtD0ToKSLKK()
void getName(std::string &name)
void decay(EvtParticle *p)
void checkSpinDaughter(int d1, EvtSpinType::spintype sp)
void checkSpinParent(EvtSpinType::spintype sp)
double getArg(int j)
void setProbMax(double prbmx)
EvtId getParentId()
void checkNDaug(int d1, int d2=-1)
EvtId * getDaugs()
EvtId getDaug(int i)
void setProb(double prob)
static int getStdHep(EvtId id)
Definition EvtPDL.hh:56
const EvtVector4R & getP4() const
EvtParticle * getDaug(int i)
double initializePhaseSpace(int numdaughter, EvtId *daughters, double poleSize=-1., int whichTwo1=0, int whichTwo2=1)
double get(int i) const
float charge
double double double double * s12
Definition qcdloop1.h:77
double double double double double * s23
Definition qcdloop1.h:77
double precision pisqo6 one
Definition qlconstants.h:4
const double b
Definition slope.cxx:9