15 #include "EvtGenBase/EvtCPUtil.hh"
16 #include "EvtGenBase/EvtPDL.hh"
17 #include "EvtGenBase/EvtTensor4C.hh"
19 #include <generators/evtgen/EvtGenModelRegister.h>
20 #include "generators/evtgen/models/EvtB0toKsKK.h"
31 EvtB0toKsKK::~EvtB0toKsKK() {}
40 return new EvtB0toKsKK;
52 if ((getParentId() != EvtPDL::getId(
"B0") &&
53 getParentId() != EvtPDL::getId(
"anti-B0")) ||
54 (getDaug(0) != EvtPDL::getId(
"K_S0")) ||
55 (getDaug(1) != EvtPDL::getId(
"K+") &&
56 getDaug(1) != EvtPDL::getId(
"K-")) ||
57 (getDaug(2) != EvtPDL::getId(
"K+") &&
58 getDaug(2) != EvtPDL::getId(
"K-"))) {
59 std::cout <<
"ERROR: Invalid decay" << std::endl;
60 std::cout <<
"USAGE: K_S0 K+ K-" << std::endl;
65 (1.0 + getArg(3)) * EvtComplex(getArg(1) * cos(getArg(2) * M_PI / 180.0),
66 getArg(1) * sin(getArg(2) * M_PI / 180.0));
68 (1.0 + getArg(7)) * EvtComplex(getArg(5) * cos(getArg(6) * M_PI / 180.0),
69 getArg(5) * sin(getArg(6) * M_PI / 180.0));
71 (1.0 + getArg(11)) * EvtComplex(getArg(9) * cos(getArg(10) * M_PI / 180.0),
72 getArg(9) * sin(getArg(10) * M_PI / 180.0));
74 (1.0 + getArg(15)) * EvtComplex(getArg(13) * cos(getArg(14) * M_PI / 180.0),
75 getArg(13) * sin(getArg(14) * M_PI / 180.0));
77 (1.0 + getArg(19)) * EvtComplex(getArg(17) * cos(getArg(18) * M_PI / 180.0),
78 getArg(17) * sin(getArg(18) * M_PI / 180.0));
80 (1.0 + getArg(24)) * EvtComplex(getArg(22) * cos(getArg(23) * M_PI / 180.0),
81 getArg(22) * sin(getArg(23) * M_PI / 180.0));
83 (1.0 + getArg(29)) * EvtComplex(getArg(27) * cos(getArg(28) * M_PI / 180.0),
84 getArg(27) * sin(getArg(28) * M_PI / 180.0));
87 (1.0 - getArg(3)) * EvtComplex(getArg(1) * cos(getArg(2) * M_PI / 180.0),
88 getArg(1) * sin(getArg(2) * M_PI / 180.0));
90 (1.0 - getArg(7)) * EvtComplex(getArg(5) * cos(getArg(6) * M_PI / 180.0),
91 getArg(5) * sin(getArg(6) * M_PI / 180.0));
93 (1.0 - getArg(11)) * EvtComplex(getArg(9) * cos(getArg(10) * M_PI / 180.0),
94 getArg(9) * sin(getArg(10) * M_PI / 180.0));
96 (1.0 - getArg(15)) * EvtComplex(getArg(13) * cos(getArg(14) * M_PI / 180.0),
97 getArg(13) * sin(getArg(14) * M_PI / 180.0));
99 (1.0 - getArg(19)) * EvtComplex(getArg(17) * cos(getArg(18) * M_PI / 180.0),
100 getArg(17) * sin(getArg(18) * M_PI / 180.0));
102 (1.0 - getArg(24)) * EvtComplex(getArg(22) * cos(getArg(23) * M_PI / 180.0),
103 getArg(22) * sin(getArg(23) * M_PI / 180.0));
105 (1.0 - getArg(29)) * EvtComplex(getArg(27) * cos(getArg(28) * M_PI / 180.0),
106 getArg(27) * sin(getArg(28) * M_PI / 180.0));
112 std::cout << setiosflags(std::ios::left) << std::endl
113 <<
"B0 Channel " << std::setw(20) <<
"Relative amplitude" << std::setw(20) <<
"Relative phase" << std::endl
114 <<
"f0ks " << std::setw(20) << real(
a_f0ks_) << std::setw(20) << imag(
a_f0ks_) << std::endl
115 <<
"phiks " << std::setw(20) << real(
a_phiks_) << std::setw(20) << imag(
a_phiks_) << std::endl
116 <<
"fxks " << std::setw(20) << real(
a_fxks_) << std::setw(20) << imag(
a_fxks_) << std::endl
117 <<
"chic0ks " << std::setw(20) << real(
a_chic0ks_) << std::setw(20) << imag(
a_chic0ks_) << std::endl
118 <<
"kpkmnr " << std::setw(20) << real(
a_kpkmnr_) << std::setw(20) << imag(
a_kpkmnr_) << std::endl
119 <<
"kskpnr " << std::setw(20) << real(
a_kskpnr_) << std::setw(20) << imag(
a_kskpnr_) << std::endl
120 <<
"kskmnr " << std::setw(20) << real(
a_kskmnr_) << std::setw(20) << imag(
a_kskmnr_) << std::endl
123 std::cout << setiosflags(std::ios::left) << std::endl
124 <<
"B0B Channel " << std::setw(20) <<
"Relative amplitude" << std::setw(20) <<
"Relative phase" << std::endl
125 <<
"f0ks " << std::setw(20) << real(
abar_f0ks_) << std::setw(20) << imag(
abar_f0ks_) << std::endl
127 <<
"fxks " << std::setw(20) << real(
abar_fxks_) << std::setw(20) << imag(
abar_fxks_) << std::endl
145 setProbMax(1100000.0);
151 static EvtId B0 = EvtPDL::getId(
"B0");
152 static EvtId B0B = EvtPDL::getId(
"anti-B0");
159 EvtCPUtil::getInstance()->OtherB(
p, t, other_b, 0.5);
164 p->initializePhaseSpace(getNDaug(), getDaugs());
165 EvtVector4R p4ks, p4kp, p4km;
166 if (
p->getId() == B0) {
167 p4ks =
p->getDaug(0)->getP4();
168 p4kp =
p->getDaug(1)->getP4();
169 p4km =
p->getDaug(2)->getP4();
176 p4ks =
p->getDaug(0)->getP4();
177 p4kp =
p->getDaug(2)->getP4();
178 p4km =
p->getDaug(1)->getP4();
219 (1.0 + getArg(3)) * EvtComplex(getArg(1) * cos(getArg(2) * M_PI / 180.0),
220 getArg(1) * sin(getArg(2) * M_PI / 180.0));
222 (1.0 + getArg(7)) * EvtComplex(getArg(5) * cos(getArg(6) * M_PI / 180.0),
223 getArg(5) * sin(getArg(6) * M_PI / 180.0));
225 (1.0 + getArg(11)) * EvtComplex(getArg(9) * cos(getArg(10) * M_PI / 180.0),
226 getArg(9) * sin(getArg(10) * M_PI / 180.0));
228 (1.0 + getArg(15)) * EvtComplex(getArg(13) * cos(getArg(14) * M_PI / 180.0),
229 getArg(13) * sin(getArg(14) * M_PI / 180.0));
231 (1.0 + getArg(19)) * EvtComplex(getArg(17) * cos(getArg(18) * M_PI / 180.0),
232 getArg(17) * sin(getArg(18) * M_PI / 180.0));
234 (1.0 + getArg(24)) * EvtComplex(getArg(22) * cos(getArg(23) * M_PI / 180.0),
235 getArg(22) * sin(getArg(23) * M_PI / 180.0));
237 (1.0 + getArg(29)) * EvtComplex(getArg(27) * cos(getArg(28) * M_PI / 180.0),
238 getArg(27) * sin(getArg(28) * M_PI / 180.0));
241 (1.0 - getArg(3)) * EvtComplex(getArg(1) * cos(getArg(2) * M_PI / 180.0),
242 getArg(1) * sin(getArg(2) * M_PI / 180.0));
244 (1.0 - getArg(7)) * EvtComplex(getArg(5) * cos(getArg(6) * M_PI / 180.0),
245 getArg(5) * sin(getArg(6) * M_PI / 180.0));
247 (1.0 - getArg(11)) * EvtComplex(getArg(9) * cos(getArg(10) * M_PI / 180.0),
248 getArg(9) * sin(getArg(10) * M_PI / 180.0));
250 (1.0 - getArg(15)) * EvtComplex(getArg(13) * cos(getArg(14) * M_PI / 180.0),
251 getArg(13) * sin(getArg(14) * M_PI / 180.0));
253 (1.0 - getArg(19)) * EvtComplex(getArg(17) * cos(getArg(18) * M_PI / 180.0),
254 getArg(17) * sin(getArg(18) * M_PI / 180.0));
256 (1.0 - getArg(24)) * EvtComplex(getArg(22) * cos(getArg(23) * M_PI / 180.0),
257 getArg(22) * sin(getArg(23) * M_PI / 180.0));
259 (1.0 - getArg(29)) * EvtComplex(getArg(27) * cos(getArg(28) * M_PI / 180.0),
260 getArg(27) * sin(getArg(28) * M_PI / 180.0));
263 const double pCP_f0ks = getArg(4) * M_PI / 180.0;
264 const double pCP_phiks = getArg(8) * M_PI / 180.0;
265 const double pCP_fxks = getArg(12) * M_PI / 180.0;
266 const double pCP_chic0ks = getArg(16) * M_PI / 180.0;
267 const double pCP_kpkmnr = getArg(20) * M_PI / 180.0;
268 const double pCP_kskpnr = getArg(25) * M_PI / 180.0;
269 const double pCP_kskmnr = getArg(30) * M_PI / 180.0;
273 EvtComplex(cos(+2.0 * pCP_f0ks), sin(+2.0 * pCP_f0ks));
275 EvtComplex(cos(+2.0 * pCP_phiks), sin(+2.0 * pCP_phiks));
277 EvtComplex(cos(+2.0 * pCP_fxks), sin(+2.0 * pCP_fxks));
279 EvtComplex(cos(+2.0 * pCP_chic0ks), sin(+2.0 * pCP_chic0ks));
281 EvtComplex(cos(+2.0 * pCP_kpkmnr), sin(+2.0 * pCP_kpkmnr));
283 EvtComplex(cos(+2.0 * pCP_kskpnr), sin(+2.0 * pCP_kskpnr));
285 EvtComplex(cos(+2.0 * pCP_kskmnr), sin(+2.0 * pCP_kskmnr));
287 if (other_b == B0B) {
289 EvtComplex(cos(-2.0 * pCP_f0ks), sin(-2.0 * pCP_f0ks));
291 EvtComplex(cos(-2.0 * pCP_phiks), sin(-2.0 * pCP_phiks));
293 EvtComplex(cos(-2.0 * pCP_fxks), sin(-2.0 * pCP_fxks));
295 EvtComplex(cos(-2.0 * pCP_chic0ks), sin(-2.0 * pCP_chic0ks));
297 EvtComplex(cos(-2.0 * pCP_kpkmnr), sin(-2.0 * pCP_kpkmnr));
299 EvtComplex(cos(-2.0 * pCP_kskpnr), sin(-2.0 * pCP_kskpnr));
301 EvtComplex(cos(-2.0 * pCP_kskmnr), sin(-2.0 * pCP_kskmnr));
305 EvtComplex Amp_f0ks =
A_f0ks(p4ks, p4kp, p4km);
306 EvtComplex Amp_phiks =
A_phiks(p4ks, p4kp, p4km);
307 EvtComplex Amp_fxks =
A_fxks(p4ks, p4kp, p4km);
308 EvtComplex Amp_chic0ks =
A_chic0ks(p4ks, p4kp, p4km);
313 EvtComplex Ampbar_f0ks =
A_f0ks(p4ks, p4km, p4kp);
314 EvtComplex Ampbar_phiks =
A_phiks(p4ks, p4km, p4kp);
315 EvtComplex Ampbar_fxks =
A_fxks(p4ks, p4km, p4kp);
316 EvtComplex Ampbar_chic0ks =
A_chic0ks(p4ks, p4km, p4kp);
321 const EvtComplex A_B0toKsKK =
330 const EvtComplex Abar_B0toKsKK =
340 const double dm = getArg(0);
343 if (other_b == B0B) {
345 A_B0toKsKK * cos(dm * t / (2.0 * EvtConst::c)) +
346 EvtComplex(0.0, 1.0) * Abar_B0toKsKK * sin(dm * t / (2.0 * EvtConst::c));
351 EvtComplex(0.0, 1.0) * sin(dm * t / (2.0 * EvtConst::c)) +
352 Abar_B0toKsKK * cos(dm * t / (2.0 * EvtConst::c));
355 if (abs2(amp) > 50000.0)
362 EvtVector4R
EvtB0toKsKK::umu(
const EvtVector4R& p4a,
const EvtVector4R& p4b,
363 const EvtVector4R& p4c)
365 const double s = (p4a + p4b + p4c) * (p4a + p4b + p4c);
367 const EvtVector4R
umu = (p4a + p4b + p4c) / sqrt(s);
372 EvtVector4R
EvtB0toKsKK::Smu(
const EvtVector4R& p4a,
const EvtVector4R& p4b,
375 const double ma = p4a.mass();
376 const double mb = p4b.mass();
377 const double mres = (p4a + p4b).mass();
378 const double ma2 = ma * ma;
379 const double mb2 = mb * mb;
380 const double mres2 = mres * mres;
384 (sqrt(mres2 - ((ma + mb) * (ma + mb))) * sqrt(mres2 - ((ma - mb) * (ma - mb))));
386 const EvtVector4R
Smu =
387 N1 * ((p4a - p4b) - ((p4a + p4b) * (ma2 - mb2) / mres2));
392 EvtVector4R
EvtB0toKsKK::Lmu(
const EvtVector4R& p4a,
const EvtVector4R& p4b,
393 const EvtVector4R& p4c)
395 const double mc = p4c.mass();
396 const double mres = (p4a + p4b).mass();
397 const double mc2 = mc * mc;
398 const double mres2 = mres * mres;
399 const double s = (p4a + p4b + p4c) * (p4a + p4b + p4c);
403 (sqrt(s - ((mres + mc) * (mres + mc))) * sqrt(s - ((mres - mc) * (mres - mc))));
405 const EvtVector4R
Lmu =
406 N2 * ((p4c - (p4a + p4b)) - ((p4c + (p4a + p4b)) * (mc2 - mres2) / s));
412 const EvtVector4R& p4b,
413 const EvtVector4R& p4c)
415 const double s = (p4a + p4b + p4c) * (p4a + p4b + p4c);
417 const EvtVector4R
umu = (p4a + p4b + p4c) / sqrt(s);
420 EvtTensor4C::g() - EvtGenFunctions::directProd(
umu,
umu);
426 const EvtVector4R& p4c)
428 const double mres = (p4a + p4b).mass();
429 const EvtVector4R wmu = (p4a + p4b) / mres;
431 const EvtTensor4C
Tmunu =
433 (EvtGenFunctions::directProd(
Smu(p4a, p4b, p4c),
Smu(p4a, p4b, p4c)) +
434 (1.0 / 3.0 * (EvtTensor4C::g() - EvtGenFunctions::directProd(wmu, wmu))));
440 const EvtTensor4C& t2)
443 for (
unsigned int i = 0; i < 4; i++)
444 for (
unsigned int j = 0; j < 4; j++) {
445 const EvtComplex element =
446 t1.get(i, 0) * t2.get(0, j) +
447 t1.get(i, 1) * t2.get(1, j) +
448 t1.get(i, 2) * t2.get(2, j) +
449 t1.get(i, 3) * t2.get(3, j);
450 t.set(i, j, element);
457 const EvtTensor4C tmp =
Multiply(t, EvtTensor4C::g());
458 const EvtTensor4C t_raised =
Multiply(EvtTensor4C::g(), tmp);
465 vector.set(1, -vector.get(1));
466 vector.set(2, -vector.get(2));
467 vector.set(3, -vector.get(3));
471 const EvtVector4R& p4c)
473 const EvtTensor4C
Mmunu =
475 (EvtGenFunctions::directProd(
Lmu(p4a, p4b, p4c),
Lmu(p4a, p4b, p4c)) +
484 const double d = 1.0;
486 const double z =
q *
q * d * d;
491 bwbf = sqrt(2.0 * z / (1.0 + z));
493 bwbf = sqrt(13.0 * z * z / (((z - 3.0) * (z - 3.0)) + (9.0 * z)));
495 std::cout <<
"ERROR: BWBF not implemented for L>2" << std::endl;
502 const unsigned int& L)
505 const double d = 1.0;
507 const double z =
q *
q * d * d;
508 const double z0 = q0 * q0 * d * d;
513 bwbf = sqrt((1.0 + z0) / (1.0 + z));
515 bwbf = sqrt((((z0 - 3.0) * (z0 - 3.0)) + (9.0 * z0)) / (((z - 3.0) * (z - 3.0)) + (9.0 * z)));
517 std::cout <<
"ERROR: BWBF not implemented for L>2" << std::endl;
524 const double& Gamma0,
525 const double& q,
const double& q0,
526 const unsigned int& L)
528 const double s = m * m;
529 const double s0 = m0 * m0;
532 Gamma0 * m0 / m * pow(
q / q0, (2 * L) + 1) *
535 const EvtComplex
BreitWigner = 1.0 / EvtComplex(s0 - s, -m0 * Gamma);
541 const EvtVector4R& boost)
543 const double bx = boost.get(1) / boost.get(0);
544 const double by = boost.get(2) / boost.get(0);
545 const double bz = boost.get(3) / boost.get(0);
547 const double bx2 = bx * bx;
548 const double by2 = by * by;
549 const double bz2 = bz * bz;
551 const double b2 = bx2 + by2 + bz2;
556 const double gamma = 1.0 / sqrt(1 - b2);
558 const double gb2 = (gamma - 1.0) / b2;
560 const double gb2xy = gb2 * bx * by;
561 const double gb2xz = gb2 * bx * bz;
562 const double gb2yz = gb2 * by * bz;
564 const double gbx = gamma * bx;
565 const double gby = gamma * by;
566 const double gbz = gamma * bz;
569 (gamma * p4.get(0)) - (gbx * p4.get(1)) - (gby * p4.get(2)) - (gbz * p4.get(3));
571 (-gbx * p4.get(0)) + ((1.0 + (gb2 * bx2)) * p4.get(1)) +
572 (gb2xy * p4.get(2)) + (gb2xz * p4.get(3));
574 (-gby * p4.get(0)) + (gb2xy * p4.get(1)) +
575 ((1.0 + (gb2 * by2)) * p4.get(2)) + (gb2yz * p4.get(3));
577 (-gbz * p4.get(0)) + (gb2xz * p4.get(1)) +
578 (gb2yz * p4.get(2)) + ((1.0 + (gb2 * bz2)) * p4.get(3));
580 const EvtVector4R p4_b(e_b, x_b, y_b, z_b);
585 double EvtB0toKsKK::p(
const double& mab,
const double& M,
const double& mc)
587 const double sab = mab * mab;
588 const double M_p_mc2 = (M + mc) * (M + mc);
589 const double M_m_mc2 = (M - mc) * (M - mc);
591 const double p = sqrt((sab - M_p_mc2) * (sab - M_m_mc2)) / (2.0 * mab);
596 double EvtB0toKsKK::q(
const double& mab,
const double& ma,
const double& mb)
598 const double mab2 = mab * mab;
599 const double ma_p_mb2 = (ma + mb) * (ma + mb);
600 const double ma_m_mb2 = (ma - mb) * (ma - mb);
602 const double q = sqrt((mab2 - ma_p_mb2) * (mab2 - ma_m_mb2)) / (2.0 * mab);
609 const double k2 = 1.0 - (4.0 * m_h * m_h / s);
612 k = EvtComplex(0.0, sqrt(fabs(k2)));
621 const double g_pipi = 0.165;
622 const double g_kk = 4.21 * g_pipi;
624 static EvtId pip = EvtPDL::getId(
"pi+");
625 static EvtId kp = EvtPDL::getId(
"K+");
627 const double s = m * m;
628 const double s0 = m0 * m0;
630 const EvtComplex rho_pipi = 2.0 *
Flatte_k(s, EvtPDL::getMeanMass(pip)) / m;
631 const EvtComplex rho_kk = 2.0 *
Flatte_k(s, EvtPDL::getMeanMass(kp)) / m;
633 const EvtComplex Gamma =
634 EvtComplex(0.0, 1.0) * ((g_pipi * rho_pipi) + (g_kk * rho_kk));
636 const EvtComplex
Flatte = 1.0 / (s0 - s - Gamma);
642 const EvtVector4R& p4kp,
643 const EvtVector4R& p4km)
645 static EvtId f0 = EvtPDL::getId(
"f_0");
647 const double f0_m = EvtPDL::getMeanMass(f0);
649 const EvtVector4R p4kpkm = p4kp + p4km;
650 const double mkpkm = p4kpkm.mass();
653 const EvtComplex H_f0ks = 1.0;
656 const EvtComplex D_f0ks = 1.0;
659 const EvtComplex L_f0ks =
Flatte(mkpkm, f0_m);
662 const EvtComplex
A_f0ks = D_f0ks * H_f0ks * L_f0ks;
668 double s13_min(
const double& s12,
const double& M,
669 const double& m1,
const double& m2,
const double& m3)
671 const double E2star = (s12 - (m2 * m2) + (m1 * m1)) / 2.0 / sqrt(s12);
672 const double E3star = (M * M - s12 - (m3 * m3)) / 2.0 / sqrt(s12);
674 const double s23_min =
675 ((E2star + E3star) * (E2star + E3star)) -
676 ((sqrt((E2star * E2star) - (m1 * m1)) + sqrt((E3star * E3star) - (m3 * m3))) *
677 (sqrt((E2star * E2star) - (m1 * m1)) + sqrt((E3star * E3star) - (m3 * m3))));
683 double s13_max(
const double& s12,
const double& M,
684 const double& m1,
const double& m2,
const double& m3)
686 const double E2star = (s12 - (m2 * m2) + (m1 * m1)) / 2.0 / sqrt(s12);
687 const double E3star = (M * M - s12 - (m3 * m3)) / 2.0 / sqrt(s12);
689 const double s23_max =
690 ((E2star + E3star) * (E2star + E3star)) -
691 ((sqrt((E2star * E2star) - (m1 * m1)) - sqrt((E3star * E3star) - (m3 * m3))) *
692 (sqrt((E2star * E2star) - (m1 * m1)) - sqrt((E3star * E3star) - (m3 * m3))));
698 const EvtVector4R& p4kp,
699 const EvtVector4R& p4km)
701 static EvtId phi = EvtPDL::getId(
"phi");
703 const double phi_m = EvtPDL::getMeanMass(phi);
704 const double phi_w = EvtPDL::getWidth(phi);
706 const EvtVector4R p4kpkm = p4kp + p4km;
707 const double mkpkm = p4kpkm.mass();
710 const EvtVector4R S_mu =
Smu(p4kp, p4km, p4ks);
711 const EvtVector4R L_mu =
Lmu(p4kp, p4km, p4ks);
713 const EvtTensor4C g_munu_tilde =
gmunu_tilde(p4kp, p4km, p4ks);
714 const EvtTensor4C g__munu_tilde =
RaiseIndices(g_munu_tilde);
716 EvtComplex H_phiks = 0.0;
717 for (
unsigned int i = 0; i < 4; i++)
718 for (
unsigned int j = 0; j < 4; j++)
719 H_phiks += g__munu_tilde.get(i, j) * S_mu.get(i) * L_mu.get(j);
722 const EvtVector4R p4kp_kpkm =
Boost(p4kp, p4kpkm);
724 const EvtComplex D_phiks =
BWBF(p4kp_kpkm.d3mag(), 1);
727 const double q0_phi =
q(phi_m, p4kp.mass(), p4km.mass());
728 const EvtComplex L_phiks =
729 BreitWigner(mkpkm, phi_m, phi_w, p4kp_kpkm.d3mag(), q0_phi, 1);
732 const EvtComplex
A_phiks = D_phiks * H_phiks * L_phiks;
773 const EvtVector4R& p4kp,
774 const EvtVector4R& p4km)
776 static EvtId fx = EvtPDL::getId(
"f_0(1500)");
778 const double fx_m = EvtPDL::getMeanMass(fx);
779 const double fx_w = EvtPDL::getWidth(fx);
781 const EvtVector4R p4kpkm = p4kp + p4km;
782 const double mkpkm = p4kpkm.mass();
785 const EvtComplex H_fxks = 1.0;
788 const EvtComplex D_fxks = 1.0;
791 const EvtVector4R p4kp_kpkm =
Boost(p4kp, p4kpkm);
792 const double q0_fx =
q(fx_m, p4kp.mass(), p4km.mass());
793 const EvtComplex L_fxks =
794 BreitWigner(mkpkm, fx_m, fx_w, p4kp_kpkm.d3mag(), q0_fx, 0);
797 const EvtComplex
A_fxks = D_fxks * H_fxks * L_fxks;
803 const EvtVector4R& p4kp,
804 const EvtVector4R& p4km)
806 static EvtId chic0 = EvtPDL::getId(
"chi_c0");
808 const double chic0_m = EvtPDL::getMeanMass(chic0);
809 const double chic0_w = EvtPDL::getWidth(chic0);
811 const EvtVector4R p4kpkm = p4kp + p4km;
812 const double mkpkm = p4kpkm.mass();
815 const EvtComplex H_chic0ks = 1.0;
818 const EvtComplex D_chic0ks = 1.0;
821 const EvtVector4R p4kp_kpkm =
Boost(p4kp, p4kpkm);
822 const double q0_chic0 =
q(chic0_m, p4kp.mass(), p4km.mass());
823 const EvtComplex L_chic0ks =
824 BreitWigner(mkpkm, chic0_m, chic0_w, p4kp_kpkm.d3mag(), q0_chic0, 0);
827 const EvtComplex
A_chic0ks = D_chic0ks * H_chic0ks * L_chic0ks;
833 const EvtVector4R& p4k2,
834 const double& alpha_kk)
836 const EvtVector4R p4kk = p4k1 + p4k2;
837 const double m2kk = p4kk.mass2();
839 const EvtComplex
A_kknr = exp(-alpha_kk * m2kk);