Bug Summary

File:generators/evtgen/models/besiii/src/EvtDToKSpipi0pi0.cc
Warning:line 456, column 9
Value stored to 'temp_PDF' is never read

Annotated Source Code

Press '?' to see keyboard shortcuts

clang -cc1 -cc1 -triple x86_64-unknown-linux-gnu -O3 -analyze -disable-free -clear-ast-before-backend -disable-llvm-verifier -discard-value-names -main-file-name EvtDToKSpipi0pi0.cc -analyzer-checker=core -analyzer-checker=apiModeling -analyzer-checker=unix -analyzer-checker=deadcode -analyzer-checker=cplusplus -analyzer-checker=security.insecureAPI.UncheckedReturn -analyzer-checker=security.insecureAPI.getpw -analyzer-checker=security.insecureAPI.gets -analyzer-checker=security.insecureAPI.mktemp -analyzer-checker=security.insecureAPI.mkstemp -analyzer-checker=security.insecureAPI.vfork -analyzer-checker=nullability.NullPassedToNonnull -analyzer-checker=nullability.NullReturnedFromNonnull -analyzer-output plist -w -setup-static-analyzer -mrelocation-model pic -pic-level 2 -fhalf-no-semantic-interposition -mframe-pointer=none -fmath-errno -ffp-contract=on -fno-rounding-math -mconstructor-aliases -funwind-tables=2 -target-cpu x86-64 -tune-cpu generic -debugger-tuning=gdb -fdebug-compilation-dir=/data/b2soft/buildbot/development/build -fcoverage-compilation-dir=/data/b2soft/buildbot/development/build -resource-dir /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/lib/clang/21 -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/c++ -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/c++/x86_64-redhat-linux -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/c++/backward -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/include -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/python3.12 -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/include/CLHEP -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/Geant4 -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/include/root -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/include/belle_legacy -I include/ -D _PACKAGE_="generators" -D G4UI_USE_TCSH -D RaveDllExport= -D HAS_SQLITE -D HAS_CALLGRIND -I include -I /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/libxml2 -internal-isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/bin/../lib64/gcc/x86_64-redhat-linux/15.2.0/../../../../include/c++ -internal-isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/bin/../lib64/gcc/x86_64-redhat-linux/15.2.0/../../../../include/c++/x86_64-redhat-linux -internal-isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/bin/../lib64/gcc/x86_64-redhat-linux/15.2.0/../../../../include/c++/backward -internal-isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/lib/clang/21/include -internal-isystem /usr/local/include -internal-isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/bin/../lib64/gcc/x86_64-redhat-linux/15.2.0/../../../../x86_64-redhat-linux/include -internal-externc-isystem /include -internal-externc-isystem /usr/include -Wno-missing-braces -Wno-unused-command-line-argument -std=c++20 -fdeprecated-macro -ferror-limit 19 -fgnuc-version=4.2.1 -fno-implicit-modules -fskip-odr-check-in-gmf -fcxx-exceptions -fexceptions -vectorize-loops -vectorize-slp -analyzer-output=html -faddrsig -D__GCC_HAVE_DWARF2_CFI_ASM=1 -o /scan_build/2026-05-31-004316-385593-1 -x c++ generators/evtgen/models/besiii/src/EvtDToKSpipi0pi0.cc
1// Model: EvtDToKSpipi0pi0
2// This file is an amplitude model for D+ -> K_S0 pi- pi0 pi0.
3// The model is from the BESIII Collaboration in JHEP09 (2023) 077. DOI:  https://doi.org/10.1007/JHEP09(2023)077
4//
5// Permission to include these files in basf2 was generously granted by the BESIII Collaboration.
6//
7// Please cite the original reference for any public/published results where this model was used.
8
9#include <EvtGenBase/EvtPatches.hh>
10#include <EvtGenBase/EvtParticle.hh>
11#include <EvtGenBase/EvtGenKine.hh>
12#include <EvtGenBase/EvtPDL.hh>
13#include <EvtGenBase/EvtReport.hh>
14#include <EvtGenBase/EvtComplex.hh>
15#include <EvtGenBase/EvtDecayTable.hh>
16#include <stdlib.h>
17#include <stdio.h>
18#include <iostream>
19#include <cmath>
20
21#include <generators/evtgen/EvtGenModelRegister.h>
22#include <generators/evtgen/models/besiii/EvtDToKSpipi0pi0.h>
23
24using namespace std;
25
26namespace Belle2 {
27
28 /** register the model in EvtGen */
29 B2_EVTGEN_REGISTER_MODEL(EvtDToKSpipi0pi0)namespace { Belle2::EvtGenModelRegister::Factory<EvtDToKSpipi0pi0
> EvtGenModelFactory_EvtDToKSpipi0pi0; };
30
31 EvtDToKSpipi0pi0::~EvtDToKSpipi0pi0() {}
32
33 std::string EvtDToKSpipi0pi0::getName()
34 {
35 return "DToKSpipi0pi0";
36 }
37
38 EvtDecayBase* EvtDToKSpipi0pi0::clone()
39 {
40 return new EvtDToKSpipi0pi0;
41 }
42
43 void EvtDToKSpipi0pi0::init()
44 {
45
46 checkNArg(0);
47 checkNDaug(4);
48 checkSpinParent(EvtSpinType::SCALAR);
49
50 mK1270 = 1.272; mK1400 = 1.403;
51 GK1270 = 0.09; GK1400 = 0.174;
52 mKstr0 = 0.89555; mrho = 0.77511;
53 GKstr0 = 0.0473; Grho = 0.1485;
54 msigma = 0.526; ma1 = 1.230;
55 Gsigma = 0.535; Ga1 = 0.420;
56 mD = 1.86965;
57 math_pi = 3.1415926;
58
59 rho[0] = 3.0177e+00;//5
60 phi[0] = 7.5497e-01;
61
62 rho[1] = 1;//a1(1260)
63 phi[1] = 0;
64
65 rho[2] = 2.3527e-01;//7S_1400
66 phi[2] = -2.9985e+00;
67
68 rho[3] = 5.5731e-01;//7D_1400
69 phi[3] = 4.2940e+00;
70
71 rho[4] = 9.0803e-01;//11_sigma
72 phi[4] = 4.7731e+00;
73
74 rho[5] = 2.5345e-01;//23S
75 phi[5] = -3.3205e+00;
76
77 rho[6] = 6.0504e-02;//23P
78 phi[6] = -1.6803e+00;
79
80 rho[7] = 7.6978e-01;//25S
81 phi[7] = -5.5937e+00;
82
83 modetype[0] = 403;
84 modetype[1] = 100;
85 modetype[2] = 200;
86 modetype[3] = 200;
87 modetype[4] = 204;
88 modetype[5] = 2;
89 modetype[6] = 2;
90 modetype[7] = 2;
91
92 int GG[4][4] = { {1, 0, 0, 0}, {0, -1, 0, 0}, {0, 0, -1, 0}, {0, 0, 0, -1} };
93 int EE[4][4][4][4] = {
94 { {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0} },
95 {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 1}, {0, 0, -1, 0}},
96 {{0, 0, 0, 0}, {0, 0, 0, -1}, {0, 0, 0, 0}, {0, 1, 0, 0} },
97 {{0, 0, 0, 0}, {0, 0, 1, 0}, {0, -1, 0, 0}, {0, 0, 0, 0} }
98 },
99 { {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, -1}, {0, 0, 1, 0} },
100 {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0} },
101 {{0, 0, 0, 1}, {0, 0, 0, 0}, {0, 0, 0, 0}, {-1, 0, 0, 0}},
102 {{0, 0, -1, 0}, {0, 0, 0, 0}, {1, 0, 0, 0}, {0, 0, 0, 0} }
103 },
104 { {{0, 0, 0, 0}, {0, 0, 0, 1}, {0, 0, 0, 0}, {0, -1, 0, 0}},
105 {{0, 0, 0, -1}, {0, 0, 0, 0}, {0, 0, 0, 0}, {1, 0, 0, 0} },
106 {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0} },
107 {{0, 1, 0, 0}, {-1, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0} }
108 },
109 { {{0, 0, 0, 0}, {0, 0, -1, 0}, {0, 1, 0, 0}, {0, 0, 0, 0} },
110 {{0, 0, 1, 0}, {0, 0, 0, 0}, {-1, 0, 0, 0}, {0, 0, 0, 0} },
111 {{0, -1, 0, 0}, {1, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0} },
112 {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0} }
113 }
114 };
115 for (int i = 0; i < 4; i++) {
116 for (int j = 0; j < 4; j++) {
117 G[i][j] = GG[i][j];
118 for (int k = 0; k < 4; k++) {
119 for (int l = 0; l < 4; l++) {
120 E[i][j][k][l] = EE[i][j][k][l];
121 }
122 }
123 }
124 }
125 }
126
127 void EvtDToKSpipi0pi0::initProbMax()
128 {
129 setProbMax(2325.0);
130 }
131
132 void EvtDToKSpipi0pi0::decay(EvtParticle* p)
133 {
134 //-----------for max value------------------
135 /* double maxprob = 0.0;
136 for(int ir=0;ir<=60000000;ir++){
137 p->initializePhaseSpace(getNDaug(),getDaugs());
138 EvtVector4R Ks0 = p->getDaug(0)->getP4();
139 EvtVector4R pi1 = p->getDaug(1)->getP4();
140 EvtVector4R pi2 = p->getDaug(2)->getP4();
141 EvtVector4R pi3 = p->getDaug(3)->getP4();
142 double value;
143 double Ks[4],Pip[4],Pi01[4],Pi02[4];
144 mother_c=EvtPDL::getStdHep(p->getId());
145 if(mother_c==411){
146 Ks[0] = Ks0.get(0); Pip[0] = pi1.get(0); Pi01[0] = pi2.get(0); Pi02[0] = pi3.get(0);
147 Ks[1] = Ks0.get(1); Pip[1] = pi1.get(1); Pi01[1] = pi2.get(1); Pi02[1] = pi3.get(1);
148 Ks[2] = Ks0.get(2); Pip[2] = pi1.get(2); Pi01[2] = pi2.get(2); Pi02[2] = pi3.get(2);
149 Ks[3] = Ks0.get(3); Pip[3] = pi1.get(3); Pi01[3] = pi2.get(3); Pi02[3] = pi3.get(3);
150 }else if(mother_c==-411){
151 Ks[0] = Ks0.get(0); Pip[0] = pi1.get(0); Pi01[0] = pi2.get(0); Pi02[0] = pi3.get(0);
152 Ks[1] = -1.0*Ks0.get(1); Pip[1] = -1.0*pi1.get(1); Pi01[1] = -1.0*pi2.get(1); Pi02[1] = -1.0*pi3.get(1);
153 Ks[2] = -1.0*Ks0.get(2); Pip[2] = -1.0*pi1.get(2); Pi01[2] = -1.0*pi2.get(2); Pi02[2] = -1.0*pi3.get(2);
154 Ks[3] = -1.0*Ks0.get(3); Pip[3] = -1.0*pi1.get(3); Pi01[3] = -1.0*pi2.get(3); Pi02[3] = -1.0*pi3.get(3);}
155 calPDF(Ks, Pip, Pi01, Pi02, value);
156 if(value>maxprob) {
157 maxprob=value;
158 }
159 }
160 return;*/
161 //-----------------------------------------------
162 p->initializePhaseSpace(getNDaug(), getDaugs());
163 EvtVector4R Ks0 = p->getDaug(0)->getP4();
164 EvtVector4R pip = p->getDaug(1)->getP4();
165 EvtVector4R pi01 = p->getDaug(2)->getP4();
166 EvtVector4R pi02 = p->getDaug(3)->getP4();
167
168 mother_c = EvtPDL::getStdHep(p->getId());
169 double Ks[4], Pip[4], Pi01[4], Pi02[4];
170 if (mother_c == 411) {
171 Ks[0] = Ks0.get(0); Pip[0] = pip.get(0); Pi01[0] = pi01.get(0); Pi02[0] = pi02.get(0);
172 Ks[1] = Ks0.get(1); Pip[1] = pip.get(1); Pi01[1] = pi01.get(1); Pi02[1] = pi02.get(1);
173 Ks[2] = Ks0.get(2); Pip[2] = pip.get(2); Pi01[2] = pi01.get(2); Pi02[2] = pi02.get(2);
174 Ks[3] = Ks0.get(3); Pip[3] = pip.get(3); Pi01[3] = pi01.get(3); Pi02[3] = pi02.get(3);
175 } else if (mother_c == -411) {
176 Ks[0] = Ks0.get(0); Pip[0] = pip.get(0); Pi01[0] = pi01.get(0); Pi02[0] = pi02.get(0);
177 Ks[1] = -1.0 * Ks0.get(1); Pip[1] = -1.0 * pip.get(1); Pi01[1] = -1.0 * pi01.get(1); Pi02[1] = -1.0 * pi02.get(1);
178 Ks[2] = -1.0 * Ks0.get(2); Pip[2] = -1.0 * pip.get(2); Pi01[2] = -1.0 * pi01.get(2); Pi02[2] = -1.0 * pi02.get(2);
179 Ks[3] = -1.0 * Ks0.get(3); Pip[3] = -1.0 * pip.get(3); Pi01[3] = -1.0 * pi01.get(3); Pi02[3] = -1.0 * pi02.get(3);
180 }
181
182 //Ks[0] = 0.656878382777544; Pip[0] = 0.465468253141211; Pi01[0] = 0.364737466814715; Pi02[0] = 0.399546086976042;
183 //Ks[1] = -0.426050514556759; Pip[1] = 0.440786345532845; Pi01[1] = -0.195437195063332; Pi02[1] = 0.352799494203111;
184 //Ks[2] = -0.015331547333183; Pip[2] = -0.003143271463748; Pi01[2] = 0.195716648451034; Pi02[2] = -0.106603408016079;
185 //Ks[3] = -0.046026135420350; Pip[3] = -0.053650975933667; Pi01[3] = 0.195739708428659; Pi02[3] = 0.074743717862502;
186
187 double value;
188 calPDF(Ks, Pip, Pi01, Pi02, value);
189 setProb(value);
190 return;
191 }
192
193 double EvtDToKSpipi0pi0::calPDF(const double Ks[], const double Pip[], const double Pi01[], const double Pi02[], double& Result)
194 {
195 double cof[2], amp_tmp1[2], amp_tmp2[2], amp_tmp[2], amp_PDF[2], PDF[2];
196 double flag[3], mass_R[2], width_R[2];
197 double sa[3], sb[3], sc[3], B[3] = {};
198 double t1D[4], t1V[4], t1A[4];
199 double pS[4], pV[4], pD[4], pA[4];
200 double pro1[2], pro2[2], proKPi_S[2], pro[2];
201
202 double rD2 = 25.0;
203 double rRes2 = 9.0;
204 double mass1[8] = {mrho, mrho, mKstr0, mKstr0, msigma, mKstr0, mKstr0, mKstr0 };
205 double mass2[8] = {mrho, ma1, mK1400, mK1400, ma1, mrho, mrho, mrho };
206 double width1[8] = {Grho, Grho, GKstr0, GKstr0, Gsigma, GKstr0, GKstr0, GKstr0 };
207 double width2[8] = {Grho, Ga1, GK1400, GK1400, Ga1, Grho, Grho, Grho };
208 double g0[8] = { 0, 1, 1, 1, 1, 1, 1, 1 };
209 double g1[8] = { 0, 1, 1, 1, 1, 1, 1, 0 };
210 double g2[8] = { 0, 0, 0, 2, 1, 0, 1, 0 };
211 double temp_PDF = 0;
212 PDF[0] = 0;
213 PDF[1] = 0;
214
215 for (int i = 0; i < 8; i++) { //0719 i=0, 7
216 flag[0] = g0[i]; flag[1] = g1[i]; flag[2] = g2[i];
217 mass_R[0] = mass1[i]; mass_R[1] = mass2[i];
218 width_R[0] = width1[i]; width_R[1] = width2[i];
219
220 amp_tmp[0] = 0;
221 amp_tmp[1] = 0;
222 amp_tmp1[0] = 0;
223 amp_tmp1[1] = 0;
224 amp_tmp2[0] = 0;
225 amp_tmp2[1] = 0;
226
227 cof[0] = rho[i] * cos(phi[i]);
228 cof[1] = rho[i] * sin(phi[i]);
229
230 if (modetype[i] == 100) { //D->Ks (rho+ pi0 ),D->Ks a1(1260)+
231
232 temp_PDF = 0;
233 double t2A[4][4];
234
235 for (int ii = 0; ii != 4; ii++) {
236 pV[ii] = Pip[ii] + Pi02[ii];
237 pA[ii] = pV[ii] + Pi01[ii];
238 pD[ii] = pA[ii] + Ks[ii];
239 }
240 sa[0] = SCADot(pV, pV);
241 sb[0] = SCADot(Pip, Pip);
242 sc[0] = SCADot(Pi02, Pi02);
243 sa[1] = SCADot(pA, pA);
244 sb[1] = sa[0];
245 sc[1] = SCADot(Pi01, Pi01);
246 sa[2] = SCADot(pD, pD);
247 sb[2] = sa[1];
248 sc[2] = SCADot(Ks, Ks);
249 if (flag[0] == 1) {
250 propagatorGS(mass_R[0]*mass_R[0], mass_R[0], width_R[0], sa[0], sb[0], sc[0], rRes2, pro1);
251 } else if (flag[0] == 0) {pro1[0] = 1; pro1[1] = 0;}
252
253 if (flag[1] == 1) {
254 propagatorRBW_a1(mass_R[1]*mass_R[1], mass_R[1], width_R[1], sa[1], sb[1], sc[1], rRes2, flag[2], pro2);
255 } else if (flag[1] == 0) {pro2[0] = 1; pro2[1] = 0;}
256 B[0] = Barrier(mass_R[0] * mass_R[0], 1, sa[0], sb[0], sc[0], rRes2);
257 B[2] = Barrier(mD * mD, 1, sa[2], sb[2], sc[2], rD2);
258 calt1(Pip, Pi02, t1V);
259 calt1(pA, Ks, t1D);
260 if (flag[2] == 0) {
261 for (int ii = 0; ii != 4; ii++) {
262 for (int j = 0; j != 4; j++) {
263 temp_PDF += t1D[ii] * (G[ii][j] - pA[ii] * pA[j] / sa[1]) * t1V[j] * G[ii][ii] * G[j][j];
264 }
265 }
266 B[1] = 1;
267 } else if (flag[2] == 2) {
268 calt2(pV, Pi01, t2A);
269 for (int ii = 0; ii != 4; ii++) {
270 for (int j = 0; j != 4; j++) {
271 temp_PDF += t1D[ii] * t2A[ii][j] * t1V[j] * G[ii][ii] * G[j][j];
272 }
273 }
274 B[1] = Barrier(mass_R[1] * mass_R[1], 2, sa[1], sb[1], sc[1], rRes2);
275 }
276 Com_Multi(pro1, pro2, pro);
277 amp_tmp1[0] = temp_PDF * B[0] * B[1] * B[2] * pro[0];
278 amp_tmp1[1] = temp_PDF * B[0] * B[1] * B[2] * pro[1];
279
280
281 temp_PDF = 0;
282 for (int ii = 0; ii != 4; ii++) {
283 pV[ii] = Pip[ii] + Pi01[ii];
284 pA[ii] = pV[ii] + Pi02[ii];
285 pD[ii] = pA[ii] + Ks[ii];
286 }
287 sa[0] = SCADot(pV, pV);
288 sb[0] = SCADot(Pip, Pip);
289 sc[0] = SCADot(Pi01, Pi01);
290 sa[1] = SCADot(pA, pA);
291 sb[1] = sa[0];
292 sc[1] = SCADot(Pi02, Pi02);
293 sa[2] = SCADot(pD, pD);
294 sb[2] = sa[1];
295 sc[2] = SCADot(Ks, Ks);
296 if (flag[0] == 1) {
297 propagatorGS(mass_R[0]*mass_R[0], mass_R[0], width_R[0], sa[0], sb[0], sc[0], rRes2, pro1);
298 } else if (flag[0] == 0) {pro1[0] = 1; pro1[1] = 0;}
299
300 if (flag[1] == 1) {
301 propagatorRBW_a1(mass_R[1]*mass_R[1], mass_R[1], width_R[1], sa[1], sb[1], sc[1], rRes2, flag[2], pro2);
302 } else if (flag[1] == 0) {pro2[0] = 1; pro2[1] = 0;}
303 B[0] = Barrier(mass_R[0] * mass_R[0], 1, sa[0], sb[0], sc[0], rRes2);
304 B[2] = Barrier(mD * mD, 1, sa[2], sb[2], sc[2], rD2);
305 calt1(Pip, Pi01, t1V);
306 calt1(pA, Ks, t1D);
307 if (flag[2] == 0) {
308 for (int ii = 0; ii != 4; ii++) {
309 for (int j = 0; j != 4; j++) {
310 temp_PDF += t1D[ii] * (G[ii][j] - pA[ii] * pA[j] / sa[1]) * t1V[j] * G[ii][ii] * G[j][j];
311 }
312 }
313 B[1] = 1;
314 } else if (flag[2] == 2) {
315 calt2(pV, Pi02, t2A);
316 for (int ii = 0; ii != 4; ii++) {
317 for (int j = 0; j != 4; j++) {
318 temp_PDF += t1D[ii] * t2A[ii][j] * t1V[j] * G[ii][ii] * G[j][j];
319 }
320 }
321 B[1] = Barrier(mass_R[1] * mass_R[1], 2, sa[1], sb[1], sc[1], rRes2);
322 }
323 Com_Multi(pro1, pro2, pro);
324 amp_tmp2[0] = temp_PDF * B[0] * B[1] * B[2] * pro[0];
325 amp_tmp2[1] = temp_PDF * B[0] * B[1] * B[2] * pro[1];
326 } else if (modetype[i] == 200) { //D->(K*0 pi0) pi+
327
328 temp_PDF = 0;
329 double t2A[4][4];
330 for (int ii = 0; ii != 4; ii++) {
331 pV[ii] = Ks[ii] + Pi02[ii];
332 pA[ii] = pV[ii] + Pi01[ii];
333 pD[ii] = pA[ii] + Pip[ii];
334 }
335 sa[0] = SCADot(pV, pV);
336 sb[0] = SCADot(Ks, Ks);
337 sc[0] = SCADot(Pi02, Pi02);
338 sa[1] = SCADot(pA, pA);
339 sb[1] = sa[0];
340 sc[1] = SCADot(Pi01, Pi01);
341 sa[2] = SCADot(pD, pD);
342 sb[2] = sa[1];
343 sc[2] = SCADot(Pip, Pip);
344 if (flag[0] == 1) {
345 propagatorRBW(mass_R[0]*mass_R[0], mass_R[0], width_R[0], sa[0], sb[0], sc[0], rRes2, 1, pro1);
346 } else if (flag[0] == 0) {pro1[0] = 1; pro1[1] = 0;}
347
348 if (flag[1] == 1) {
349 propagatorRBW(mass_R[1]*mass_R[1], mass_R[1], width_R[1], sa[1], sb[1], sc[1], rRes2, flag[2], pro2);
350 } else if (flag[1] == 0) {pro2[0] = 1; pro2[1] = 0;}
351 B[0] = Barrier(mass_R[0] * mass_R[0], 1, sa[0], sb[0], sc[0], rRes2);
352 B[2] = Barrier(mD * mD, 1, sa[2], sb[2], sc[2], rD2);
353 calt1(Ks, Pi02, t1V);
354 calt1(pA, Pip, t1D);
355 if (flag[2] == 0) {
356 for (int ii = 0; ii != 4; ii++) {
357 for (int j = 0; j != 4; j++) {
358 temp_PDF += t1D[ii] * (G[ii][j] - pA[ii] * pA[j] / sa[1]) * t1V[j] * G[ii][ii] * G[j][j];
359 }
360 }
361 B[1] = 1;
362 } else if (flag[2] == 2) {
363 calt2(pV, Pi01, t2A);
364 for (int ii = 0; ii != 4; ii++) {
365 for (int j = 0; j != 4; j++) {
366 temp_PDF += t1D[ii] * t2A[ii][j] * t1V[j] * G[ii][ii] * G[j][j];
367 }
368 }
369 B[1] = Barrier(mass_R[1] * mass_R[1], 2, sa[1], sb[1], sc[1], rRes2);
370 }
371 Com_Multi(pro1, pro2, pro);
372 amp_tmp1[0] = temp_PDF * B[0] * B[1] * B[2] * pro[0];
373 amp_tmp1[1] = temp_PDF * B[0] * B[1] * B[2] * pro[1];
374 temp_PDF = 0;
375 for (int ii = 0; ii != 4; ii++) {
376 pV[ii] = Ks[ii] + Pi01[ii];
377 pA[ii] = pV[ii] + Pi02[ii];
378 pD[ii] = pA[ii] + Pip[ii];
379 }
380 sa[0] = SCADot(pV, pV);
381 sb[0] = SCADot(Ks, Ks);
382 sc[0] = SCADot(Pi01, Pi01);
383 sa[1] = SCADot(pA, pA);
384 sb[1] = sa[0];
385 sc[1] = SCADot(Pi02, Pi02);
386 sa[2] = SCADot(pD, pD);
387 sb[2] = sa[1];
388 sc[2] = SCADot(Pip, Pip);
389 if (flag[0] == 1) {
390 propagatorRBW(mass_R[0]*mass_R[0], mass_R[0], width_R[0], sa[0], sb[0], sc[0], rRes2, 1, pro1);
391 } else if (flag[0] == 0) {pro1[0] = 1; pro1[1] = 0;}
392
393 if (flag[1] == 1) {
394 propagatorRBW(mass_R[1]*mass_R[1], mass_R[1], width_R[1], sa[1], sb[1], sc[1], rRes2, flag[2], pro2);
395 } else if (flag[1] == 0) {pro2[0] = 1; pro2[1] = 0;}
396
397 B[0] = Barrier(mass_R[0] * mass_R[0], 1, sa[0], sb[0], sc[0], rRes2);
398 B[2] = Barrier(mD * mD, 1, sa[2], sb[2], sc[2], rD2);
399 calt1(Ks, Pi01, t1V);
400 calt1(pA, Pip, t1D);
401 if (flag[2] == 0) {
402 for (int ii = 0; ii != 4; ii++) {
403 for (int j = 0; j != 4; j++) {
404 temp_PDF += t1D[ii] * (G[ii][j] - pA[ii] * pA[j] / sa[1]) * t1V[j] * G[ii][ii] * G[j][j];
405 }
406 }
407 B[1] = 1;
408 } else if (flag[2] == 2) {
409 calt2(pV, Pi02, t2A);
410 for (int ii = 0; ii != 4; ii++) {
411 for (int j = 0; j != 4; j++) {
412 temp_PDF += t1D[ii] * t2A[ii][j] * t1V[j] * G[ii][ii] * G[j][j];
413 }
414 }
415 B[1] = Barrier(mass_R[1] * mass_R[1], 2, sa[1], sb[1], sc[1], rRes2);
416 }
417
418 Com_Multi(pro1, pro2, pro);
419 amp_tmp2[0] = temp_PDF * B[0] * B[1] * B[2] * pro[0];
420 amp_tmp2[1] = temp_PDF * B[0] * B[1] * B[2] * pro[1];
421 } else if (modetype[i] == 204) { //D->Ks a1(1260)+, a1(1260)+->sigma pi+
422 temp_PDF = 0;
423
424 for (int ii = 0; ii != 4; ii++) {
425 pS[ii] = Pi01[ii] + Pi02[ii];
426 pA[ii] = pS[ii] + Pip[ii];
427 pD[ii] = pA[ii] + Ks[ii];
428 }
429 sa[0] = SCADot(pS, pS);
430 sb[0] = SCADot(Pi01, Pi01);
431 sc[0] = SCADot(Pi02, Pi02);
432 sa[1] = SCADot(pA, pA);
433 sb[1] = sa[0];
434 sc[1] = SCADot(Pip, Pip);
435 sa[2] = SCADot(pD, pD);
436 sb[2] = sa[1];
437 sc[2] = SCADot(Ks, Ks);
438 B[1] = Barrier(mass_R[1] * mass_R[1], 1, sa[1], sb[1], sc[1], rRes2);
439 B[2] = Barrier(mD * mD, 1, sa[2], sb[2], sc[2], rD2);
440 calt1(pA, Ks, t1D);
441 calt1(pS, Pip, t1A);
442 for (int ii = 0; ii != 4; ii++) {
443 temp_PDF += t1D[ii] * t1A[ii] * (G[ii][ii]);
444 }
445 if (flag[0] == 1) { //sigma500
446 propagatorsigma500(sa[0], sb[0], sc[0], pro1);
447 } else if (flag[0] == 0) {pro1[0] = 1; pro1[1] = 0;}
448
449 if (flag[1] == 1) { //a1
450 propagatorRBW_a1(mass_R[1]*mass_R[1], mass_R[1], width_R[1], sa[1], sb[1], sc[1], rRes2, 1, pro2);
451 } else if (flag[1] == 0) {pro2[0] = 1; pro2[1] = 0;}
452
453 Com_Multi(pro1, pro2, pro);
454 amp_tmp1[0] = temp_PDF * B[1] * B[2] * pro[0];
455 amp_tmp1[1] = temp_PDF * B[1] * B[2] * pro[1];
456 temp_PDF = 0;
Value stored to 'temp_PDF' is never read
457 amp_tmp2[0] = amp_tmp1[0];
458 amp_tmp2[1] = amp_tmp1[1];
459 } else if (modetype[i] == 120) { //D->(Ks pi0)S rho+
460 temp_PDF = 0;
461
462 for (int ii = 0; ii != 4; ii++) {
463 pS[ii] = Ks[ii] + Pi02[ii];
464 pV[ii] = Pip[ii] + Pi01[ii];
465 pD[ii] = pS[ii] + pV[ii];
466 }
467 sa[0] = SCADot(pS, pS);
468 sb[0] = SCADot(Ks, Ks);
469 sc[0] = SCADot(Pi02, Pi02);
470 sa[1] = SCADot(pV, pV);
471 sb[1] = SCADot(Pip, Pip);
472 sc[1] = SCADot(Pi01, Pi01);
473 sa[2] = SCADot(pD, pD);
474 sb[2] = sa[0];
475 sc[2] = sa[1];
476 if (flag[0] == 1) {
477 propagatorGS(mass_R[1]*mass_R[1], mass_R[1], width_R[0], sa[1], sb[1], sc[1], rRes2, pro1);
478 } else if (flag[0] == 0) { pro1[0] = 1; pro1[1] = 0; }
479 KPiSLASS(sa[0], sb[0], sc[0], proKPi_S);
480 if (flag[1] == 1) {
481 pro2[0] = proKPi_S[0]; pro2[1] = proKPi_S[1];
482 } else if (flag[1] == 0) {pro2[0] = 1; pro2[1] = 0;}
483
484 B[1] = Barrier(mass_R[1] * mass_R[1], 1, sa[1], sb[1], sc[1], rRes2);
485 B[2] = Barrier(mD * mD, 1, sa[2], sb[2], sc[2], rD2);
486 calt1(Pip, Pi01, t1V);
487 calt1(pS, pV, t1D);
488 for (int ii = 0; ii != 4; ii++) {
489 temp_PDF += G[ii][ii] * t1D[ii] * t1V[ii];
490 }
491 Com_Multi(pro1, pro2, pro);
492 amp_tmp1[0] = temp_PDF * B[1] * B[2] * pro[0];
493 amp_tmp1[1] = temp_PDF * B[1] * B[2] * pro[1];
494
495 temp_PDF = 0;
496 for (int ii = 0; ii != 4; ii++) {
497 pS[ii] = Ks[ii] + Pi01[ii];
498 pV[ii] = Pip[ii] + Pi02[ii];
499 pD[ii] = pS[ii] + pV[ii];
500 }
501 sa[0] = SCADot(pS, pS);
502 sb[0] = SCADot(Ks, Ks);
503 sc[0] = SCADot(Pi01, Pi01);
504 sa[1] = SCADot(pV, pV);
505 sb[1] = SCADot(Pip, Pip);
506 sc[1] = SCADot(Pi02, Pi02);
507 sa[2] = SCADot(pD, pD);
508 sb[2] = sa[0];
509 sc[2] = sa[1];
510 if (flag[0] == 1) {
511 propagatorGS(mass_R[1]*mass_R[1], mass_R[1], width_R[0], sa[1], sb[1], sc[1], rRes2, pro1);
512 } else if (flag[0] == 0) { pro1[0] = 1; pro1[1] = 0; }
513 KPiSLASS(sa[0], sb[0], sc[0], proKPi_S);
514 if (flag[1] == 1) {
515 pro2[0] = proKPi_S[0]; pro2[1] = proKPi_S[1];
516 } else if (flag[1] == 0) {pro2[0] = 1; pro2[1] = 0;}
517 B[1] = Barrier(mass_R[1] * mass_R[1], 1, sa[1], sb[1], sc[1], rRes2);
518 B[2] = Barrier(mD * mD, 1, sa[2], sb[2], sc[2], rD2);
519 calt1(Pip, Pi02, t1V);
520 calt1(pS, pV, t1D);
521 for (int ii = 0; ii != 4; ii++) {
522 temp_PDF += G[ii][ii] * t1D[ii] * t1V[ii];
523 }
524 Com_Multi(pro1, pro2, pro);
525 amp_tmp2[0] = temp_PDF * B[1] * B[2] * pro[0];
526 amp_tmp2[1] = temp_PDF * B[1] * B[2] * pro[1];
527 } else if (modetype[i] == 2) {
528 double t1V1[4], t1V2[4], t2D[4][4];
529 temp_PDF = 0;
530
531 double pV1[4], pV2[4];
532 for (int ii = 0; ii != 4; ii++) {
533 pV1[ii] = Ks[ii] + Pi01[ii];
534 pV2[ii] = Pip[ii] + Pi02[ii];
535 pD[ii] = pV1[ii] + pV2[ii];
536 }
537 sa[0] = SCADot(pV1, pV1);
538 sb[0] = SCADot(Ks, Ks);
539 sc[0] = SCADot(Pi01, Pi01);
540 sa[1] = SCADot(pV2, pV2);
541 sb[1] = SCADot(Pip, Pip);
542 sc[1] = SCADot(Pi02, Pi02);
543 sa[2] = SCADot(pD, pD);
544 sb[2] = sa[0];
545 sc[2] = sa[1];
546 if (flag[0] == 1) {
547 propagatorRBW(mass_R[0]*mass_R[0], mass_R[0], width_R[0], sa[0], sb[0], sc[0], rRes2, 1, pro1);
548 } else if (flag[0] == 0) {pro1[0] = 1; pro1[1] = 0;}
549 if (flag[1] == 1) {
550 propagatorGS(mass_R[1]*mass_R[1], mass_R[1], width_R[1], sa[1], sb[1], sc[1], rRes2, pro2);
551 } else if (flag[1] == 0) {pro2[0] = 1; pro2[1] = 0;}
552 B[0] = Barrier(mass_R[0] * mass_R[0], 1, sa[0], sb[0], sc[0], rRes2);
553 B[1] = Barrier(mass_R[1] * mass_R[1], 1, sa[1], sb[1], sc[1], rRes2);
554 calt1(Ks, Pi01, t1V1);
555 calt1(Pip, Pi02, t1V2);
556 if (flag[2] == 0) {
557 for (int ii = 0; ii != 4; ii++) {
558 temp_PDF += (G[ii][ii]) * t1V1[ii] * t1V2[ii];
559 }
560 B[2] = 1;
561 }
562 if (flag[2] == 1) {
563 calt1(pV1, pV2, t1D);
564 for (int ii = 0; ii != 4; ii++) {
565 for (int j = 0; j != 4; j++) {
566 for (int k = 0; k != 4; k++) {
567 for (int l = 0; l != 4; l++) {
568 temp_PDF += E[ii][j][k][l] * pD[ii] * t1D[j] * t1V1[k] * t1V2[l] *
569 (G[ii][ii]) * (G[j][j]) * (G[l][l]) * (G[k][k]);
570 }
571 }
572 }
573 }
574 B[2] = Barrier(mD * mD, 1, sa[2], sb[2], sc[2], rD2);
575 }
576 if (flag[2] == 2) {
577 calt2(pV1, pV2, t2D);
578 for (int ii = 0; ii != 4; ii++) {
579 for (int j = 0; j != 4; j++) {
580 temp_PDF += t2D[ii][j] * t1V1[ii] * t1V2[j] * (G[ii][ii]) * (G[j][j]);
581 }
582 }
583 B[2] = Barrier(mD * mD, 2, sa[2], sb[2], sc[2], rD2);
584 }
585 Com_Multi(pro1, pro2, pro);
586 amp_tmp1[0] = temp_PDF * B[0] * B[1] * B[2] * pro[0];
587 amp_tmp1[1] = temp_PDF * B[0] * B[1] * B[2] * pro[1];
588 temp_PDF = 0;
589 for (int ii = 0; ii != 4; ii++) {
590 pV1[ii] = Ks[ii] + Pi02[ii];
591 pV2[ii] = Pip[ii] + Pi01[ii];
592 pD[ii] = pV1[ii] + pV2[ii];
593 }
594 sa[0] = SCADot(pV1, pV1);
595 sb[0] = SCADot(Ks, Ks);
596 sc[0] = SCADot(Pi02, Pi02);
597 sa[1] = SCADot(pV2, pV2);
598 sb[1] = SCADot(Pip, Pip);
599 sc[1] = SCADot(Pi01, Pi01);
600 sa[2] = SCADot(pD, pD);
601 sb[2] = sa[0];
602 sc[2] = sa[1];
603 if (flag[0] == 1) {
604 propagatorRBW(mass_R[0]*mass_R[0], mass_R[0], width_R[0], sa[0], sb[0], sc[0], rRes2, 1, pro1);
605 } else if (flag[0] == 0) {pro1[0] = 1; pro1[1] = 0;}
606 if (flag[1] == 1) {
607 propagatorGS(mass_R[1]*mass_R[1], mass_R[1], width_R[1], sa[1], sb[1], sc[1], rRes2, pro2);
608 } else if (flag[1] == 0) {pro2[0] = 1; pro2[1] = 0;}
609 B[0] = Barrier(mass_R[0] * mass_R[0], 1, sa[0], sb[0], sc[0], rRes2);
610 B[1] = Barrier(mass_R[1] * mass_R[1], 1, sa[1], sb[1], sc[1], rRes2);
611 calt1(Ks, Pi02, t1V1);
612 calt1(Pip, Pi01, t1V2);
613 if (flag[2] == 0) {
614 for (int ii = 0; ii != 4; ii++) {
615 temp_PDF += (G[ii][ii]) * t1V1[ii] * t1V2[ii];
616 }
617 B[2] = 1;
618 }
619 if (flag[2] == 1) {
620 calt1(pV1, pV2, t1D);
621 for (int ii = 0; ii != 4; ii++) {
622 for (int j = 0; j != 4; j++) {
623 for (int k = 0; k != 4; k++) {
624 for (int l = 0; l != 4; l++) {
625 temp_PDF += E[ii][j][k][l] * pD[ii] * t1D[j] * t1V1[k] * t1V2[l] *
626 (G[ii][ii]) * (G[j][j]) * (G[l][l]) * (G[k][k]);
627 }
628 }
629 }
630 }
631 B[2] = Barrier(mD * mD, 1, sa[2], sb[2], sc[2], rD2);
632 }
633 if (flag[2] == 2) {
634 calt2(pV1, pV2, t2D);
635 for (int ii = 0; ii != 4; ii++) {
636 for (int j = 0; j != 4; j++) {
637 temp_PDF += t2D[ii][j] * t1V1[ii] * t1V2[j] * (G[ii][ii]) * (G[j][j]);
638 }
639 }
640 B[2] = Barrier(mD * mD, 2, sa[2], sb[2], sc[2], rD2);
641 }
642 Com_Multi(pro1, pro2, pro);
643 amp_tmp2[0] = temp_PDF * B[0] * B[1] * B[2] * pro[0];
644 amp_tmp2[1] = temp_PDF * B[0] * B[1] * B[2] * pro[1];
645 } else if (modetype[i] == 43) {
646 double KPi[4];
647 for (int ii = 0; ii != 4; ii++) {
648 KPi[ii] = Ks[ii] + Pi01[ii];
649 }
650 sa[0] = SCADot(KPi, KPi);
651 sb[0] = SCADot(Ks, Ks);
652 sc[0] = SCADot(Pi01, Pi01);
653 KPiSLASS(sa[0], sb[0], sc[0], proKPi_S);
654 if (flag[0] == 1) {
655 pro1[0] = proKPi_S[0]; pro1[1] = proKPi_S[1];
656 } else if (flag[0] == 0) {pro1[0] = 1; pro1[1] = 0;}
657
658 amp_tmp1[0] = pro1[0];
659 amp_tmp1[1] = pro1[1];
660
661 for (int ii = 0; ii != 4; ii++) {
662 KPi[ii] = Ks[ii] + Pi02[ii];
663 }
664 sa[0] = SCADot(KPi, KPi);
665 sb[0] = SCADot(Ks, Ks);
666 sc[0] = SCADot(Pi02, Pi02);
667 KPiSLASS(sa[0], sb[0], sc[0], proKPi_S);
668 if (flag[0] == 1) {
669 pro1[0] = proKPi_S[0]; pro1[1] = proKPi_S[1];
670 } else if (flag[0] == 0) {pro1[0] = 1; pro1[1] = 0;}
671
672 amp_tmp2[0] = pro1[0];
673 amp_tmp2[1] = pro1[1];
674 } else if (modetype[i] == 403) {
675 temp_PDF = 0;
676 double pP[4];
677 for (int ii = 0; ii != 4; ii++) {
678 pV[ii] = Pip[ii] + Pi02[ii];
679 pP[ii] = pV[ii] + Pi01[ii];
680 pD[ii] = pP[ii] + Ks[ii];
681 }
682 sa[0] = SCADot(pV, pV);
683 sb[0] = SCADot(Pip, Pip);
684 sc[0] = SCADot(Pi02, Pi02);
685 sa[1] = SCADot(pP, pP);
686 sb[1] = sa[0];
687 sc[1] = SCADot(Pi01, Pi01);
688 sa[2] = SCADot(pD, pD);
689 sb[2] = sa[1];
690 sc[2] = SCADot(Ks, Ks);
691 B[0] = Barrier(mass_R[0] * mass_R[0], 1, sa[0], sb[0], sc[0], rRes2);
692 B[1] = Barrier(mass_R[1] * mass_R[1], 1, sa[1], sb[1], sc[1], rRes2);
693 propagatorGS(mass_R[0]*mass_R[0], mass_R[0], width_R[0], sa[0], sb[0], sc[0], rRes2, pro);
694 calt1(Pip, Pi02, t1V);
695 for (int ii = 0; ii != 4; ii++) {
696 temp_PDF += Pi01[ii] * t1V[ii] * (G[ii][ii]);
697 }
698 amp_tmp1[0] = temp_PDF * B[0] * B[1] * pro[0];
699 amp_tmp1[1] = temp_PDF * B[0] * B[1] * pro[1];
700 temp_PDF = 0;
701 for (int ii = 0; ii != 4; ii++) {
702 pV[ii] = Pip[ii] + Pi01[ii];
703 pP[ii] = pV[ii] + Pi02[ii];
704 pD[ii] = pP[ii] + Ks[ii];
705 }
706 sa[0] = SCADot(pV, pV);
707 sb[0] = SCADot(Pip, Pip);
708 sc[0] = SCADot(Pi01, Pi01);
709 sa[1] = SCADot(pP, pP);
710 sb[1] = sa[0];
711 sc[1] = SCADot(Pi02, Pi02);
712 sa[2] = SCADot(pD, pD);
713 sb[2] = sa[1];
714 sc[2] = SCADot(Ks, Ks);
715 B[0] = Barrier(mass_R[0] * mass_R[0], 1, sa[0], sb[0], sc[0], rRes2);
716 B[1] = Barrier(mass_R[1] * mass_R[1], 1, sa[1], sb[1], sc[1], rRes2);
717 propagatorGS(mass_R[0]*mass_R[0], mass_R[0], width_R[0], sa[0], sb[0], sc[0], rRes2, pro);
718 calt1(Pip, Pi01, t1V);
719 for (int ii = 0; ii != 4; ii++) {
720 temp_PDF += Pi02[ii] * t1V[ii] * (G[ii][ii]);
721 }
722 amp_tmp2[0] = temp_PDF * B[0] * B[1] * pro[0];
723 amp_tmp2[1] = temp_PDF * B[0] * B[1] * pro[1];
724 } else if (modetype[i] == 220) {
725 temp_PDF = 0;
726 for (int ii = 0; ii != 4; ii++) {
727 pS[ii] = Pip[ii] + Pi02[ii];
728 pV[ii] = Ks[ii] + Pi01[ii];
729 pD[ii] = pS[ii] + pV[ii];
730 }
731 sa[0] = SCADot(pS, pS);
732 sb[0] = SCADot(Pip, Pip);
733 sc[0] = SCADot(Pi02, Pi02);
734 sa[1] = SCADot(pV, pV);
735 sb[1] = SCADot(Ks, Ks);
736 sc[1] = SCADot(Pi01, Pi01);
737 sa[2] = SCADot(pD, pD);
738 sb[2] = sa[0];
739 sc[2] = sa[1];
740 if (flag[0] == 1) {
741 propagatorRBW(mass_R[0]*mass_R[0], mass_R[0], width_R[0], sa[1], sb[1], sc[1], rRes2, 1, pro);
742 }
743 if (flag[0] == 0) { pro[0] = 1; pro[1] = 0; }
744 B[1] = Barrier(mass_R[0] * mass_R[0], 1, sa[1], sb[1], sc[1], rRes2);
745 B[2] = Barrier(mD * mD, 1, sa[2], sb[2], sc[2], rD2);
746 calt1(Ks, Pi01, t1V);
747 calt1(pS, pV, t1D);
748 for (int ii = 0; ii != 4; ii++) {
749 temp_PDF += G[ii][ii] * t1D[ii] * t1V[ii];
750 }
751 amp_tmp1[0] = temp_PDF * B[1] * B[2] * pro[0];
752 amp_tmp1[1] = temp_PDF * B[1] * B[2] * pro[1];
753 temp_PDF = 0;
754 for (int ii = 0; ii != 4; ii++) {
755 pS[ii] = Pip[ii] + Pi01[ii];
756 pV[ii] = Ks[ii] + Pi02[ii];
757 pD[ii] = pS[ii] + pV[ii];
758 }
759 sa[0] = SCADot(pS, pS);
760 sb[0] = SCADot(Pip, Pip);
761 sc[0] = SCADot(Pi01, Pi01);
762 sa[1] = SCADot(pV, pV);
763 sb[1] = SCADot(Ks, Ks);
764 sc[1] = SCADot(Pi02, Pi02);
765 sa[2] = SCADot(pD, pD);
766 sb[2] = sa[0];
767 sc[2] = sa[1];
768 if (flag[0] == 1) {
769 propagatorRBW(mass_R[0]*mass_R[0], mass_R[0], width_R[0], sa[1], sb[1], sc[1], rRes2, 1, pro);
770 }
771 if (flag[0] == 0) { pro[0] = 1; pro[1] = 0; }
772 B[1] = Barrier(mass_R[0] * mass_R[0], 1, sa[1], sb[1], sc[1], rRes2);
773 B[2] = Barrier(mD * mD, 1, sa[2], sb[2], sc[2], rD2);
774 calt1(Ks, Pi02, t1V);
775 calt1(pS, pV, t1D);
776 for (int ii = 0; ii != 4; ii++) {
777 temp_PDF += G[ii][ii] * t1D[ii] * t1V[ii];
778 }
779 amp_tmp2[0] = temp_PDF * B[1] * B[2] * pro[0];
780 amp_tmp2[1] = temp_PDF * B[1] * B[2] * pro[1];
781 }
782 amp_tmp[0] = amp_tmp1[0] + amp_tmp2[0];
783 amp_tmp[1] = amp_tmp1[1] + amp_tmp2[1];
784 Com_Multi(amp_tmp, cof, amp_PDF);
785
786 PDF[0] += amp_PDF[0];
787 PDF[1] += amp_PDF[1];
788 }
789 double value = PDF[0] * PDF[0] + PDF[1] * PDF[1];
790
791 Result = value;
792
793 return value;
794 }
795
796 void EvtDToKSpipi0pi0::Com_Multi(const double a1[2], const double a2[2], double res[2])
797 {
798 res[0] = a1[0] * a2[0] - a1[1] * a2[1];
799 res[1] = a1[1] * a2[0] + a1[0] * a2[1];
800 }
801 void EvtDToKSpipi0pi0::Com_Divide(const double a1[2], const double a2[2], double res[2])
802 {
803 double tmp = a2[0] * a2[0] + a2[1] * a2[1];
804 res[0] = (a1[0] * a2[0] + a1[1] * a2[1]) / tmp;
805 res[1] = (a1[1] * a2[0] - a1[0] * a2[1]) / tmp;
806 }
807 double EvtDToKSpipi0pi0::SCADot(const double a1[4], const double a2[4])
808 {
809 double _cal = a1[0] * a2[0] - a1[1] * a2[1] - a1[2] * a2[2] - a1[3] * a2[3];
810 return _cal;
811 }
812 double EvtDToKSpipi0pi0::Barrier(const double mass2, const int l, const double sa, const double sb, const double sc,
813 const double r2)
814 {
815 double F;
816 double tmp = sa + sb - sc;
817 double q = fabs(0.25 * tmp * tmp / sa - sb);
818 double tmp2 = mass2 + sb - sc;
819 double q0 = fabs(0.25 * tmp2 * tmp2 / mass2 - sb);
820 double z = q * r2;
821 double z0 = q0 * r2;
822 if (l == 1) {
823 F = sqrt((1.0 + z0) / (1.0 + z));
824 } else if (l == 2) {
825 double z2 = z * z; double z02 = z0 * z0;
826 F = sqrt((9.0 + 3.0 * z0 + z02) / (9.0 + 3.0 * z + z2));
827 } else {
828 F = 1.0;
829 }
830 return F;
831 }
832 void EvtDToKSpipi0pi0::calt1(const double daug1[4], const double daug2[4], double t1[4])
833 {
834 double p, pq, tmp;
835 double pa[4], qa[4];
836 for (int i = 0; i < 4; i++) {
837 pa[i] = daug1[i] + daug2[i];
838 qa[i] = daug1[i] - daug2[i];
839 }
840 p = SCADot(pa, pa);
841 pq = SCADot(pa, qa);
842 tmp = pq / p;
843 for (int i = 0; i < 4; i++) {
844 t1[i] = qa[i] - tmp * pa[i];
845 }
846 }
847 void EvtDToKSpipi0pi0::calt2(const double daug1[4], const double daug2[4], double t2[4][4])
848 {
849 double p, r;
850 double pa[4], t1[4];
851 calt1(daug1, daug2, t1);
852 r = SCADot(t1, t1) / 3.0;
853 for (int i = 0; i < 4; i++) {
854 pa[i] = daug1[i] + daug2[i];
855 }
856 p = SCADot(pa, pa);
857 for (int i = 0; i < 4; i++) {
858 for (int j = 0; j < 4; j++) {
859 t2[i][j] = t1[i] * t1[j] - r * (G[i][j] - pa[i] * pa[j] / p);
860 }
861 }
862 }
863 double EvtDToKSpipi0pi0::wid(const double mass2, const double mass, const double sa, const double sb, const double sc,
864 const double r2, const int l)
865 {
866 double widm = 0.;
867 double m = sqrt(sa);
868 double tmp = sb - sc;
869 double tmp1 = sa + tmp;
870 double q = fabs(0.25 * tmp1 * tmp1 / sa - sb);
871 double tmp2 = mass2 + tmp;
872 double q0 = fabs(0.25 * tmp2 * tmp2 / mass2 - sb);
873 double z = q * r2;
874 double z0 = q0 * r2;
875 double t = q / q0;
876 if (l == 0) {widm = sqrt(t) * mass / m;}
877 else if (l == 1) {widm = t * sqrt(t) * mass / m * (1 + z0) / (1 + z);}
878 else if (l == 2) {widm = t * t * sqrt(t) * mass / m * (9 + 3 * z0 + z0 * z0) / (9 + 3 * z + z * z);}
879 return widm;
880 }
881 double EvtDToKSpipi0pi0::widl1(const double mass2, const double mass, const double sa, const double sb, const double sc,
882 const double r2)
883 {
884 double widm = 0.;
885 double m = sqrt(sa);
886 double tmp = sb - sc;
887 double tmp1 = sa + tmp;
888 double q = fabs(0.25 * tmp1 * tmp1 / sa - sb);
889 double tmp2 = mass2 + tmp;
890 double q0 = fabs(0.25 * tmp2 * tmp2 / mass2 - sb);
891 double z = q * r2;
892 double z0 = q0 * r2;
893 double F = (1 + z0) / (1 + z);
894 double t = q / q0;
895 widm = t * sqrt(t) * mass / m * F;
896 return widm;
897 }
898 void EvtDToKSpipi0pi0::propagatorRBW(const double mass2, const double mass, const double width, const double sa, const double sb,
899 const double sc, const double r2, const int l, double prop[2])
900 {
901 double a[2], b[2];
902 a[0] = 1;
903 a[1] = 0;
904 b[0] = mass2 - sa;
905 b[1] = -mass * width * wid(mass2, mass, sa, sb, sc, r2, l);
906 Com_Divide(a, b, prop);
907 }
908 void EvtDToKSpipi0pi0::propagatorRBWl1(const double mass2, const double mass, const double width, const double sa, const double sb,
909 const double sc, const double r2, double prop[2])
910 {
911 double a[2], b[2];
912 a[0] = 1;
913 a[1] = 0;
914 b[0] = mass2 - sa;
915 b[1] = -mass * width * widl1(mass2, mass, sa, sb, sc, r2);
916 Com_Divide(a, b, prop);
917 }
918 void EvtDToKSpipi0pi0::propagatorRBW_a1(const double mass2, const double mass, const double width, const double sa, const double sb,
919 const double sc, const double r2, const int l, double prop[2])
920 {
921 (void)width;
922 (void)sb;
923 (void)sc;
924 (void)r2;
925 (void)l;
926 double a[2], b[2];
927 int iii = int(sqrt(sa) * 1000) - 1;
928 double a1width[3000] = {
929 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
930 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
931 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
932 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
933 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
934 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
935 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
936 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
937 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
938 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
939 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
940 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
941 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
942 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
943 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
944 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
945 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
946 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
947 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
948 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
949 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
950 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
951 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
952 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
953 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
954 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
955 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
956 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
957 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
958 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
959 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
960 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
961 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
962 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
963 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
964 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
965 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
966 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
967 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
968 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
969 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
970 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
971 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
972 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
973 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
974 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
975 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
976 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
977 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
978 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
979 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
980 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
981 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000,
982 0.000000, 0.000001, 0.000001, 0.000001, 0.000002, 0.000002, 0.000002, 0.000003,
983 0.000004, 0.000004, 0.000005, 0.000006, 0.000007, 0.000008, 0.000009, 0.000010,
984 0.000011, 0.000012, 0.000014, 0.000015, 0.000017, 0.000019, 0.000021, 0.000023,
985 0.000025, 0.000027, 0.000029, 0.000032, 0.000035, 0.000038, 0.000041, 0.000044,
986 0.000047, 0.000050, 0.000054, 0.000058, 0.000062, 0.000066, 0.000070, 0.000075,
987 0.000079, 0.000084, 0.000089, 0.000094, 0.000100, 0.000105, 0.000111, 0.000117,
988 0.000124, 0.000130, 0.000137, 0.000143, 0.000151, 0.000158, 0.000165, 0.000173,
989 0.000182, 0.000190, 0.000199, 0.000207, 0.000216, 0.000225, 0.000235, 0.000245,
990 0.000256, 0.000266, 0.000277, 0.000288, 0.000300, 0.000311, 0.000322, 0.000335,
991 0.000347, 0.000360, 0.000373, 0.000385, 0.000400, 0.000415, 0.000429, 0.000442,
992 0.000457, 0.000473, 0.000488, 0.000504, 0.000520, 0.000539, 0.000555, 0.000572,
993 0.000590, 0.000608, 0.000626, 0.000646, 0.000664, 0.000684, 0.000704, 0.000725,
994 0.000745, 0.000766, 0.000787, 0.000809, 0.000828, 0.000854, 0.000878, 0.000901,
995 0.000927, 0.000952, 0.000973, 0.001001, 0.001027, 0.001048, 0.001080, 0.001104,
996 0.001132, 0.001159, 0.001189, 0.001219, 0.001245, 0.001277, 0.001308, 0.001338,
997 0.001370, 0.001404, 0.001433, 0.001468, 0.001498, 0.001533, 0.001570, 0.001600,
998 0.001638, 0.001678, 0.001711, 0.001745, 0.001780, 0.001825, 0.001857, 0.001898,
999 0.001941, 0.001972, 0.002017, 0.002065, 0.002104, 0.002146, 0.002189, 0.002234,
1000 0.002277, 0.002319, 0.002369, 0.002410, 0.002461, 0.002511, 0.002557, 0.002605,
1001 0.002661, 0.002704, 0.002762, 0.002807, 0.002855, 0.002910, 0.002965, 0.003020,
1002 0.003074, 0.003127, 0.003178, 0.003228, 0.003288, 0.003351, 0.003409, 0.003471,
1003 0.003532, 0.003598, 0.003660, 0.003720, 0.003793, 0.003854, 0.003910, 0.003972,
1004 0.004050, 0.004108, 0.004181, 0.004254, 0.004309, 0.004380, 0.004464, 0.004533,
1005 0.004603, 0.004679, 0.004756, 0.004811, 0.004898, 0.004974, 0.005048, 0.005142,
1006 0.005215, 0.005279, 0.005363, 0.005449, 0.005533, 0.005604, 0.005695, 0.005783,
1007 0.005869, 0.005971, 0.006060, 0.006142, 0.006247, 0.006332, 0.006409, 0.006502,
1008 0.006594, 0.006713, 0.006784, 0.006889, 0.006995, 0.007079, 0.007190, 0.007303,
1009 0.007381, 0.007487, 0.007592, 0.007710, 0.007801, 0.007910, 0.008032, 0.008149,
1010 0.008247, 0.008378, 0.008462, 0.008559, 0.008706, 0.008843, 0.008943, 0.009091,
1011 0.009207, 0.009308, 0.009448, 0.009555, 0.009698, 0.009810, 0.009936, 0.010020,
1012 0.010186, 0.010320, 0.010474, 0.010611, 0.010742, 0.010840, 0.011011, 0.011167,
1013 0.011281, 0.011395, 0.011541, 0.011714, 0.011853, 0.012046, 0.012169, 0.012277,
1014 0.012460, 0.012617, 0.012805, 0.012922, 0.013072, 0.013234, 0.013389, 0.013561,
1015 0.013704, 0.013917, 0.014025, 0.014239, 0.014425, 0.014600, 0.014716, 0.014958,
1016 0.015114, 0.015325, 0.015488, 0.015630, 0.015797, 0.016035, 0.016206, 0.016404,
1017 0.016591, 0.016842, 0.016964, 0.017199, 0.017392, 0.017557, 0.017798, 0.017987,
1018 0.018178, 0.018337, 0.018631, 0.018829, 0.019008, 0.019221, 0.019467, 0.019698,
1019 0.019941, 0.020166, 0.020379, 0.020585, 0.020806, 0.021040, 0.021309, 0.021482,
1020 0.021764, 0.022046, 0.022306, 0.022478, 0.022736, 0.023049, 0.023222, 0.023536,
1021 0.023756, 0.024090, 0.024337, 0.024557, 0.024873, 0.025099, 0.025392, 0.025682,
1022 0.025995, 0.026291, 0.026498, 0.026927, 0.027119, 0.027377, 0.027804, 0.028135,
1023 0.028279, 0.028682, 0.028871, 0.029355, 0.029531, 0.029956, 0.030243, 0.030592,
1024 0.030873, 0.031246, 0.031494, 0.031771, 0.032167, 0.032515, 0.032881, 0.033211,
1025 0.033653, 0.033988, 0.034394, 0.034639, 0.035055, 0.035569, 0.035879, 0.036211,
1026 0.036611, 0.036932, 0.037489, 0.037779, 0.038284, 0.038723, 0.039001, 0.039574,
1027 0.039854, 0.040274, 0.040881, 0.041117, 0.041644, 0.042055, 0.042531, 0.043030,
1028 0.043354, 0.043832, 0.044277, 0.044956, 0.045284, 0.045828, 0.046440, 0.046800,
1029 0.047518, 0.047727, 0.048258, 0.048850, 0.049316, 0.049992, 0.050486, 0.050987,
1030 0.051410, 0.051928, 0.052613, 0.053110, 0.053824, 0.054351, 0.055078, 0.055654,
1031 0.056030, 0.056763, 0.057245, 0.057832, 0.058569, 0.059292, 0.060048, 0.060569,
1032 0.061056, 0.061869, 0.062612, 0.063186, 0.063886, 0.064655, 0.065198, 0.065815,
1033 0.066649, 0.067577, 0.068012, 0.068967, 0.069630, 0.070181, 0.070786, 0.071989,
1034 0.072764, 0.073466, 0.074461, 0.075093, 0.075994, 0.076834, 0.077455, 0.078709,
1035 0.079581, 0.080408, 0.080884, 0.081965, 0.082882, 0.083658, 0.084824, 0.085513,
1036 0.086662, 0.087602, 0.088678, 0.089492, 0.090641, 0.091369, 0.092494, 0.093484,
1037 0.094615, 0.095385, 0.096168, 0.097668, 0.098611, 0.099630, 0.100772, 0.102020,
1038 0.103145, 0.104110, 0.105071, 0.106604, 0.107791, 0.108451, 0.109509, 0.111356,
1039 0.112026, 0.113921, 0.114507, 0.116071, 0.117027, 0.118213, 0.120164, 0.120701,
1040 0.122121, 0.123894, 0.124937, 0.126134, 0.127391, 0.128882, 0.130056, 0.131649,
1041 0.133046, 0.134275, 0.135119, 0.137072, 0.138476, 0.139612, 0.140388, 0.142734,
1042 0.143576, 0.145445, 0.147414, 0.148856, 0.149891, 0.150963, 0.152477, 0.153717,
1043 0.155275, 0.156859, 0.158462, 0.159257, 0.161865, 0.163182, 0.164465, 0.165538,
1044 0.167003, 0.169257, 0.171211, 0.172093, 0.173261, 0.174639, 0.176510, 0.177684,
1045 0.179077, 0.181041, 0.182446, 0.184769, 0.184926, 0.186741, 0.188844, 0.190884,
1046 0.191714, 0.192254, 0.193921, 0.195917, 0.196766, 0.199052, 0.200603, 0.201808,
1047 0.202699, 0.204636, 0.205712, 0.206849, 0.208741, 0.209424, 0.211698, 0.212753,
1048 0.215516, 0.215857, 0.217790, 0.217774, 0.220454, 0.221821, 0.223466, 0.224494,
1049 0.225632, 0.227231, 0.229456, 0.229581, 0.231537, 0.232263, 0.233834, 0.234725,
1050 0.237079, 0.238015, 0.239400, 0.240193, 0.241693, 0.243787, 0.244317, 0.244971,
1051 0.246711, 0.248615, 0.249387, 0.250905, 0.252702, 0.253535, 0.254385, 0.255375,
1052 0.256671, 0.258405, 0.259741, 0.260875, 0.262131, 0.262920, 0.264860, 0.265893,
1053 0.266016, 0.267727, 0.270039, 0.270689, 0.271047, 0.272313, 0.272474, 0.274724,
1054 0.275813, 0.275937, 0.278793, 0.278783, 0.281407, 0.281351, 0.282481, 0.284226,
1055 0.284113, 0.284999, 0.285655, 0.288361, 0.287856, 0.288893, 0.290211, 0.291708,
1056 0.291985, 0.294298, 0.294849, 0.296796, 0.296197, 0.296851, 0.298011, 0.300368,
1057 0.299982, 0.302378, 0.304363, 0.303711, 0.304729, 0.306789, 0.306378, 0.307372,
1058 0.308720, 0.309509, 0.309712, 0.310782, 0.311699, 0.312668, 0.312755, 0.313675,
1059 0.315311, 0.316640, 0.317217, 0.317403, 0.318478, 0.319916, 0.321803, 0.322678,
1060 0.323237, 0.324343, 0.324433, 0.324493, 0.324969, 0.325894, 0.328563, 0.328721,
1061 0.328954, 0.330640, 0.328164, 0.331267, 0.331695, 0.333772, 0.333619, 0.334351,
1062 0.334605, 0.336434, 0.337510, 0.336535, 0.337362, 0.338799, 0.340732, 0.339896,
1063 0.342707, 0.343471, 0.342318, 0.342431, 0.344543, 0.345611, 0.345786, 0.346590,
1064 0.346610, 0.347761, 0.348914, 0.349558, 0.350577, 0.352128, 0.350982, 0.354134,
1065 0.352773, 0.353213, 0.352972, 0.354927, 0.355784, 0.355778, 0.355801, 0.357040,
1066 0.358013, 0.358432, 0.360045, 0.359743, 0.360238, 0.359850, 0.362184, 0.361580,
1067 0.363430, 0.362333, 0.364397, 0.364472, 0.364370, 0.365303, 0.366644, 0.367777,
1068 0.368604, 0.367631, 0.368324, 0.369782, 0.371121, 0.370653, 0.370040, 0.371649,
1069 0.370201, 0.373362, 0.373900, 0.374159, 0.374916, 0.374503, 0.376703, 0.372802,
1070 0.376191, 0.379596, 0.377325, 0.376363, 0.379369, 0.379791, 0.378703, 0.380177,
1071 0.381762, 0.381335, 0.381374, 0.384668, 0.381763, 0.382746, 0.384723, 0.385089,
1072 0.386229, 0.386702, 0.387749, 0.384423, 0.384714, 0.384181, 0.388489, 0.388618,
1073 0.388179, 0.390092, 0.389871, 0.390496, 0.391181, 0.390679, 0.392614, 0.392269,
1074 0.393899, 0.393466, 0.391421, 0.391090, 0.395586, 0.391776, 0.396882, 0.393254,
1075 0.394400, 0.395749, 0.398063, 0.397138, 0.397585, 0.397288, 0.397847, 0.395375,
1076 0.400170, 0.400007, 0.401191, 0.398513, 0.401922, 0.400477, 0.404257, 0.403271,
1077 0.400677, 0.403913, 0.403172, 0.404727, 0.403406, 0.404404, 0.405265, 0.406389,
1078 0.405738, 0.402173, 0.407831, 0.405895, 0.409172, 0.408934, 0.405915, 0.408486,
1079 0.407320, 0.407437, 0.405444, 0.408400, 0.410909, 0.412427, 0.409881, 0.411021,
1080 0.413001, 0.410369, 0.414702, 0.413372, 0.413095, 0.410972, 0.416346, 0.416095,
1081 0.414132, 0.414344, 0.416952, 0.415197, 0.417583, 0.416582, 0.416622, 0.416895,
1082 0.416576, 0.415551, 0.417925, 0.414838, 0.417051, 0.416831, 0.420000, 0.419132,
1083 0.418173, 0.417645, 0.419679, 0.419866, 0.419581, 0.421531, 0.420878, 0.422737,
1084 0.421872, 0.421304, 0.425486, 0.424434, 0.420842, 0.426753, 0.422761, 0.422178,
1085 0.422372, 0.424173, 0.425582, 0.425080, 0.425831, 0.423551, 0.422949, 0.425784,
1086 0.427977, 0.427948, 0.426368, 0.425138, 0.425351, 0.428643, 0.428148, 0.427488,
1087 0.431704, 0.430167, 0.429655, 0.429584, 0.425458, 0.430728, 0.429845, 0.431145,
1088 0.429180, 0.428874, 0.430720, 0.430024, 0.432034, 0.431359, 0.431535, 0.432995,
1089 0.432425, 0.432454, 0.433140, 0.432574, 0.433814, 0.433348, 0.432886, 0.435472,
1090 0.436517, 0.432681, 0.436999, 0.435182, 0.434834, 0.435478, 0.438255, 0.436650,
1091 0.434464, 0.438530, 0.434077, 0.436471, 0.434012, 0.436822, 0.437505, 0.440135,
1092 0.438322, 0.438032, 0.439001, 0.440270, 0.438661, 0.439233, 0.439274, 0.437945,
1093 0.443080, 0.439191, 0.438233, 0.440415, 0.441063, 0.440926, 0.440929, 0.439731,
1094 0.443584, 0.439729, 0.441597, 0.442615, 0.444637, 0.443180, 0.440789, 0.440261,
1095 0.442202, 0.445081, 0.445484, 0.445415, 0.445532, 0.442806, 0.444188, 0.441073,
1096 0.444299, 0.445897, 0.445279, 0.442830, 0.445506, 0.445272, 0.447267, 0.443522,
1097 0.445519, 0.446459, 0.446753, 0.446377, 0.446129, 0.446383, 0.448556, 0.446593,
1098 0.445293, 0.449199, 0.447590, 0.445968, 0.447482, 0.448474, 0.449890, 0.450004,
1099 0.447765, 0.449274, 0.450652, 0.448210, 0.449360, 0.449577, 0.448575, 0.452112,
1100 0.448780, 0.451393, 0.450200, 0.452018, 0.451182, 0.452050, 0.451748, 0.451377,
1101 0.451402, 0.448810, 0.452311, 0.452909, 0.452491, 0.452418, 0.454190, 0.454420,
1102 0.452121, 0.452307, 0.456857, 0.453506, 0.454058, 0.457203, 0.454394, 0.453596,
1103 0.452240, 0.453692, 0.456516, 0.453753, 0.455541, 0.452702, 0.456481, 0.452226,
1104 0.454280, 0.454855, 0.456297, 0.456482, 0.454154, 0.455387, 0.454748, 0.455764,
1105 0.457282, 0.455487, 0.454822, 0.454257, 0.457678, 0.454225, 0.458689, 0.456123,
1106 0.457011, 0.457386, 0.458351, 0.458638, 0.456164, 0.455884, 0.458525, 0.457575,
1107 0.458340, 0.458912, 0.457836, 0.461734, 0.457545, 0.460755, 0.460960, 0.459226,
1108 0.458613, 0.461078, 0.460958, 0.460337, 0.460237, 0.461190, 0.460760, 0.457911,
1109 0.461310, 0.459657, 0.461960, 0.461040, 0.459578, 0.461650, 0.461550, 0.461251,
1110 0.461054, 0.463082, 0.461732, 0.461324, 0.462547, 0.461261, 0.461629, 0.464067,
1111 0.462430, 0.462525, 0.464232, 0.462921, 0.463202, 0.465558, 0.462914, 0.461698,
1112 0.463963, 0.463040, 0.464275, 0.461940, 0.462913, 0.465261, 0.461500, 0.463679,
1113 0.463354, 0.465205, 0.464529, 0.462220, 0.464279, 0.463427, 0.465387, 0.465288,
1114 0.464839, 0.464926, 0.466100, 0.465531, 0.466187, 0.464647, 0.466285, 0.465461,
1115 0.464134, 0.466783, 0.466763, 0.466183, 0.467089, 0.464497, 0.466080, 0.466109,
1116 0.468166, 0.466984, 0.465335, 0.466721, 0.466856, 0.465113, 0.468377, 0.467904,
1117 0.464546, 0.468787, 0.465648, 0.469841, 0.469477, 0.466311, 0.468700, 0.465183,
1118 0.466559, 0.470433, 0.468563, 0.468109, 0.466980, 0.467567, 0.467670, 0.466991,
1119 0.467992, 0.468784, 0.469406, 0.469652, 0.468527, 0.470460, 0.467308, 0.470693,
1120 0.469539, 0.468000, 0.469295, 0.467038, 0.471908, 0.468829, 0.470663, 0.469266,
1121 0.468975, 0.470222, 0.468649, 0.469507, 0.472307, 0.471611, 0.470419, 0.471181,
1122 0.471140, 0.473187, 0.471086, 0.469801, 0.472234, 0.472131, 0.468996, 0.470229,
1123 0.471597, 0.469625, 0.472230, 0.470164, 0.468404, 0.472264, 0.471336, 0.471597,
1124 0.472280, 0.471256, 0.473151, 0.471863, 0.474458, 0.471956, 0.473099, 0.473956,
1125 0.471725, 0.472809, 0.473065, 0.473180, 0.470611, 0.473614, 0.474263, 0.472792,
1126 0.473543, 0.472656, 0.469728, 0.473431, 0.474538, 0.475322, 0.474962, 0.473598,
1127 0.474114, 0.473486, 0.472934, 0.473252, 0.477149, 0.471719, 0.476383, 0.473076,
1128 0.473952, 0.473104, 0.472459, 0.474433, 0.474494, 0.473588, 0.473839, 0.478113,
1129 0.472435, 0.475571, 0.475194, 0.475626, 0.474617, 0.474520, 0.474472, 0.476437,
1130 0.474512, 0.474497, 0.474628, 0.476203, 0.475698, 0.473907, 0.477144, 0.479000,
1131 0.475553, 0.477481, 0.473998, 0.476672, 0.477115, 0.477114, 0.476282, 0.476152,
1132 0.477009, 0.479854, 0.474354, 0.477645, 0.477517, 0.477111, 0.474843, 0.476173,
1133 0.477321, 0.477384, 0.477880, 0.475726, 0.476004, 0.478204, 0.475586, 0.477973,
1134 0.477935, 0.480640, 0.478234, 0.476349, 0.477493, 0.476994, 0.479815, 0.477771,
1135 0.476333, 0.476325, 0.478245, 0.477284, 0.479238, 0.478339, 0.478966, 0.478012,
1136 0.479304, 0.480148, 0.476125, 0.481267, 0.479801, 0.476720, 0.478898, 0.479284,
1137 0.479153, 0.480157, 0.478681, 0.479712, 0.478993, 0.479943, 0.478349, 0.478930,
1138 0.478052, 0.477173, 0.479244, 0.480454, 0.479128, 0.480530, 0.477843, 0.478369,
1139 0.478561, 0.478639, 0.479191, 0.481763, 0.481321, 0.480979, 0.479702, 0.479777,
1140 0.479384, 0.477571, 0.481880, 0.478615, 0.481303, 0.478783, 0.479384, 0.480517,
1141 0.481928, 0.481199, 0.479041, 0.479188, 0.481491, 0.482840, 0.478766, 0.481941,
1142 0.481298, 0.478105, 0.482933, 0.479744, 0.483361, 0.482332, 0.482556, 0.482057,
1143 0.483616, 0.480599, 0.482245, 0.481091, 0.480871, 0.481938, 0.480678, 0.481851,
1144 0.482902, 0.482158, 0.480187, 0.481772, 0.484967, 0.483094, 0.482133, 0.483929,
1145 0.483354, 0.483382, 0.483964, 0.479941, 0.481375, 0.480255, 0.482184, 0.482541,
1146 0.482032, 0.483484, 0.479492, 0.483305, 0.481070, 0.483573, 0.485689, 0.485767,
1147 0.484221, 0.481365, 0.482440, 0.481507, 0.483418, 0.480978
1148 };
1149 double width_a1 = a1width[iii];
1150 a[0] = 1;
1151 a[1] = 0;
1152 b[0] = mass2 - sa;
1153 b[1] = -mass * width_a1;
1154 Com_Divide(a, b, prop);
1155 }
1156 void EvtDToKSpipi0pi0::propagatorGS(const double mass2, const double mass, const double width, const double sa, const double sb,
1157 const double sc, const double r2, double prop[2])
1158 {
1159
1160 double GS1 = 0.636619783;
1161 double GS2 = 0.01860182466;
1162 double GS3 = 0.1591549458;
1163 double GS4 = 0.00620060822;
1164 double a[2], b[2];
1165 double tmp = sb - sc;
1166 double tmp1 = sa + tmp;
1167 double q2 = fabs(0.25 * tmp1 * tmp1 / sa - sb);
1168 double tmp2 = mass2 + tmp;
1169 double q02 = fabs(0.25 * tmp2 * tmp2 / mass2 - sb);
1170
1171 double q = sqrt(q2);
1172 double q0 = sqrt(q02);
1173 double m = sqrt(sa);
1174 double q03 = q0 * q02;
1175 double tmp3 = log(mass + 2 * q0) + 1.2760418309; // log(mass_2Pion) = 1.2760418309;
1176
1177 double h = GS1 * q / m * (log(m + 2 * q) + 1.2760418309);
1178 double h0 = GS1 * q0 / mass * tmp3;
1179 double dh = h0 * (0.125 / q02 - 0.5 / mass2) + GS3 / mass2;
1180 double d = GS2 / q02 * tmp3 + GS3 * mass / q0 - GS4 * mass / q03;
1181 double f = mass2 / q03 * (q2 * (h - h0) + (mass2 - sa) * q02 * dh);
1182
1183 a[0] = 1.0 + d * width / mass;
1184 a[1] = 0.0;
1185 b[0] = mass2 - sa + width * f;
1186 b[1] = -mass * width * widl1(mass2, mass, sa, sb, sc, r2);
1187 Com_Divide(a, b, prop);
1188 }
1189 void EvtDToKSpipi0pi0::rhoab(const double sa, const double sb, const double sc, double res[2])
1190 {
1191 double tmp = sa + sb - sc;
1192 double q = 0.25 * tmp * tmp / sa - sb;
1193 if (q >= 0) {
1194 res[0] = 2.0 * sqrt(q / sa);
1195 res[1] = 0.0;
1196 } else {
1197 res[0] = 0.0;
1198 res[1] = 2.0 * sqrt(-q / sa);
1199 }
1200 }
1201 void EvtDToKSpipi0pi0::rho4Pi(const double sa, double res[2])
1202 {
1203 double temp = 1.0 - 0.3116765584 / sa; // 0.3116765584=0.13957*0.13957*16
1204 if (temp >= 0) {
1205 res[0] = sqrt(temp) / (1.0 + exp(9.8 - 3.5 * sa));
1206 res[1] = 0.0;
1207 } else {
1208 res[0] = 0.0;
1209 res[1] = sqrt(-temp) / (1.0 + exp(9.8 - 3.5 * sa));
1210 }
1211 }
1212 void EvtDToKSpipi0pi0::propagatorsigma500(const double sa, const double sb, const double sc, double prop[2])
1213 {
1214 double f = 0.5843 + 1.6663 * sa;
1215 const double M = 0.9264;
1216 const double mass2 = 0.85821696; // M*M
1217 const double mpi2d2 = 0.00973989245;
1218 double g1 = f * (sa - mpi2d2) / (mass2 - mpi2d2) * exp((mass2 - sa) / 1.082);
1219 double rho1s[2], rho1M[2], rho2s[2], rho2M[2], rho1[2], rho2[2];
1220 rhoab(sa, sb, sc, rho1s);
1221 rhoab(mass2, sb, sc, rho1M);
1222 rho4Pi(sa, rho2s);
1223 rho4Pi(mass2, rho2M);
1224 Com_Divide(rho1s, rho1M, rho1);
1225 Com_Divide(rho2s, rho2M, rho2);
1226 double a[2], b[2];
1227 a[0] = 1.0;
1228 a[1] = 0.0;
1229 b[0] = mass2 - sa + M * (g1 * rho1[1] + 0.0024 * rho2[1]);
1230 b[1] = -M * (g1 * rho1[0] + 0.0024 * rho2[0]);
1231 Com_Divide(a, b, prop);
1232 }
1233 void EvtDToKSpipi0pi0::KPiSLASS(const double sa, const double sb, const double sc, double prop[2])
1234 {
1235 const double m1430 = 1.441;
1236 const double sa0 = 2.076481; // m1430*m1430;
1237 const double w1430 = 0.193;
1238 const double Lass1 = 0.25 / sa0;
1239 double tmp = sb - sc;
1240 double tmp1 = sa0 + tmp;
1241 double q0 = fabs(Lass1 * tmp1 * tmp1 - sb);
1242 double tmp2 = sa + tmp;
1243 double qs = fabs(0.25 * tmp2 * tmp2 / sa - sb);
1244 double q = sqrt(qs);
1245 double width = w1430 * q * m1430 / sqrt(sa * q0);
1246 double temp_R = atan(m1430 * width / (sa0 - sa));
1247 if (temp_R < 0) temp_R += math_pi;
1248 double deltaR = -109.7 * math_pi / 180.0 + temp_R;
1249 double temp_F = atan(0.226 * q / (2.0 - 3.8194 * qs)); // 2.0*0.113 = 0.226; -33.8*0.113 = -3.8194
1250 if (temp_F < 0) temp_F += math_pi;
1251 double deltaF = 0.1 * math_pi / 180.0 + temp_F;
1252 double deltaS = deltaR + 2.0 * deltaF;
1253 double t1 = 0.96 * sin(deltaF);
1254 double t2 = sin(deltaR);
1255 double CF[2], CS[2];
1256 CF[0] = cos(deltaF);
1257 CF[1] = sin(deltaF);
1258 CS[0] = cos(deltaS);
1259 CS[1] = sin(deltaS);
1260 prop[0] = t1 * CF[0] + t2 * CS[0];
1261 prop[1] = t1 * CF[1] + t2 * CS[1];
1262
1263 }
1264
1265} // Belle 2 namespace