AmpGen to GooFit conversion#
[1]:
from __future__ import annotations
from particle import SpinType
from decaylanguage.modeling.goofit import GooFitChain, SF_4Body
Pretty colors
[2]:
from plumbum import colors
colors.use_color = 2
[3]:
lines, all_states = GooFitChain.read_ampgen("../models/DtoKpipipi_v2.txt")
Look at an example line#
[4]:
lines[0]
[4]:
[5]:
for seen_factor in {p.spindetails() for p in lines}:
my_lines = [p for p in lines if p.spindetails() == seen_factor]
print(colors.bold | seen_factor, ":", *my_lines[0].spinfactors)
for line in my_lines:
colors.blue.print(" ", line)
DtoA1P1_A1toV2P2Dwave_V2toP3P4 : SF_4Body.DtoAP1_AtoVP2Dwave_VtoP3P4 SF_4Body.FF_123_4_L1
D0{K(1)(1270)-[D;GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}
D0{a(1)(1260)+[D;GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-}
DtoV1S2_V1toP1P2_S2toP3P4 : SF_4Body.DtoVS_VtoP1P2_StoP3P4 SF_4Body.FF_12_34_L1
D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.pole.1]{pi+,pi-}}
D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.prod.0]{pi+,pi-}}
D0{rho(770)0{pi+,pi-},KPi10[FOCUS.I32]{K-,pi+}}
D0{rho(770)0{pi+,pi-},KPi10[FOCUS.Kpi]{K-,pi+}}
DtoV1V2_V1toP1P2_V2toP3P4 : SF_4Body.DtoV1V2_V1toP1P2_V2toP3P4_S
D0{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}
D0{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}
DtoT1P1_T1toV2P2_V2toP3P4 : SF_4Body.DtoTP1_TtoVP2_VtoP3P4 SF_4Body.FF_123_4_L2
D0{K(2)*(1430)-{K*(892)~0{K-,pi+},pi-},pi+}
DtoV1V2_V1toP1P2_V2toP3P4_P : SF_4Body.DtoV1V2_V1toP1P2_V2toP3P4_P SF_4Body.FF_12_34_L1
D0[P]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}
D0[P]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}
Dtos1P1_s1toV2P2_V2toP3P4 : SF_4Body.DtoPP1_PtoVP2_VtoP3P4
D0{K(1460)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}
DtoS1S2_S1toP1P2_S2toP3P4 : SF_4Body.ONE
D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}
D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}
D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}
D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}
D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}
D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}
Dtos1P1_s1toS2P2_S2toP3P4 : SF_4Body.DtoPP1_PtoSP2_StoP3P4
D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.0]{pi+,pi-},K-},pi+}
D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.1]{pi+,pi-},K-},pi+}
D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.prod.1]{pi+,pi-},K-},pi+}
DtoA1P1_A1toV2P2_V2toP3P4 : SF_4Body.DtoAP1_AtoVP2_VtoP3P4 SF_4Body.FF_123_4_L1
D0{K(1)(1270)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}
D0{K(1)(1270)-[GSpline.EFF]{omega(782){pi+,pi-},K-},pi+}
D0{K(1)(1270)-[GSpline.EFF]{rho(1450)0{pi+,pi-},K-},pi+}
D0{K(1)(1270)-[GSpline.EFF]{rho(770)0{pi+,pi-},K-},pi+}
D0{K(1)(1400)-{K*(892)~0{K-,pi+},pi-},pi+}
D0{a(1)(1260)+[GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-}
DtoV1V2_V1toP1P2_V2toP3P4_D : SF_4Body.DtoV1V2_V1toP1P2_V2toP3P4_D SF_4Body.FF_12_34_L2
D0[D]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}
D0[D]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}
DtoA1P1_A1toS2P2_S2toP3P4 : SF_4Body.DtoAP1_AtoSP2_StoP3P4 SF_4Body.FF_123_4_L1
D0{K(1)(1270)-[GSpline.EFF]{KPi20[FOCUS.Kpi]{K-,pi+},pi-},pi+}
D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.0]{pi+,pi-},pi+},K-}
D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.1]{pi+,pi-},pi+},K-}
D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.prod.0]{pi+,pi-},pi+},K-}
[6]:
for spintype, c in zip(SpinType, colors, strict=False):
ps = [
c | format(str(p), "11")
for p in GooFitChain.all_particles
if p.spin_type == spintype
]
print(c & colors.bold | f"{spintype.name:>12}:", *ps)
Scalar: PiPi1 KPi10 PiPi0 KPi00 PiPi3 PiPi2 KPi20
PseudoScalar: pi+ K(1460)- D0 pi- K-
Vector: K*(892)~0 rho(1450)0 omega(782) rho(770)0
Axial: K(1)(1400)- K(1)(1270)- a(1)(1260)+
Tensor: K(2)*(1430)-
PseudoTensor:
Unknown:
NonDefined:
[7]:
for n, line in enumerate(lines):
print(
colors.bold | f"{n:2}",
f"{line!s:<70}",
colors.bold & colors.blue | "spinfactors:",
colors.blue | str(len(line.spinfactors)),
colors.bold & colors.magenta | "L:",
colors.magenta | "{0} [{1[0]}-{1[1]}]".format(line.L, line.L_range()),
)
0 D0[D]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}} spinfactors: 2 L: 2 [0-2.0]
1 D0[D]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}} spinfactors: 2 L: 2 [0-2.0]
2 D0[P]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}} spinfactors: 2 L: 1 [0-2.0]
3 D0[P]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}} spinfactors: 2 L: 1 [0-2.0]
4 D0{K(1)(1270)-[D;GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+} spinfactors: 2 L: 1.0 [1.0-1.0]
5 D0{K(1)(1270)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+} spinfactors: 2 L: 1.0 [1.0-1.0]
6 D0{K(1)(1270)-[GSpline.EFF]{KPi20[FOCUS.Kpi]{K-,pi+},pi-},pi+} spinfactors: 2 L: 1.0 [1.0-1.0]
7 D0{K(1)(1270)-[GSpline.EFF]{omega(782){pi+,pi-},K-},pi+} spinfactors: 2 L: 1.0 [1.0-1.0]
8 D0{K(1)(1270)-[GSpline.EFF]{rho(1450)0{pi+,pi-},K-},pi+} spinfactors: 2 L: 1.0 [1.0-1.0]
9 D0{K(1)(1270)-[GSpline.EFF]{rho(770)0{pi+,pi-},K-},pi+} spinfactors: 2 L: 1.0 [1.0-1.0]
10 D0{K(1)(1400)-{K*(892)~0{K-,pi+},pi-},pi+} spinfactors: 2 L: 1.0 [1.0-1.0]
11 D0{K(1460)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+} spinfactors: 1 L: 0 [0-0.0]
12 D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.0]{pi+,pi-},K-},pi+} spinfactors: 1 L: 0 [0-0.0]
13 D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.1]{pi+,pi-},K-},pi+} spinfactors: 1 L: 0 [0-0.0]
14 D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.prod.1]{pi+,pi-},K-},pi+} spinfactors: 1 L: 0 [0-0.0]
15 D0{K(2)*(1430)-{K*(892)~0{K-,pi+},pi-},pi+} spinfactors: 2 L: 2.0 [2.0-2.0]
16 D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.pole.1]{pi+,pi-}} spinfactors: 2 L: 1.0 [1.0-1.0]
17 D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.prod.0]{pi+,pi-}} spinfactors: 2 L: 1.0 [1.0-1.0]
18 D0{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}} spinfactors: 1 L: 0 [0-2.0]
19 D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}} spinfactors: 1 L: 0 [0-0.0]
20 D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}} spinfactors: 1 L: 0 [0-0.0]
21 D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}} spinfactors: 1 L: 0 [0-0.0]
22 D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}} spinfactors: 1 L: 0 [0-0.0]
23 D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}} spinfactors: 1 L: 0 [0-0.0]
24 D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}} spinfactors: 1 L: 0 [0-0.0]
25 D0{a(1)(1260)+[D;GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-} spinfactors: 2 L: 1.0 [1.0-1.0]
26 D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.0]{pi+,pi-},pi+},K-} spinfactors: 2 L: 1.0 [1.0-1.0]
27 D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.1]{pi+,pi-},pi+},K-} spinfactors: 2 L: 1.0 [1.0-1.0]
28 D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.prod.0]{pi+,pi-},pi+},K-} spinfactors: 2 L: 1.0 [1.0-1.0]
29 D0{a(1)(1260)+[GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-} spinfactors: 2 L: 1.0 [1.0-1.0]
30 D0{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}} spinfactors: 1 L: 0 [0-2.0]
31 D0{rho(770)0{pi+,pi-},KPi10[FOCUS.I32]{K-,pi+}} spinfactors: 2 L: 1.0 [1.0-1.0]
32 D0{rho(770)0{pi+,pi-},KPi10[FOCUS.Kpi]{K-,pi+}} spinfactors: 2 L: 1.0 [1.0-1.0]
We can make the GooFit Intro code:#
[8]:
colors.bold.print("All discovered spin configurations:")
for line in sorted({line.spindetails() for line in lines}):
print(line)
All discovered spin configurations:
DtoA1P1_A1toS2P2_S2toP3P4
DtoA1P1_A1toV2P2Dwave_V2toP3P4
DtoA1P1_A1toV2P2_V2toP3P4
DtoS1S2_S1toP1P2_S2toP3P4
DtoT1P1_T1toV2P2_V2toP3P4
DtoV1S2_V1toP1P2_S2toP3P4
DtoV1V2_V1toP1P2_V2toP3P4
DtoV1V2_V1toP1P2_V2toP3P4_D
DtoV1V2_V1toP1P2_V2toP3P4_P
Dtos1P1_s1toS2P2_S2toP3P4
Dtos1P1_s1toV2P2_V2toP3P4
[9]:
colors.bold.print("All known spin configurations:")
for e in SF_4Body:
print(e.name)
All known spin configurations:
DtoPP1_PtoSP2_StoP3P4
DtoPP1_PtoVP2_VtoP3P4
DtoV1V2_V1toP1P2_V2toP3P4_S
DtoV1V2_V1toP1P2_V2toP3P4_P
DtoV1V2_V1toP1P2_V2toP3P4_D
DtoAP1_AtoVP2_VtoP3P4
DtoAP1_AtoVP2Dwave_VtoP3P4
DtoVS_VtoP1P2_StoP3P4
DtoV1P1_V1toV2P2_V2toP3P4
DtoAP1_AtoSP2_StoP3P4
DtoTP1_TtoVP2_VtoP3P4
FF_12_34_L1
FF_12_34_L2
FF_123_4_L1
FF_123_4_L2
ONE
[10]:
print(GooFitChain.make_intro(all_states))
// Event type: D0 -> K- (0) pi+ (1) pi+ (2) pi- (3)
std::vector<std::vector<Lineshape*>> line_factor_list;
std::vector<std::vector<SpinFactor*>> spin_factor_list;
std::vector<Amplitude*> amplitudes_list;
constexpr fptype PI_PLUS { 139.57039 };
constexpr fptype PI_MINUS { 139.57039 };
constexpr fptype D_0 { 1864.84 };
constexpr fptype K_MINUS { 493.677 };
Variable PiPi1_M { "PiPi1_M" , 9990 };
Variable PiPi1_W { "PiPi1_W" , 9.9e+09 };
Variable Kst_892_0_bar_M { "Kst_892_0_bar_M" , 895.55 };
Variable Kst_892_0_bar_W { "Kst_892_0_bar_W" , 47.3 };
Variable KPi1_0_M { "KPi1_0_M" , 9990 };
Variable KPi1_0_W { "KPi1_0_W" , 9.9e+09 };
Variable rho_1450_0_M { "rho_1450_0_M" , 1465 };
Variable rho_1450_0_W { "rho_1450_0_W" , 400 };
Variable PiPi_0_M { "PiPi_0_M" , 9990 };
Variable PiPi_0_W { "PiPi_0_W" , 9.9e+09 };
Variable K_1_1400_minus_M { "K_1_1400_minus_M" , 1403 };
Variable K_1_1400_minus_W { "K_1_1400_minus_W" , 174 };
Variable omega_782_M { "omega_782_M" , 782.66 };
Variable omega_782_W { "omega_782_W" , 8.68 };
Variable KPi0_0_M { "KPi0_0_M" , 9990 };
Variable KPi0_0_W { "KPi0_0_W" , 9.9e+09 };
Variable K_1460_minus_M { "K_1460_minus_M" , 1482.4 };
Variable K_1460_minus_W { "K_1460_minus_W" , 335.6 };
Variable PiPi3_M { "PiPi3_M" , 9990 };
Variable PiPi3_W { "PiPi3_W" , 9.9e+09 };
Variable K_1_1270_minus_M { "K_1_1270_minus_M" , 1253 };
Variable K_1_1270_minus_W { "K_1_1270_minus_W" , 90 };
Variable rho_770_0_M { "rho_770_0_M" , 775.26 };
Variable rho_770_0_W { "rho_770_0_W" , 147.4 };
Variable a_1_1260_plus_M { "a_1_1260_plus_M" , 1230 };
Variable a_1_1260_plus_W { "a_1_1260_plus_W" , 420 };
Variable PiPi2_M { "PiPi2_M" , 9990 };
Variable PiPi2_W { "PiPi2_W" , 9.9e+09 };
Variable K_2st_1430_minus_M { "K_2st_1430_minus_M" , 1427.3 };
Variable K_2st_1430_minus_W { "K_2st_1430_minus_W" , 100 };
Variable KPi2_0_M { "KPi2_0_M" , 9990 };
Variable KPi2_0_W { "KPi2_0_W" , 9.9e+09 };
DK3P_DI.meson_radius = 5;
DK3P_DI.particle_masses = {D_0, K_MINUS, PI_PLUS, PI_PLUS, PI_MINUS};
[11]:
colors.green.print(" // Parameters")
print(GooFitChain.make_pars())
// Parameters
Variable D0_radius {"D0_radius", 0.0037559 };
Variable IS_p1_4pi {"IS_p1_4pi", 0.0 };
Variable IS_p1_EtaEta {"IS_p1_EtaEta", -0.39899 };
Variable IS_p1_EtapEta {"IS_p1_EtapEta", -0.34639 };
Variable IS_p1_KK {"IS_p1_KK", -0.55377 };
Variable IS_p1_mass {"IS_p1_mass", 0.651 };
Variable IS_p1_pipi {"IS_p1_pipi", 0.22889 };
Variable IS_p2_4pi {"IS_p2_4pi", 0.0 };
Variable IS_p2_EtaEta {"IS_p2_EtaEta", 0.39065 };
Variable IS_p2_EtapEta {"IS_p2_EtapEta", 0.31503 };
Variable IS_p2_KK {"IS_p2_KK", 0.55095 };
Variable IS_p2_mass {"IS_p2_mass", 1.2036 };
Variable IS_p2_pipi {"IS_p2_pipi", 0.94128 };
Variable IS_p3_4pi {"IS_p3_4pi", 0.55639 };
Variable IS_p3_EtaEta {"IS_p3_EtaEta", 0.1834 };
Variable IS_p3_EtapEta {"IS_p3_EtapEta", 0.18681 };
Variable IS_p3_KK {"IS_p3_KK", 0.23888 };
Variable IS_p3_mass {"IS_p3_mass", 1.55817 };
Variable IS_p3_pipi {"IS_p3_pipi", 0.36856 };
Variable IS_p4_4pi {"IS_p4_4pi", 0.85679 };
Variable IS_p4_EtaEta {"IS_p4_EtaEta", 0.19906 };
Variable IS_p4_EtapEta {"IS_p4_EtapEta", -0.00984 };
Variable IS_p4_KK {"IS_p4_KK", 0.40907 };
Variable IS_p4_mass {"IS_p4_mass", 1.21 };
Variable IS_p4_pipi {"IS_p4_pipi", 0.3365 };
Variable IS_p5_4pi {"IS_p5_4pi", -0.79658 };
Variable IS_p5_EtaEta {"IS_p5_EtaEta", -0.00355 };
Variable IS_p5_EtapEta {"IS_p5_EtapEta", 0.22358 };
Variable IS_p5_KK {"IS_p5_KK", -0.17558 };
Variable IS_p5_mass {"IS_p5_mass", 1.82206 };
Variable IS_p5_pipi {"IS_p5_pipi", 0.18171 };
Variable K_1_1270bar_minus_Spline_Gamma__0 {"K(1)(1270)bar-::Spline::Gamma::0", 6.62044e-09 };
Variable K_1_1270bar_minus_Spline_Gamma_1 {"K(1)(1270)bar-::Spline::Gamma::1", 8.1034e-05 };
Variable K_1_1270bar_minus_Spline_Gamma_1_0 {"K(1)(1270)bar-::Spline::Gamma::10", 0.142193 };
Variable K_1_1270bar_minus_Spline_Gamma_11 {"K(1)(1270)bar-::Spline::Gamma::11", 0.182565 };
Variable K_1_1270bar_minus_Spline_Gamma_12 {"K(1)(1270)bar-::Spline::Gamma::12", 0.231309 };
Variable K_1_1270bar_minus_Spline_Gamma_13 {"K(1)(1270)bar-::Spline::Gamma::13", 0.295139 };
Variable K_1_1270bar_minus_Spline_Gamma_14 {"K(1)(1270)bar-::Spline::Gamma::14", 0.383415 };
Variable K_1_1270bar_minus_Spline_Gamma_15 {"K(1)(1270)bar-::Spline::Gamma::15", 0.507206 };
Variable K_1_1270bar_minus_Spline_Gamma_16 {"K(1)(1270)bar-::Spline::Gamma::16", 0.673751 };
Variable K_1_1270bar_minus_Spline_Gamma_17 {"K(1)(1270)bar-::Spline::Gamma::17", 0.918665 };
Variable K_1_1270bar_minus_Spline_Gamma_18 {"K(1)(1270)bar-::Spline::Gamma::18", 1.18142 };
Variable K_1_1270bar_minus_Spline_Gamma_19 {"K(1)(1270)bar-::Spline::Gamma::19", 1.41125 };
Variable K_1_1270bar_minus_Spline_Gamma_2 {"K(1)(1270)bar-::Spline::Gamma::2", 0.000601093 };
Variable K_1_1270bar_minus_Spline_Gamma_2_0 {"K(1)(1270)bar-::Spline::Gamma::20", 1.61709 };
Variable K_1_1270bar_minus_Spline_Gamma_21 {"K(1)(1270)bar-::Spline::Gamma::21", 1.80236 };
Variable K_1_1270bar_minus_Spline_Gamma_22 {"K(1)(1270)bar-::Spline::Gamma::22", 1.97044 };
Variable K_1_1270bar_minus_Spline_Gamma_23 {"K(1)(1270)bar-::Spline::Gamma::23", 2.12449 };
Variable K_1_1270bar_minus_Spline_Gamma_24 {"K(1)(1270)bar-::Spline::Gamma::24", 2.26729 };
Variable K_1_1270bar_minus_Spline_Gamma_25 {"K(1)(1270)bar-::Spline::Gamma::25", 2.40124 };
Variable K_1_1270bar_minus_Spline_Gamma_26 {"K(1)(1270)bar-::Spline::Gamma::26", 2.52841 };
Variable K_1_1270bar_minus_Spline_Gamma_27 {"K(1)(1270)bar-::Spline::Gamma::27", 2.65063 };
Variable K_1_1270bar_minus_Spline_Gamma_28 {"K(1)(1270)bar-::Spline::Gamma::28", 2.76952 };
Variable K_1_1270bar_minus_Spline_Gamma_29 {"K(1)(1270)bar-::Spline::Gamma::29", 2.88658 };
Variable K_1_1270bar_minus_Spline_Gamma_3 {"K(1)(1270)bar-::Spline::Gamma::3", 0.00199416 };
Variable K_1_1270bar_minus_Spline_Gamma_3_0 {"K(1)(1270)bar-::Spline::Gamma::30", 3.0032 };
Variable K_1_1270bar_minus_Spline_Gamma_31 {"K(1)(1270)bar-::Spline::Gamma::31", 3.12073 };
Variable K_1_1270bar_minus_Spline_Gamma_32 {"K(1)(1270)bar-::Spline::Gamma::32", 3.24046 };
Variable K_1_1270bar_minus_Spline_Gamma_33 {"K(1)(1270)bar-::Spline::Gamma::33", 3.36373 };
Variable K_1_1270bar_minus_Spline_Gamma_34 {"K(1)(1270)bar-::Spline::Gamma::34", 3.49188 };
Variable K_1_1270bar_minus_Spline_Gamma_35 {"K(1)(1270)bar-::Spline::Gamma::35", 3.62634 };
Variable K_1_1270bar_minus_Spline_Gamma_36 {"K(1)(1270)bar-::Spline::Gamma::36", 3.76865 };
Variable K_1_1270bar_minus_Spline_Gamma_37 {"K(1)(1270)bar-::Spline::Gamma::37", 3.92048 };
Variable K_1_1270bar_minus_Spline_Gamma_38 {"K(1)(1270)bar-::Spline::Gamma::38", 4.08365 };
Variable K_1_1270bar_minus_Spline_Gamma_39 {"K(1)(1270)bar-::Spline::Gamma::39", 4.26015 };
Variable K_1_1270bar_minus_Spline_Gamma_4 {"K(1)(1270)bar-::Spline::Gamma::4", 0.00478532 };
Variable K_1_1270bar_minus_Spline_Gamma_5 {"K(1)(1270)bar-::Spline::Gamma::5", 0.0097556 };
Variable K_1_1270bar_minus_Spline_Gamma_6 {"K(1)(1270)bar-::Spline::Gamma::6", 0.0184164 };
Variable K_1_1270bar_minus_Spline_Gamma_7 {"K(1)(1270)bar-::Spline::Gamma::7", 0.0348721 };
Variable K_1_1270bar_minus_Spline_Gamma_8 {"K(1)(1270)bar-::Spline::Gamma::8", 0.0672518 };
Variable K_1_1270bar_minus_Spline_Gamma_9 {"K(1)(1270)bar-::Spline::Gamma::9", 0.105115 };
Variable K_1_1270bar_minus_mass {"K(1)(1270)bar-_mass", 1289.81, 0.557988 };
Variable K_1_1270bar_minus_radius {"K(1)(1270)bar-_radius", 0.0017 };
Variable K_1_1270bar_minus_width {"K(1)(1270)bar-_width", 116.114, 1.6492 };
Variable K_1460bar_minus_Spline_Gamma__0 {"K(1460)bar-::Spline::Gamma::0", 8.29869e-06 };
Variable K_1460bar_minus_Spline_Gamma_1 {"K(1460)bar-::Spline::Gamma::1", 0.00482634 };
Variable K_1460bar_minus_Spline_Gamma_1_0 {"K(1460)bar-::Spline::Gamma::10", 0.162536 };
Variable K_1460bar_minus_Spline_Gamma_11 {"K(1460)bar-::Spline::Gamma::11", 0.195329 };
Variable K_1460bar_minus_Spline_Gamma_12 {"K(1460)bar-::Spline::Gamma::12", 0.232727 };
Variable K_1460bar_minus_Spline_Gamma_13 {"K(1460)bar-::Spline::Gamma::13", 0.273932 };
Variable K_1460bar_minus_Spline_Gamma_14 {"K(1460)bar-::Spline::Gamma::14", 0.318217 };
Variable K_1460bar_minus_Spline_Gamma_15 {"K(1460)bar-::Spline::Gamma::15", 0.364931 };
Variable K_1460bar_minus_Spline_Gamma_16 {"K(1460)bar-::Spline::Gamma::16", 0.413518 };
Variable K_1460bar_minus_Spline_Gamma_17 {"K(1460)bar-::Spline::Gamma::17", 0.463526 };
Variable K_1460bar_minus_Spline_Gamma_18 {"K(1460)bar-::Spline::Gamma::18", 0.514629 };
Variable K_1460bar_minus_Spline_Gamma_19 {"K(1460)bar-::Spline::Gamma::19", 0.566657 };
Variable K_1460bar_minus_Spline_Gamma_2 {"K(1460)bar-::Spline::Gamma::2", 0.0176071 };
Variable K_1460bar_minus_Spline_Gamma_2_0 {"K(1460)bar-::Spline::Gamma::20", 0.619642 };
Variable K_1460bar_minus_Spline_Gamma_21 {"K(1460)bar-::Spline::Gamma::21", 0.67389 };
Variable K_1460bar_minus_Spline_Gamma_22 {"K(1460)bar-::Spline::Gamma::22", 0.730087 };
Variable K_1460bar_minus_Spline_Gamma_23 {"K(1460)bar-::Spline::Gamma::23", 0.789463 };
Variable K_1460bar_minus_Spline_Gamma_24 {"K(1460)bar-::Spline::Gamma::24", 0.853918 };
Variable K_1460bar_minus_Spline_Gamma_25 {"K(1460)bar-::Spline::Gamma::25", 0.925529 };
Variable K_1460bar_minus_Spline_Gamma_26 {"K(1460)bar-::Spline::Gamma::26", 0.999944 };
Variable K_1460bar_minus_Spline_Gamma_27 {"K(1460)bar-::Spline::Gamma::27", 1.05939 };
Variable K_1460bar_minus_Spline_Gamma_28 {"K(1460)bar-::Spline::Gamma::28", 1.11703 };
Variable K_1460bar_minus_Spline_Gamma_29 {"K(1460)bar-::Spline::Gamma::29", 1.17704 };
Variable K_1460bar_minus_Spline_Gamma_3 {"K(1460)bar-::Spline::Gamma::3", 0.0349086 };
Variable K_1460bar_minus_Spline_Gamma_3_0 {"K(1460)bar-::Spline::Gamma::30", 1.2408 };
Variable K_1460bar_minus_Spline_Gamma_31 {"K(1460)bar-::Spline::Gamma::31", 1.30884 };
Variable K_1460bar_minus_Spline_Gamma_32 {"K(1460)bar-::Spline::Gamma::32", 1.37941 };
Variable K_1460bar_minus_Spline_Gamma_33 {"K(1460)bar-::Spline::Gamma::33", 1.4487 };
Variable K_1460bar_minus_Spline_Gamma_34 {"K(1460)bar-::Spline::Gamma::34", 1.51729 };
Variable K_1460bar_minus_Spline_Gamma_35 {"K(1460)bar-::Spline::Gamma::35", 1.58528 };
Variable K_1460bar_minus_Spline_Gamma_36 {"K(1460)bar-::Spline::Gamma::36", 1.6526 };
Variable K_1460bar_minus_Spline_Gamma_37 {"K(1460)bar-::Spline::Gamma::37", 1.71901 };
Variable K_1460bar_minus_Spline_Gamma_38 {"K(1460)bar-::Spline::Gamma::38", 1.78418 };
Variable K_1460bar_minus_Spline_Gamma_39 {"K(1460)bar-::Spline::Gamma::39", 1.8476 };
Variable K_1460bar_minus_Spline_Gamma_4 {"K(1460)bar-::Spline::Gamma::4", 0.0532389 };
Variable K_1460bar_minus_Spline_Gamma_5 {"K(1460)bar-::Spline::Gamma::5", 0.0703992 };
Variable K_1460bar_minus_Spline_Gamma_6 {"K(1460)bar-::Spline::Gamma::6", 0.0855373 };
Variable K_1460bar_minus_Spline_Gamma_7 {"K(1460)bar-::Spline::Gamma::7", 0.0991271 };
Variable K_1460bar_minus_Spline_Gamma_8 {"K(1460)bar-::Spline::Gamma::8", 0.114209 };
Variable K_1460bar_minus_Spline_Gamma_9 {"K(1460)bar-::Spline::Gamma::9", 0.135205 };
Variable K_1460bar_minus_mass {"K(1460)bar-_mass", 1482.4, 3.57585 };
Variable K_1460bar_minus_radius {"K(1460)bar-_radius", 0.0017 };
Variable K_1460bar_minus_width {"K(1460)bar-_width", 335.595, 6.19588 };
Variable PiPi00_s0_prod {"PiPi00_s0_prod", -0.196872 };
Variable PiPi10_s0_prod {"PiPi10_s0_prod", -0.950027 };
Variable PiPi20_s0_prod {"PiPi20_s0_prod", -0.165753 };
Variable PiPi30_s0_prod {"PiPi30_s0_prod", -0.0676736 };
Variable a_1_1260_plus_Spline_Gamma__0 {"a(1)(1260)+::Spline::Gamma::0", 1.23936e-06 };
Variable a_1_1260_plus_Spline_Gamma_1 {"a(1)(1260)+::Spline::Gamma::1", 0.000223871 };
Variable a_1_1260_plus_Spline_Gamma_1_0 {"a(1)(1260)+::Spline::Gamma::10", 0.0832521 };
Variable a_1_1260_plus_Spline_Gamma_11 {"a(1)(1260)+::Spline::Gamma::11", 0.115406 };
Variable a_1_1260_plus_Spline_Gamma_12 {"a(1)(1260)+::Spline::Gamma::12", 0.159329 };
Variable a_1_1260_plus_Spline_Gamma_13 {"a(1)(1260)+::Spline::Gamma::13", 0.218726 };
Variable a_1_1260_plus_Spline_Gamma_14 {"a(1)(1260)+::Spline::Gamma::14", 0.295241 };
Variable a_1_1260_plus_Spline_Gamma_15 {"a(1)(1260)+::Spline::Gamma::15", 0.384295 };
Variable a_1_1260_plus_Spline_Gamma_16 {"a(1)(1260)+::Spline::Gamma::16", 0.475641 };
Variable a_1_1260_plus_Spline_Gamma_17 {"a(1)(1260)+::Spline::Gamma::17", 0.560491 };
Variable a_1_1260_plus_Spline_Gamma_18 {"a(1)(1260)+::Spline::Gamma::18", 0.635169 };
Variable a_1_1260_plus_Spline_Gamma_19 {"a(1)(1260)+::Spline::Gamma::19", 0.699435 };
Variable a_1_1260_plus_Spline_Gamma_2 {"a(1)(1260)+::Spline::Gamma::2", 0.00119329 };
Variable a_1_1260_plus_Spline_Gamma_2_0 {"a(1)(1260)+::Spline::Gamma::20", 0.754352 };
Variable a_1_1260_plus_Spline_Gamma_21 {"a(1)(1260)+::Spline::Gamma::21", 0.801255 };
Variable a_1_1260_plus_Spline_Gamma_22 {"a(1)(1260)+::Spline::Gamma::22", 0.841402 };
Variable a_1_1260_plus_Spline_Gamma_23 {"a(1)(1260)+::Spline::Gamma::23", 0.875894 };
Variable a_1_1260_plus_Spline_Gamma_24 {"a(1)(1260)+::Spline::Gamma::24", 0.905694 };
Variable a_1_1260_plus_Spline_Gamma_25 {"a(1)(1260)+::Spline::Gamma::25", 0.931693 };
Variable a_1_1260_plus_Spline_Gamma_26 {"a(1)(1260)+::Spline::Gamma::26", 0.954723 };
Variable a_1_1260_plus_Spline_Gamma_27 {"a(1)(1260)+::Spline::Gamma::27", 0.97539 };
Variable a_1_1260_plus_Spline_Gamma_28 {"a(1)(1260)+::Spline::Gamma::28", 0.994175 };
Variable a_1_1260_plus_Spline_Gamma_29 {"a(1)(1260)+::Spline::Gamma::29", 1.01148 };
Variable a_1_1260_plus_Spline_Gamma_3 {"a(1)(1260)+::Spline::Gamma::3", 0.00326416 };
Variable a_1_1260_plus_Spline_Gamma_3_0 {"a(1)(1260)+::Spline::Gamma::30", 1.02765 };
Variable a_1_1260_plus_Spline_Gamma_31 {"a(1)(1260)+::Spline::Gamma::31", 1.04297 };
Variable a_1_1260_plus_Spline_Gamma_32 {"a(1)(1260)+::Spline::Gamma::32", 1.05768 };
Variable a_1_1260_plus_Spline_Gamma_33 {"a(1)(1260)+::Spline::Gamma::33", 1.07198 };
Variable a_1_1260_plus_Spline_Gamma_34 {"a(1)(1260)+::Spline::Gamma::34", 1.08602 };
Variable a_1_1260_plus_Spline_Gamma_35 {"a(1)(1260)+::Spline::Gamma::35", 1.09995 };
Variable a_1_1260_plus_Spline_Gamma_36 {"a(1)(1260)+::Spline::Gamma::36", 1.11387 };
Variable a_1_1260_plus_Spline_Gamma_37 {"a(1)(1260)+::Spline::Gamma::37", 1.12789 };
Variable a_1_1260_plus_Spline_Gamma_38 {"a(1)(1260)+::Spline::Gamma::38", 1.14208 };
Variable a_1_1260_plus_Spline_Gamma_39 {"a(1)(1260)+::Spline::Gamma::39", 1.15651 };
Variable a_1_1260_plus_Spline_Gamma_4 {"a(1)(1260)+::Spline::Gamma::4", 0.00671647 };
Variable a_1_1260_plus_Spline_Gamma_5 {"a(1)(1260)+::Spline::Gamma::5", 0.0118496 };
Variable a_1_1260_plus_Spline_Gamma_6 {"a(1)(1260)+::Spline::Gamma::6", 0.0190462 };
Variable a_1_1260_plus_Spline_Gamma_7 {"a(1)(1260)+::Spline::Gamma::7", 0.0288353 };
Variable a_1_1260_plus_Spline_Gamma_8 {"a(1)(1260)+::Spline::Gamma::8", 0.0419745 };
Variable a_1_1260_plus_Spline_Gamma_9 {"a(1)(1260)+::Spline::Gamma::9", 0.0595699 };
Variable a_1_1260_plus_mass {"a(1)(1260)+_mass", 1195.05, 1.04514 };
Variable a_1_1260_plus_radius {"a(1)(1260)+_radius", 0.0017 };
Variable a_1_1260_plus_width {"a(1)(1260)+_width", 422.013, 2.0958 };
Variable f_scatt_0 {"f_scatt0", 0.23399 };
Variable f_scatt1 {"f_scatt1", 0.15044 };
Variable f_scatt2 {"f_scatt2", -0.20545 };
Variable f_scatt3 {"f_scatt3", 0.32825 };
Variable f_scatt4 {"f_scatt4", 0.35412 };
Variable s0_prod {"s0_prod", -1.0 };
Variable s0_scatt {"s0_scatt", -3.92637 };
Variable sA {"sA", 1.0 };
std::vector<Variable> K_1_1270bar_minus_SplineArr {{
K_1_1270bar_minus_Spline_Gamma__0,
K_1_1270bar_minus_Spline_Gamma_1,
K_1_1270bar_minus_Spline_Gamma_2,
K_1_1270bar_minus_Spline_Gamma_3,
K_1_1270bar_minus_Spline_Gamma_4,
K_1_1270bar_minus_Spline_Gamma_5,
K_1_1270bar_minus_Spline_Gamma_6,
K_1_1270bar_minus_Spline_Gamma_7,
K_1_1270bar_minus_Spline_Gamma_8,
K_1_1270bar_minus_Spline_Gamma_9,
K_1_1270bar_minus_Spline_Gamma_1_0,
K_1_1270bar_minus_Spline_Gamma_11,
K_1_1270bar_minus_Spline_Gamma_12,
K_1_1270bar_minus_Spline_Gamma_13,
K_1_1270bar_minus_Spline_Gamma_14,
K_1_1270bar_minus_Spline_Gamma_15,
K_1_1270bar_minus_Spline_Gamma_16,
K_1_1270bar_minus_Spline_Gamma_17,
K_1_1270bar_minus_Spline_Gamma_18,
K_1_1270bar_minus_Spline_Gamma_19,
K_1_1270bar_minus_Spline_Gamma_2_0,
K_1_1270bar_minus_Spline_Gamma_21,
K_1_1270bar_minus_Spline_Gamma_22,
K_1_1270bar_minus_Spline_Gamma_23,
K_1_1270bar_minus_Spline_Gamma_24,
K_1_1270bar_minus_Spline_Gamma_25,
K_1_1270bar_minus_Spline_Gamma_26,
K_1_1270bar_minus_Spline_Gamma_27,
K_1_1270bar_minus_Spline_Gamma_28,
K_1_1270bar_minus_Spline_Gamma_29,
K_1_1270bar_minus_Spline_Gamma_3_0,
K_1_1270bar_minus_Spline_Gamma_31,
K_1_1270bar_minus_Spline_Gamma_32,
K_1_1270bar_minus_Spline_Gamma_33,
K_1_1270bar_minus_Spline_Gamma_34,
K_1_1270bar_minus_Spline_Gamma_35,
K_1_1270bar_minus_Spline_Gamma_36,
K_1_1270bar_minus_Spline_Gamma_37,
K_1_1270bar_minus_Spline_Gamma_38,
K_1_1270bar_minus_Spline_Gamma_39
}};
std::vector<Variable> K_1460bar_minus_SplineArr {{
K_1460bar_minus_Spline_Gamma__0,
K_1460bar_minus_Spline_Gamma_1,
K_1460bar_minus_Spline_Gamma_2,
K_1460bar_minus_Spline_Gamma_3,
K_1460bar_minus_Spline_Gamma_4,
K_1460bar_minus_Spline_Gamma_5,
K_1460bar_minus_Spline_Gamma_6,
K_1460bar_minus_Spline_Gamma_7,
K_1460bar_minus_Spline_Gamma_8,
K_1460bar_minus_Spline_Gamma_9,
K_1460bar_minus_Spline_Gamma_1_0,
K_1460bar_minus_Spline_Gamma_11,
K_1460bar_minus_Spline_Gamma_12,
K_1460bar_minus_Spline_Gamma_13,
K_1460bar_minus_Spline_Gamma_14,
K_1460bar_minus_Spline_Gamma_15,
K_1460bar_minus_Spline_Gamma_16,
K_1460bar_minus_Spline_Gamma_17,
K_1460bar_minus_Spline_Gamma_18,
K_1460bar_minus_Spline_Gamma_19,
K_1460bar_minus_Spline_Gamma_2_0,
K_1460bar_minus_Spline_Gamma_21,
K_1460bar_minus_Spline_Gamma_22,
K_1460bar_minus_Spline_Gamma_23,
K_1460bar_minus_Spline_Gamma_24,
K_1460bar_minus_Spline_Gamma_25,
K_1460bar_minus_Spline_Gamma_26,
K_1460bar_minus_Spline_Gamma_27,
K_1460bar_minus_Spline_Gamma_28,
K_1460bar_minus_Spline_Gamma_29,
K_1460bar_minus_Spline_Gamma_3_0,
K_1460bar_minus_Spline_Gamma_31,
K_1460bar_minus_Spline_Gamma_32,
K_1460bar_minus_Spline_Gamma_33,
K_1460bar_minus_Spline_Gamma_34,
K_1460bar_minus_Spline_Gamma_35,
K_1460bar_minus_Spline_Gamma_36,
K_1460bar_minus_Spline_Gamma_37,
K_1460bar_minus_Spline_Gamma_38,
K_1460bar_minus_Spline_Gamma_39
}};
std::vector<Variable> a_1_1260_plus_SplineArr {{
a_1_1260_plus_Spline_Gamma__0,
a_1_1260_plus_Spline_Gamma_1,
a_1_1260_plus_Spline_Gamma_2,
a_1_1260_plus_Spline_Gamma_3,
a_1_1260_plus_Spline_Gamma_4,
a_1_1260_plus_Spline_Gamma_5,
a_1_1260_plus_Spline_Gamma_6,
a_1_1260_plus_Spline_Gamma_7,
a_1_1260_plus_Spline_Gamma_8,
a_1_1260_plus_Spline_Gamma_9,
a_1_1260_plus_Spline_Gamma_1_0,
a_1_1260_plus_Spline_Gamma_11,
a_1_1260_plus_Spline_Gamma_12,
a_1_1260_plus_Spline_Gamma_13,
a_1_1260_plus_Spline_Gamma_14,
a_1_1260_plus_Spline_Gamma_15,
a_1_1260_plus_Spline_Gamma_16,
a_1_1260_plus_Spline_Gamma_17,
a_1_1260_plus_Spline_Gamma_18,
a_1_1260_plus_Spline_Gamma_19,
a_1_1260_plus_Spline_Gamma_2_0,
a_1_1260_plus_Spline_Gamma_21,
a_1_1260_plus_Spline_Gamma_22,
a_1_1260_plus_Spline_Gamma_23,
a_1_1260_plus_Spline_Gamma_24,
a_1_1260_plus_Spline_Gamma_25,
a_1_1260_plus_Spline_Gamma_26,
a_1_1260_plus_Spline_Gamma_27,
a_1_1260_plus_Spline_Gamma_28,
a_1_1260_plus_Spline_Gamma_29,
a_1_1260_plus_Spline_Gamma_3_0,
a_1_1260_plus_Spline_Gamma_31,
a_1_1260_plus_Spline_Gamma_32,
a_1_1260_plus_Spline_Gamma_33,
a_1_1260_plus_Spline_Gamma_34,
a_1_1260_plus_Spline_Gamma_35,
a_1_1260_plus_Spline_Gamma_36,
a_1_1260_plus_Spline_Gamma_37,
a_1_1260_plus_Spline_Gamma_38,
a_1_1260_plus_Spline_Gamma_39
}};
std::vector<Variable> f_scatt {{
f_scatt_0,
f_scatt1,
f_scatt2,
f_scatt3,
f_scatt4
}};
std::vector<Variable> IS_poles {{
IS_p1_pipi,
IS_p1_KK,
IS_p1_4pi,
IS_p1_EtaEta,
IS_p1_EtapEta,
IS_p1_mass,
IS_p2_pipi,
IS_p2_KK,
IS_p2_4pi,
IS_p2_EtaEta,
IS_p2_EtapEta,
IS_p2_mass,
IS_p3_pipi,
IS_p3_KK,
IS_p3_4pi,
IS_p3_EtaEta,
IS_p3_EtapEta,
IS_p3_mass,
IS_p4_pipi,
IS_p4_KK,
IS_p4_4pi,
IS_p4_EtaEta,
IS_p4_EtapEta,
IS_p4_mass,
IS_p5_pipi,
IS_p5_KK,
IS_p5_4pi,
IS_p5_EtaEta,
IS_p5_EtapEta,
IS_p5_mass
}};
And the lines can be turned into code, as well:#
[12]:
colors.green.print(" // Lines")
for n, line in enumerate(lines):
colors.green.print(" // Line", n)
print(line.to_goofit(all_states[1:]), end="\n\n\n")
// Lines
// Line 0
// D0[D]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_D , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::FF_12_34_L2 , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_D , 0, 2, 1, 3),
new SpinFactor("SF", SF_4Body::FF_12_34_L2 , 0, 2, 1, 3)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_34, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_24, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0[D]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}",
mkvar("D0[D]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}_r", true, 1.0, 0.0),
mkvar("D0[D]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}_i", true, 0.0, 0.0),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 1
// D0[D]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_D , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::FF_12_34_L2 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_D , 2, 3, 0, 1),
new SpinFactor("SF", SF_4Body::FF_12_34_L2 , 2, 3, 0, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("rho(1450)0", rho_1450_0_M, rho_1450_0_W, 1.0, M_24, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2),
new Lineshapes::RBW("rho(1450)0", rho_1450_0_M, rho_1450_0_W, 1.0, M_34, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0[D]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}",
mkvar("D0[D]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}_r", false, 0.625141, 0.0138651),
mkvar("D0[D]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}_i", false, -0.174115, 0.0158658),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 2
// D0[P]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_P , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_P , 0, 2, 1, 3),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 0, 2, 1, 3)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_34, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_24, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0[P]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}",
mkvar("D0[P]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}_r", false, -0.080874, -0.00287233),
mkvar("D0[P]{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}_i", false, -0.35291, -0.00284993),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 3
// D0[P]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_P , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_P , 2, 3, 0, 1),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 2, 3, 0, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("rho(1450)0", rho_1450_0_M, rho_1450_0_W, 1.0, M_24, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2),
new Lineshapes::RBW("rho(1450)0", rho_1450_0_M, rho_1450_0_W, 1.0, M_34, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0[P]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}",
mkvar("D0[P]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}_r", false, -0.0817223, 0.00501346),
mkvar("D0[P]{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}_i", false, 0.637565, 0.00491896),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 4
// D0{K(1)(1270)-[D;GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2Dwave_VtoP3P4 , 0, 1, 3, 2),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 0, 1, 3, 2),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2Dwave_VtoP3P4 , 0, 2, 3, 1),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 0, 2, 3, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 2, M_12_4, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2),
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 2, M_13_4, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(1)(1270)-[D;GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}",
mkvar("D0{K(1)(1270)-[D;GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}_r", false, -0.148416, 0.00282651),
mkvar("D0{K(1)(1270)-[D;GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}_i", false, 0.330131, 0.00147999),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 5
// D0{K(1)(1270)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 0, 1, 3, 2),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 0, 1, 3, 2),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 0, 2, 3, 1),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 0, 2, 3, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 0, M_12_4, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2),
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 0, M_13_4, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(1)(1270)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}",
mkvar("D0{K(1)(1270)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}_r", false, -0.148416, 0.00282651),
mkvar("D0{K(1)(1270)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}_i", false, 0.330131, 0.00147999),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 6
// D0{K(1)(1270)-[GSpline.EFF]{KPi20[FOCUS.Kpi]{K-,pi+},pi-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoSP2_StoP3P4 , 0, 1, 3, 2),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 0, 1, 3, 2),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoSP2_StoP3P4 , 0, 2, 3, 1),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 0, 2, 3, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 1.0, M_12_4, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::FOCUS("KPi20", Lineshapes::FOCUS::Mod::Kpi, KPi2_0_M, KPi2_0_W, 0, M_12, FF::BL2, 1.5),
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 1.0, M_13_4, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::FOCUS("KPi20", Lineshapes::FOCUS::Mod::Kpi, KPi2_0_M, KPi2_0_W, 0, M_13, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(1)(1270)-[GSpline.EFF]{KPi20[FOCUS.Kpi]{K-,pi+},pi-},pi+}",
mkvar("D0{K(1)(1270)-[GSpline.EFF]{KPi20[FOCUS.Kpi]{K-,pi+},pi-},pi+}_r", false, -0.148416, 0.00282651),
mkvar("D0{K(1)(1270)-[GSpline.EFF]{KPi20[FOCUS.Kpi]{K-,pi+},pi-},pi+}_i", false, 0.330131, 0.00147999),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 7
// D0{K(1)(1270)-[GSpline.EFF]{omega(782){pi+,pi-},K-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 2, 3, 0, 1),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 2, 3, 0, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 0, M_24_1, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("omega(782)0", omega_782_M, omega_782_W, 1.0, M_24, FF::BL2),
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 0, M_34_1, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("omega(782)0", omega_782_M, omega_782_W, 1.0, M_34, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(1)(1270)-[GSpline.EFF]{omega(782){pi+,pi-},K-},pi+}",
mkvar("D0{K(1)(1270)-[GSpline.EFF]{omega(782){pi+,pi-},K-},pi+}_r", false, -0.148416, 0.00282651),
mkvar("D0{K(1)(1270)-[GSpline.EFF]{omega(782){pi+,pi-},K-},pi+}_i", false, 0.330131, 0.00147999),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 8
// D0{K(1)(1270)-[GSpline.EFF]{rho(1450)0{pi+,pi-},K-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 2, 3, 0, 1),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 2, 3, 0, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 0, M_24_1, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("rho(1450)0", rho_1450_0_M, rho_1450_0_W, 1.0, M_24, FF::BL2),
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 0, M_34_1, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("rho(1450)0", rho_1450_0_M, rho_1450_0_W, 1.0, M_34, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(1)(1270)-[GSpline.EFF]{rho(1450)0{pi+,pi-},K-},pi+}",
mkvar("D0{K(1)(1270)-[GSpline.EFF]{rho(1450)0{pi+,pi-},K-},pi+}_r", false, -0.148416, 0.00282651),
mkvar("D0{K(1)(1270)-[GSpline.EFF]{rho(1450)0{pi+,pi-},K-},pi+}_i", false, 0.330131, 0.00147999),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 9
// D0{K(1)(1270)-[GSpline.EFF]{rho(770)0{pi+,pi-},K-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 2, 3, 0, 1),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 2, 3, 0, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 0, M_24_1, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_24, FF::BL2),
new Lineshapes::GSpline("K(1)(1270)bar-", K_1_1270_minus_M, K_1_1270_minus_W, 0, M_34_1, FF::BL2,
1.5, K_1_1270bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_34, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(1)(1270)-[GSpline.EFF]{rho(770)0{pi+,pi-},K-},pi+}",
mkvar("D0{K(1)(1270)-[GSpline.EFF]{rho(770)0{pi+,pi-},K-},pi+}_r", false, -0.148416, 0.00282651),
mkvar("D0{K(1)(1270)-[GSpline.EFF]{rho(770)0{pi+,pi-},K-},pi+}_i", false, 0.330131, 0.00147999),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 10
// D0{K(1)(1400)-{K*(892)~0{K-,pi+},pi-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 0, 1, 3, 2),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 0, 1, 3, 2),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 0, 2, 3, 1),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 0, 2, 3, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("K(1)(1400)bar-", K_1_1400_minus_M, K_1_1400_minus_W, 0, M_12_4, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2),
new Lineshapes::RBW("K(1)(1400)bar-", K_1_1400_minus_M, K_1_1400_minus_W, 0, M_13_4, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(1)(1400)-{K*(892)~0{K-,pi+},pi-},pi+}",
mkvar("D0{K(1)(1400)-{K*(892)~0{K-,pi+},pi-},pi+}_r", false, -0.125471, -0.00260923),
mkvar("D0{K(1)(1400)-{K*(892)~0{K-,pi+},pi-},pi+}_i", false, -0.0225264, -0.00280319),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 11
// D0{K(1460)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoPP1_PtoVP2_VtoP3P4 , 0, 1, 3, 2),
new SpinFactor("SF", SF_4Body::DtoPP1_PtoVP2_VtoP3P4 , 0, 2, 3, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("K(1460)bar-", K_1460_minus_M, K_1460_minus_W, 1.0, M_12_4, FF::BL2,
1.5, K_1460bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2),
new Lineshapes::GSpline("K(1460)bar-", K_1460_minus_M, K_1460_minus_W, 1.0, M_13_4, FF::BL2,
1.5, K_1460bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(1460)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}",
mkvar("D0{K(1460)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}_r", false, -0.120933, -0.0017034),
mkvar("D0{K(1460)-[GSpline.EFF]{K*(892)~0{K-,pi+},pi-},pi+}_i", false, 0.0155464, -0.00440427),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 12
// D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.0]{pi+,pi-},K-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoPP1_PtoSP2_StoP3P4 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::DtoPP1_PtoSP2_StoP3P4 , 2, 3, 0, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("K(1460)bar-", K_1460_minus_M, K_1460_minus_W, 0, M_24_1, FF::BL2,
1.5, K_1460bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::kMatrix("PiPi30", 0, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi3_M, PiPi3_W, 0, M_24, FF::BL2, 1.5),
new Lineshapes::GSpline("K(1460)bar-", K_1460_minus_M, K_1460_minus_W, 0, M_34_1, FF::BL2,
1.5, K_1460bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::kMatrix("PiPi30", 0, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi3_M, PiPi3_W, 0, M_34, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.0]{pi+,pi-},K-},pi+}",
mkvar("D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.0]{pi+,pi-},K-},pi+}_r", false, -0.120933, -0.0017034),
mkvar("D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.0]{pi+,pi-},K-},pi+}_i", false, 0.0155464, -0.00440427),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 13
// D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.1]{pi+,pi-},K-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoPP1_PtoSP2_StoP3P4 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::DtoPP1_PtoSP2_StoP3P4 , 2, 3, 0, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("K(1460)bar-", K_1460_minus_M, K_1460_minus_W, 0, M_24_1, FF::BL2,
1.5, K_1460bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::kMatrix("PiPi30", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi3_M, PiPi3_W, 0, M_24, FF::BL2, 1.5),
new Lineshapes::GSpline("K(1460)bar-", K_1460_minus_M, K_1460_minus_W, 0, M_34_1, FF::BL2,
1.5, K_1460bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::kMatrix("PiPi30", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi3_M, PiPi3_W, 0, M_34, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.1]{pi+,pi-},K-},pi+}",
mkvar("D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.1]{pi+,pi-},K-},pi+}_r", false, -0.120933, -0.0017034),
mkvar("D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.pole.1]{pi+,pi-},K-},pi+}_i", false, 0.0155464, -0.00440427),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 14
// D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.prod.1]{pi+,pi-},K-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoPP1_PtoSP2_StoP3P4 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::DtoPP1_PtoSP2_StoP3P4 , 2, 3, 0, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("K(1460)bar-", K_1460_minus_M, K_1460_minus_W, 0, M_24_1, FF::BL2,
1.5, K_1460bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::kMatrix("PiPi30", 1, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi3_M, PiPi3_W, 0, M_24, FF::BL2, 1.5),
new Lineshapes::GSpline("K(1460)bar-", K_1460_minus_M, K_1460_minus_W, 0, M_34_1, FF::BL2,
1.5, K_1460bar_minus_SplineArr, Lineshapes::spline_t(0.6,3.0,40)),
new Lineshapes::kMatrix("PiPi30", 1, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi3_M, PiPi3_W, 0, M_34, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.prod.1]{pi+,pi-},K-},pi+}",
mkvar("D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.prod.1]{pi+,pi-},K-},pi+}_r", false, -0.120933, -0.0017034),
mkvar("D0{K(1460)-[GSpline.EFF]{PiPi3[kMatrix.prod.1]{pi+,pi-},K-},pi+}_i", false, 0.0155464, -0.00440427),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 15
// D0{K(2)*(1430)-{K*(892)~0{K-,pi+},pi-},pi+}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoTP1_TtoVP2_VtoP3P4 , 0, 1, 3, 2),
new SpinFactor("SF", SF_4Body::FF_123_4_L2 , 0, 1, 3, 2),
new SpinFactor("SF", SF_4Body::DtoTP1_TtoVP2_VtoP3P4 , 0, 2, 3, 1),
new SpinFactor("SF", SF_4Body::FF_123_4_L2 , 0, 2, 3, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("K(2)*(1430)bar-", K_2st_1430_minus_M, K_2st_1430_minus_W, 1.0, M_12_4, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2),
new Lineshapes::RBW("K(2)*(1430)bar-", K_2st_1430_minus_M, K_2st_1430_minus_W, 1.0, M_13_4, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{K(2)*(1430)-{K*(892)~0{K-,pi+},pi-},pi+}",
mkvar("D0{K(2)*(1430)-{K*(892)~0{K-,pi+},pi-},pi+}_r", false, 0.0643742, -0.00291957),
mkvar("D0{K(2)*(1430)-{K*(892)~0{K-,pi+},pi-},pi+}_i", false, -0.294991, -0.00307298),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 16
// D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.pole.1]{pi+,pi-}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoVS_VtoP1P2_StoP3P4 , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::DtoVS_VtoP1P2_StoP3P4 , 0, 2, 1, 3),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 0, 2, 1, 3)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2),
new Lineshapes::kMatrix("PiPi10", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi1_M, PiPi1_W, 0, M_34, FF::BL2, 1.5),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2),
new Lineshapes::kMatrix("PiPi10", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi1_M, PiPi1_W, 0, M_24, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.pole.1]{pi+,pi-}}",
mkvar("D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.pole.1]{pi+,pi-}}_r", true, 1.0, 0.0),
mkvar("D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.pole.1]{pi+,pi-}}_i", true, 0.0, 0.0),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 17
// D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.prod.0]{pi+,pi-}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoVS_VtoP1P2_StoP3P4 , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::DtoVS_VtoP1P2_StoP3P4 , 0, 2, 1, 3),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 0, 2, 1, 3)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2),
new Lineshapes::kMatrix("PiPi10", 0, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi1_M, PiPi1_W, 0, M_34, FF::BL2, 1.5),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2),
new Lineshapes::kMatrix("PiPi10", 0, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi1_M, PiPi1_W, 0, M_24, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.prod.0]{pi+,pi-}}",
mkvar("D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.prod.0]{pi+,pi-}}_r", true, 1.0, 0.0),
mkvar("D0{K*(892)~0{K-,pi+},PiPi1[kMatrix.prod.0]{pi+,pi-}}_i", true, 0.0, 0.0),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 18
// D0{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_S , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_S , 0, 2, 1, 3)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_34, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_24, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}",
mkvar("D0{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}_r", false, 0.181293, 0.00065235),
mkvar("D0{K*(892)~0{K-,pi+},rho(770)0{pi+,pi-}}_i", false, -0.0745874, 0.000680398),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 19
// D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::ONE , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::ONE , 0, 2, 1, 3)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::I32, KPi0_0_M, KPi0_0_W, 0, M_12, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_34, FF::BL2, 1.5),
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::I32, KPi0_0_M, KPi0_0_W, 0, M_13, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_24, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}",
mkvar("D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}_r", true, 1.0, 0.0),
mkvar("D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}_i", true, 0.0, 0.0),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 20
// D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::ONE , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::ONE , 0, 2, 1, 3)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::I32, KPi0_0_M, KPi0_0_W, 0, M_12, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 0, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_34, FF::BL2, 1.5),
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::I32, KPi0_0_M, KPi0_0_W, 0, M_13, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 0, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_24, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}",
mkvar("D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}_r", true, 1.0, 0.0),
mkvar("D0{KPi00[FOCUS.I32]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}_i", true, 0.0, 0.0),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 21
// D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::ONE , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::ONE , 0, 2, 1, 3)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::KEta, KPi0_0_M, KPi0_0_W, 0, M_12, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_34, FF::BL2, 1.5),
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::KEta, KPi0_0_M, KPi0_0_W, 0, M_13, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_24, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}",
mkvar("D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}_r", true, 1.0, 0.0),
mkvar("D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}_i", true, 0.0, 0.0),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 22
// D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::ONE , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::ONE , 0, 2, 1, 3)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::KEta, KPi0_0_M, KPi0_0_W, 0, M_12, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 0, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_34, FF::BL2, 1.5),
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::KEta, KPi0_0_M, KPi0_0_W, 0, M_13, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 0, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_24, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}",
mkvar("D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}_r", true, 1.0, 0.0),
mkvar("D0{KPi00[FOCUS.KEta]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}_i", true, 0.0, 0.0),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 23
// D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::ONE , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::ONE , 0, 2, 1, 3)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::Kpi, KPi0_0_M, KPi0_0_W, 0, M_12, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_34, FF::BL2, 1.5),
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::Kpi, KPi0_0_M, KPi0_0_W, 0, M_13, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_24, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}",
mkvar("D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}_r", true, 1.0, 0.0),
mkvar("D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.pole.1]{pi+,pi-}}_i", true, 0.0, 0.0),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 24
// D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::ONE , 0, 1, 2, 3),
new SpinFactor("SF", SF_4Body::ONE , 0, 2, 1, 3)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::Kpi, KPi0_0_M, KPi0_0_W, 0, M_12, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 0, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_34, FF::BL2, 1.5),
new Lineshapes::FOCUS("KPi00", Lineshapes::FOCUS::Mod::Kpi, KPi0_0_M, KPi0_0_W, 0, M_13, FF::BL2, 1.5),
new Lineshapes::kMatrix("PiPi00", 0, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi_0_M, PiPi_0_W, 0, M_24, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}",
mkvar("D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}_r", true, 1.0, 0.0),
mkvar("D0{KPi00[FOCUS.Kpi]{K-,pi+},PiPi0[kMatrix.prod.0]{pi+,pi-}}_i", true, 0.0, 0.0),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 25
// D0{a(1)(1260)+[D;GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2Dwave_VtoP3P4 , 1, 3, 2, 0),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 1, 3, 2, 0),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2Dwave_VtoP3P4 , 2, 3, 1, 0),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 2, 3, 1, 0)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("a(1)(1260)+", a_1_1260_plus_M, a_1_1260_plus_W, 2, M_24_3, FF::BL2,
1.5, a_1_1260_plus_SplineArr, Lineshapes::spline_t(0.18412,1.9,40)),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_24, FF::BL2),
new Lineshapes::GSpline("a(1)(1260)+", a_1_1260_plus_M, a_1_1260_plus_W, 2, M_34_2, FF::BL2,
1.5, a_1_1260_plus_SplineArr, Lineshapes::spline_t(0.18412,1.9,40)),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_34, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{a(1)(1260)+[D;GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-}",
mkvar("D0{a(1)(1260)+[D;GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-}_r", false, -0.698394, -0.00833039),
mkvar("D0{a(1)(1260)+[D;GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-}_i", false, -0.417067, -0.00852572),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 26
// D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.0]{pi+,pi-},pi+},K-}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoSP2_StoP3P4 , 1, 3, 2, 0),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 1, 3, 2, 0),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoSP2_StoP3P4 , 2, 3, 1, 0),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 2, 3, 1, 0)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("a(1)(1260)+", a_1_1260_plus_M, a_1_1260_plus_W, 1.0, M_24_3, FF::BL2,
1.5, a_1_1260_plus_SplineArr, Lineshapes::spline_t(0.18412,1.9,40)),
new Lineshapes::kMatrix("PiPi20", 0, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi2_M, PiPi2_W, 0, M_24, FF::BL2, 1.5),
new Lineshapes::GSpline("a(1)(1260)+", a_1_1260_plus_M, a_1_1260_plus_W, 1.0, M_34_2, FF::BL2,
1.5, a_1_1260_plus_SplineArr, Lineshapes::spline_t(0.18412,1.9,40)),
new Lineshapes::kMatrix("PiPi20", 0, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi2_M, PiPi2_W, 0, M_34, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.0]{pi+,pi-},pi+},K-}",
mkvar("D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.0]{pi+,pi-},pi+},K-}_r", false, -0.698394, -0.00833039),
mkvar("D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.0]{pi+,pi-},pi+},K-}_i", false, -0.417067, -0.00852572),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 27
// D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.1]{pi+,pi-},pi+},K-}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoSP2_StoP3P4 , 1, 3, 2, 0),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 1, 3, 2, 0),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoSP2_StoP3P4 , 2, 3, 1, 0),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 2, 3, 1, 0)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("a(1)(1260)+", a_1_1260_plus_M, a_1_1260_plus_W, 1.0, M_24_3, FF::BL2,
1.5, a_1_1260_plus_SplineArr, Lineshapes::spline_t(0.18412,1.9,40)),
new Lineshapes::kMatrix("PiPi20", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi2_M, PiPi2_W, 0, M_24, FF::BL2, 1.5),
new Lineshapes::GSpline("a(1)(1260)+", a_1_1260_plus_M, a_1_1260_plus_W, 1.0, M_34_2, FF::BL2,
1.5, a_1_1260_plus_SplineArr, Lineshapes::spline_t(0.18412,1.9,40)),
new Lineshapes::kMatrix("PiPi20", 1, true,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi2_M, PiPi2_W, 0, M_34, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.1]{pi+,pi-},pi+},K-}",
mkvar("D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.1]{pi+,pi-},pi+},K-}_r", false, -0.698394, -0.00833039),
mkvar("D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.pole.1]{pi+,pi-},pi+},K-}_i", false, -0.417067, -0.00852572),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 28
// D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.prod.0]{pi+,pi-},pi+},K-}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoSP2_StoP3P4 , 1, 3, 2, 0),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 1, 3, 2, 0),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoSP2_StoP3P4 , 2, 3, 1, 0),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 2, 3, 1, 0)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("a(1)(1260)+", a_1_1260_plus_M, a_1_1260_plus_W, 1.0, M_24_3, FF::BL2,
1.5, a_1_1260_plus_SplineArr, Lineshapes::spline_t(0.18412,1.9,40)),
new Lineshapes::kMatrix("PiPi20", 0, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi2_M, PiPi2_W, 0, M_24, FF::BL2, 1.5),
new Lineshapes::GSpline("a(1)(1260)+", a_1_1260_plus_M, a_1_1260_plus_W, 1.0, M_34_2, FF::BL2,
1.5, a_1_1260_plus_SplineArr, Lineshapes::spline_t(0.18412,1.9,40)),
new Lineshapes::kMatrix("PiPi20", 0, false,
sA_0, sA, s0_prod, s0_scatt,
f_scatt, IS_poles,
PiPi2_M, PiPi2_W, 0, M_34, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.prod.0]{pi+,pi-},pi+},K-}",
mkvar("D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.prod.0]{pi+,pi-},pi+},K-}_r", false, -0.698394, -0.00833039),
mkvar("D0{a(1)(1260)+[GSpline.EFF]{PiPi2[kMatrix.prod.0]{pi+,pi-},pi+},K-}_i", false, -0.417067, -0.00852572),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 29
// D0{a(1)(1260)+[GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 1, 3, 2, 0),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 1, 3, 2, 0),
new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4 , 2, 3, 1, 0),
new SpinFactor("SF", SF_4Body::FF_123_4_L1 , 2, 3, 1, 0)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::GSpline("a(1)(1260)+", a_1_1260_plus_M, a_1_1260_plus_W, 0, M_24_3, FF::BL2,
1.5, a_1_1260_plus_SplineArr, Lineshapes::spline_t(0.18412,1.9,40)),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_24, FF::BL2),
new Lineshapes::GSpline("a(1)(1260)+", a_1_1260_plus_M, a_1_1260_plus_W, 0, M_34_2, FF::BL2,
1.5, a_1_1260_plus_SplineArr, Lineshapes::spline_t(0.18412,1.9,40)),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_34, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{a(1)(1260)+[GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-}",
mkvar("D0{a(1)(1260)+[GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-}_r", false, -0.698394, -0.00833039),
mkvar("D0{a(1)(1260)+[GSpline.EFF]{rho(770)0{pi+,pi-},pi+},K-}_i", false, -0.417067, -0.00852572),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 30
// D0{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_S , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_S , 2, 3, 0, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("rho(1450)0", rho_1450_0_M, rho_1450_0_W, 1.0, M_24, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_13, FF::BL2),
new Lineshapes::RBW("rho(1450)0", rho_1450_0_M, rho_1450_0_W, 1.0, M_34, FF::BL2),
new Lineshapes::RBW("K*(892)bar0", Kst_892_0_bar_M, Kst_892_0_bar_W, 1.0, M_12, FF::BL2)
});
amplitudes_list.push_back(new Amplitude{
"D0{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}",
mkvar("D0{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}_r", false, 0.0109497, -0.00487185),
mkvar("D0{rho(1450)0{pi+,pi-},K*(892)~0{K-,pi+}}_i", false, -0.161547, -0.00481352),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 31
// D0{rho(770)0{pi+,pi-},KPi10[FOCUS.I32]{K-,pi+}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoVS_VtoP1P2_StoP3P4 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::DtoVS_VtoP1P2_StoP3P4 , 2, 3, 0, 1),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 2, 3, 0, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_24, FF::BL2),
new Lineshapes::FOCUS("KPi10", Lineshapes::FOCUS::Mod::I32, KPi1_0_M, KPi1_0_W, 0, M_13, FF::BL2, 1.5),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_34, FF::BL2),
new Lineshapes::FOCUS("KPi10", Lineshapes::FOCUS::Mod::I32, KPi1_0_M, KPi1_0_W, 0, M_12, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{rho(770)0{pi+,pi-},KPi10[FOCUS.I32]{K-,pi+}}",
mkvar("D0{rho(770)0{pi+,pi-},KPi10[FOCUS.I32]{K-,pi+}}_r", false, 0.0985307, 0.00625139),
mkvar("D0{rho(770)0{pi+,pi-},KPi10[FOCUS.I32]{K-,pi+}}_i", false, 0.32325, 0.00715978),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
// Line 32
// D0{rho(770)0{pi+,pi-},KPi10[FOCUS.Kpi]{K-,pi+}}
spin_factor_list.push_back(std::vector<SpinFactor*>({
new SpinFactor("SF", SF_4Body::DtoVS_VtoP1P2_StoP3P4 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 1, 3, 0, 2),
new SpinFactor("SF", SF_4Body::DtoVS_VtoP1P2_StoP3P4 , 2, 3, 0, 1),
new SpinFactor("SF", SF_4Body::FF_12_34_L1 , 2, 3, 0, 1)
}));
line_factor_list.push_back(std::vector<Lineshape*>{
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_24, FF::BL2),
new Lineshapes::FOCUS("KPi10", Lineshapes::FOCUS::Mod::Kpi, KPi1_0_M, KPi1_0_W, 0, M_13, FF::BL2, 1.5),
new Lineshapes::RBW("rho(770)0", rho_770_0_M, rho_770_0_W, 1.0, M_34, FF::BL2),
new Lineshapes::FOCUS("KPi10", Lineshapes::FOCUS::Mod::Kpi, KPi1_0_M, KPi1_0_W, 0, M_12, FF::BL2, 1.5)
});
amplitudes_list.push_back(new Amplitude{
"D0{rho(770)0{pi+,pi-},KPi10[FOCUS.Kpi]{K-,pi+}}",
mkvar("D0{rho(770)0{pi+,pi-},KPi10[FOCUS.Kpi]{K-,pi+}}_r", false, 0.0985307, 0.00625139),
mkvar("D0{rho(770)0{pi+,pi-},KPi10[FOCUS.Kpi]{K-,pi+}}_i", false, 0.32325, 0.00715978),
line_factor_list.back(),
spin_factor_list.back(),
2});
DK3P_DI.amplitudes_B.push_back(amplitudes_list.back());
[ ]: