BOSS 8.0.0
BESIII Offline Software System
Loading...
Searching...
No Matches
EvtD0ToKSpipipi0.cc
Go to the documentation of this file.
1//--------------------------------------------------------------------------
2//
3// Environment:
4// This software is part of models developed at BES collaboration
5// based on the EvtGen framework. If you use all or part
6// of it, please give an appropriate acknowledgement.
7//
8// Copyright Information: See EvtGen/BesCopyright
9// Copyright (A) 2006 Ping Rong-Gang @IHEP
10//
11// Module: EvtD0ToKSpipipi0.cc
12//
13// Description: Model provided by user, see
14// https://indico.ihep.ac.cn/event/8797/contributions/104436/attachments/55909/64345/Ks3pi-20190703.pdf
15//
16// Modification history:
17//
18// Liaoyuan Dong Aug. 11, 2022 Module created
19//
20//------------------------------------------------------------------------
21#include "EvtD0ToKSpipipi0.hh"
22#include "../EvtGenBase/EvtGenKine.hh"
23#include "../EvtGenBase/EvtPDL.hh"
24#include "../EvtGenBase/EvtParticle.hh"
25#include "../EvtGenBase/EvtPatches.hh"
26#include "../EvtGenBase/EvtReport.hh"
27#include "../EvtGenBase/EvtVector4R.hh"
28#include "TLorentzVector.h"
29#include "TMath.h"
30#include <complex>
31#include <map>
32#include <math.h>
33#include <stdlib.h>
34#include <string>
35#include <vector>
36
37#include <cassert>
38#include <iostream>
39
40using namespace std; //::endl;
41
42template <typename T> string toString( const T& t ) {
43 stringstream oss;
44 oss << t;
45 return oss.str();
46}
47
48std::vector<double> operator+( std::vector<double>& lhs, std::vector<double>& rhs ) {
49 assert( int( lhs.size() ) == int( rhs.size() ) );
50 std::vector<double> sum;
51 sum.clear();
52
53 for ( int i = 0; i < int( lhs.size() ); i++ ) { sum.push_back( lhs[i] + rhs[i] ); }
54
55 return sum;
56}
57
58EvtD0ToKSpipipi0::~EvtD0ToKSpipipi0() {}
59
60void EvtD0ToKSpipipi0::getName( std::string& model_name ) { model_name = "D0ToKSpipipi0"; }
61
62EvtDecayBase* EvtD0ToKSpipipi0::clone() { return new EvtD0ToKSpipipi0; }
63
64void EvtD0ToKSpipipi0::init() {
65 checkNArg( 0 );
66 checkNDaug( 4 );
67 checkSpinParent( EvtSpinType::SCALAR );
68 // std::cout << "Initializing EvtD0ToKSpipipi0" << std::endl;
69
70 _nd = 4;
71 rRes = 3.0;
72 rD = 5.0;
73 math_pi = 3.1415926;
74 mpip = 0.13957;
75 mpim = 0.13957;
76 mpion = 0.13957;
77 mkaon = 0.493677;
78}
79void EvtD0ToKSpipipi0::initProbMax() { setProbMax( 29338399326.9 ); }
80
81void EvtD0ToKSpipipi0::decay( EvtParticle* p ) {
82 /*
83 double maxprob = 0.0;
84 for(int ir=0;ir<=60000000;ir++){
85 p->initializePhaseSpace(getNDaug(),getDaugs());
86 for(int i=0; i<_nd; i++){
87 _p4Lab[i]=p->getDaug(i)->getP4Lab();
88 _p4CM[i]=p->getDaug(i)->getP4();
89 }
90 double Prob = AmplitudeSquare();
91 if(Prob>maxprob) {
92 maxprob=Prob;
93 std::cout << "Max PDF = " << ir << " prob= " << Prob << std::endl;
94 }
95 }
96 std::cout << "Max!!!!!!!!!!! " << maxprob<< std::endl;
97 */
98 p->initializePhaseSpace( getNDaug(), getDaugs() );
99 for ( int i = 0; i < _nd; i++ )
100 {
101 _p4Lab[i] = p->getDaug( i )->getP4Lab();
102 _p4CM[i] = p->getDaug( i )->getP4();
103 }
104 double prob = AmplitudeSquare();
105 setProb( prob );
106 return;
107}
108
109void EvtD0ToKSpipipi0::readInputFile() {
110
111 resonance_par["kstarm_892_mass"] = 0.89166;
112 resonance_par["kstarm_892_width"] = 0.0508;
113 resonance_par["kstar0_892_mass"] = 0.89581;
114 resonance_par["kstar0_892_width"] = 0.0474;
115 resonance_par["kstarp_892_mass"] = 0.89166;
116 resonance_par["kstarp_892_width"] = 0.0508;
117 resonance_par["k1m_1270_mass"] = 1.272;
118 resonance_par["k1m_1270_width"] = 0.087;
119 resonance_par["k10_1270_mass"] = 1.272;
120 resonance_par["k10_1270_width"] = 0.087;
121 resonance_par["k1m_1400_mass"] = 1.403;
122 resonance_par["k1m_1400_width"] = 0.174;
123 resonance_par["k10_1400_mass"] = 1.403;
124 resonance_par["k10_1400_width"] = 0.174;
125 resonance_par["kstarm_1410_mass"] = 1.414;
126 resonance_par["kstarm_1410_width"] = 0.232;
127 resonance_par["kstar0_1410_mass"] = 1.414;
128 resonance_par["kstar0_1410_width"] = 0.232;
129 resonance_par["k0starm_1430_mass"] = 1.425;
130 resonance_par["k0starm_1430_width"] = 0.27;
131 resonance_par["k0star0_1430_mass"] = 1.425;
132 resonance_par["k0star0_1430_width"] = 0.27;
133 resonance_par["k2starm_1430_mass"] = 1.4256;
134 resonance_par["k2starm_1430_width"] = 0.0985;
135 resonance_par["k2star0_1430_mass"] = 1.4324;
136 resonance_par["k2star0_1430_width"] = 0.109;
137 resonance_par["km_1460_mass"] = 1.46;
138 resonance_par["km_1460_width"] = 0.26;
139 resonance_par["k0_1460_mass"] = 1.46;
140 resonance_par["k0_1460_width"] = 0.26;
141 resonance_par["k2m_1580_mass"] = 1.58;
142 resonance_par["k2m_1580_width"] = 0.11;
143 resonance_par["k20_1580_mass"] = 1.58;
144 resonance_par["k20_1580_width"] = 0.11;
145 resonance_par["km_1630_mass"] = 1.595;
146 resonance_par["km_1630_width"] = 0.016;
147 resonance_par["k0_1630_mass"] = 1.595;
148 resonance_par["k0_1630_width"] = 0.016;
149 resonance_par["k1m_1650_mass"] = 1.65;
150 resonance_par["k1m_1650_width"] = 0.15;
151 resonance_par["k10_1650_mass"] = 1.65;
152 resonance_par["k10_1650_width"] = 0.15;
153 resonance_par["kstarm_1680_mass"] = 1.717;
154 resonance_par["kstarm_1680_width"] = 0.322;
155 resonance_par["kstar0_1680_mass"] = 1.717;
156 resonance_par["kstar0_1680_width"] = 0.322;
157 resonance_par["k2m_1770_mass"] = 1.773;
158 resonance_par["k2m_1770_width"] = 0.186;
159 resonance_par["k20_1770_mass"] = 1.773;
160 resonance_par["k20_1770_width"] = 0.186;
161 resonance_par["sigma_500_mass"] = 0.9264;
162 resonance_par["sigma_500_width"] = 0.45623;
163 resonance_par["sigma_500_gf"] = 0.0024;
164 resonance_par["f0_980_mass"] = 0.965;
165 resonance_par["f0_980_width"] = 0.055;
166 resonance_par["f0_980_g1"] = 0.165;
167 resonance_par["f0_980_g2"] = 4.21;
168 resonance_par["f0_1370_mass"] = 1.37;
169 resonance_par["f0_1370_width"] = 0.26;
170 resonance_par["f2_1270_mass"] = 1.275;
171 resonance_par["f2_1270_width"] = 0.185;
172 resonance_par["rhop_770_mass"] = 0.77511;
173 resonance_par["rhop_770_width"] = 0.1491;
174 resonance_par["rhom_770_mass"] = 0.77511;
175 resonance_par["rhom_770_width"] = 0.1491;
176 resonance_par["rho0_770_mass"] = 0.77526;
177 resonance_par["rho0_770_width"] = 0.1478;
178 resonance_par["omega_782_mass"] = 0.78265;
179 resonance_par["omega_782_width"] = 0.00849;
180 resonance_par["eta_547_mass"] = 0.547862;
181 resonance_par["eta_547_width"] = 0.00131;
182 resonance_par["phi_1020_mass"] = 1.01946;
183 resonance_par["phi_1020_width"] = 0.00426;
184 resonance_par["a10_1260_mass"] = 1.366;
185 resonance_par["a10_1260_width"] = 0.542;
186 resonance_par["kstarm_phsp_mass"] = 2;
187 resonance_par["kstarm_phsp_width"] = 90;
188 resonance_par["kstar0_phsp_mass"] = 2;
189 resonance_par["kstar0_phsp_width"] = 90;
190 resonance_par["kstarp_phsp_mass"] = 2;
191 resonance_par["kstarp_phsp_width"] = 90;
192 resonance_par["rhop_phsp_mass"] = 2;
193 resonance_par["rhop_phsp_width"] = 90;
194 resonance_par["rho0_phsp_mass"] = 2;
195 resonance_par["rho0_phsp_width"] = 90;
196 resonance_par["rhom_phsp_mass"] = 2;
197 resonance_par["rhom_phsp_width"] = 90;
198 resonance_par["k1m_phsp_mass"] = 2;
199 resonance_par["k1m_phsp_width"] = 90;
200 resonance_par["k10_phsp_mass"] = 2;
201 resonance_par["k10_phsp_width"] = 90;
202 resonance_par["k1p_phsp_mass"] = 2;
203 resonance_par["k1p_phsp_width"] = 90;
204 resonance_par["a10_phsp_mass"] = 2;
205 resonance_par["a10_phsp_width"] = 90;
206
207 readInputCoeff();
208}
209
210void EvtD0ToKSpipipi0::initPar() {
211
212 g.clear();
213 for ( int i = 0; i < 4; i++ )
214 {
215 for ( int j = 0; j < 4; j++ )
216 {
217 if ( i != j ) { g.push_back( 0.0 ); }
218 else if ( i < 3 ) { g.push_back( -1.0 ); }
219 else if ( i == 3 ) { g.push_back( 1.0 ); }
220 }
221 }
222
223 epsilon.clear();
224 for ( int i = 0; i < 4; i++ )
225 {
226 for ( int j = 0; j < 4; j++ )
227 {
228 for ( int k = 0; k < 4; k++ )
229 {
230 for ( int l = 0; l < 4; l++ )
231 {
232 if ( i == j || i == k || i == l || j == k || j == l || k == l )
233 { epsilon.push_back( 0.0 ); }
234 else
235 {
236 if ( i == 0 && j == 1 && k == 2 && l == 3 ) epsilon.push_back( 1.0 );
237 if ( i == 0 && j == 1 && k == 3 && l == 2 ) epsilon.push_back( -1.0 );
238 if ( i == 0 && j == 2 && k == 1 && l == 3 ) epsilon.push_back( -1.0 );
239 if ( i == 0 && j == 2 && k == 3 && l == 1 ) epsilon.push_back( 1.0 );
240 if ( i == 0 && j == 3 && k == 1 && l == 2 ) epsilon.push_back( 1.0 );
241 if ( i == 0 && j == 3 && k == 2 && l == 1 ) epsilon.push_back( -1.0 );
242
243 if ( i == 1 && j == 0 && k == 2 && l == 3 ) epsilon.push_back( -1.0 );
244 if ( i == 1 && j == 0 && k == 3 && l == 2 ) epsilon.push_back( 1.0 );
245 if ( i == 1 && j == 2 && k == 0 && l == 3 ) epsilon.push_back( 1.0 );
246 if ( i == 1 && j == 2 && k == 3 && l == 0 ) epsilon.push_back( -1.0 );
247 if ( i == 1 && j == 3 && k == 0 && l == 2 ) epsilon.push_back( -1.0 );
248 if ( i == 1 && j == 3 && k == 2 && l == 0 ) epsilon.push_back( 1.0 );
249
250 if ( i == 2 && j == 0 && k == 1 && l == 3 ) epsilon.push_back( 1.0 );
251 if ( i == 2 && j == 0 && k == 3 && l == 1 ) epsilon.push_back( -1.0 );
252 if ( i == 2 && j == 1 && k == 0 && l == 3 ) epsilon.push_back( -1.0 );
253 if ( i == 2 && j == 1 && k == 3 && l == 0 ) epsilon.push_back( 1.0 );
254 if ( i == 2 && j == 3 && k == 0 && l == 1 ) epsilon.push_back( 1.0 );
255 if ( i == 2 && j == 3 && k == 1 && l == 0 ) epsilon.push_back( -1.0 );
256
257 if ( i == 3 && j == 0 && k == 1 && l == 2 ) epsilon.push_back( -1.0 );
258 if ( i == 3 && j == 0 && k == 2 && l == 1 ) epsilon.push_back( 1.0 );
259 if ( i == 3 && j == 1 && k == 0 && l == 2 ) epsilon.push_back( 1.0 );
260 if ( i == 3 && j == 1 && k == 2 && l == 0 ) epsilon.push_back( -1.0 );
261 if ( i == 3 && j == 2 && k == 0 && l == 1 ) epsilon.push_back( -1.0 );
262 if ( i == 3 && j == 2 && k == 1 && l == 0 ) epsilon.push_back( 1.0 );
263 }
264 }
265 }
266 }
267 }
268
269 totalAmp = ( 0., 0. );
270
271 readInputFile();
272}
273
274//============================================================//
275// Function //
276//============================================================//
277void EvtD0ToKSpipipi0::addPartialWave( EvtComplex amp, double mag, double pha ) {
278
279 EvtComplex coeff( mag * cos( pha ), mag * sin( pha ) );
280 totalAmp += ( amp * coeff );
281}
282
283void EvtD0ToKSpipipi0::addPartialWave( EvtComplex amp1, EvtComplex amp2, double mag,
284 double pha ) {
285
286 EvtComplex coeff( mag * cos( pha ), mag * sin( pha ) );
287 totalAmp += ( amp1 * amp2 * coeff );
288}
289
290void EvtD0ToKSpipipi0::addPartialWave( double t1, EvtComplex resX, EvtComplex resY, double mag,
291 double pha ) {
292 EvtComplex coeff( mag * cos( pha ), mag * sin( pha ) );
293 totalAmp += ( t1 * resX * resY * coeff );
294}
295
296// For iso-spin conjugate
297void EvtD0ToKSpipipi0::addPartialWave( double t1, EvtComplex resX1, EvtComplex resY1,
298 double t2, EvtComplex resX2, EvtComplex resY2,
299 double fact, double mag, double pha ) {
300
301 EvtComplex coeff( mag * cos( pha ), mag * sin( pha ) );
302 totalAmp += ( t1 * resX1 * resY1 + fact * t2 * resX2 * resY2 ) * coeff;
303}
304
305//============================================================//
306// Spin Factor //
307//============================================================//
308void EvtD0ToKSpipipi0::createSpinfactor( EvtVector4R Ks0, EvtVector4R Pip, EvtVector4R Pim,
309 EvtVector4R Pi0 ) {
310
311 double m2_ks0 = Ks0.mass2();
312 double m2_pip = Pip.mass2();
313 double m2_pim = Pim.mass2();
314 double m2_pi0 = Pi0.mass2();
315 double m2_ks0pip = ( Ks0 + Pip ).mass2();
316 double m2_ks0pim = ( Ks0 + Pim ).mass2();
317 double m2_ks0pi0 = ( Ks0 + Pi0 ).mass2();
318 double m2_pippim = ( Pip + Pim ).mass2();
319 double m2_pippi0 = ( Pip + Pi0 ).mass2();
320 double m2_pimpi0 = ( Pim + Pi0 ).mass2();
321 double m2_ks0pippim = ( Ks0 + Pip + Pim ).mass2();
322 double m2_ks0pippi0 = ( Ks0 + Pip + Pi0 ).mass2();
323 double m2_ks0pimpi0 = ( Ks0 + Pim + Pi0 ).mass2();
324 double m2_pippimpi0 = ( Pip + Pim + Pi0 ).mass2();
325
326 std::vector<double> ks0;
327 ks0.clear();
328 ks0.push_back( Ks0.get( 1 ) );
329 ks0.push_back( Ks0.get( 2 ) );
330 ks0.push_back( Ks0.get( 3 ) );
331 ks0.push_back( Ks0.get( 0 ) );
332 std::vector<double> pip;
333 pip.clear();
334 pip.push_back( Pip.get( 1 ) );
335 pip.push_back( Pip.get( 2 ) );
336 pip.push_back( Pip.get( 3 ) );
337 pip.push_back( Pip.get( 0 ) );
338 std::vector<double> pim;
339 pim.clear();
340 pim.push_back( Pim.get( 1 ) );
341 pim.push_back( Pim.get( 2 ) );
342 pim.push_back( Pim.get( 3 ) );
343 pim.push_back( Pim.get( 0 ) );
344 std::vector<double> pi0;
345 pi0.clear();
346 pi0.push_back( Pi0.get( 1 ) );
347 pi0.push_back( Pi0.get( 2 ) );
348 pi0.push_back( Pi0.get( 3 ) );
349 pi0.push_back( Pi0.get( 0 ) );
350
351 //------------------------------------------------------------//
352 // Cascade Decay //
353 //------------------------------------------------------------//
354 //---- D2PP_P2SP
355 spinfactor["D2PP_P2SP_pippimpi0_pimpi0_0"] = D2PP_P2SP();
356 spinfactor["D2PP_P2SP_pippimpi0_pippi0_0"] = D2PP_P2SP();
357 spinfactor["D2PP_P2SP_pippimpi0_pippim_0"] = D2PP_P2SP(); //
358 spinfactor["D2PP_P2SP_ks0pimpi0_ks0pim_0"] = D2PP_P2SP();
359 spinfactor["D2PP_P2SP_ks0pimpi0_ks0pi0_0"] = D2PP_P2SP();
360 spinfactor["D2PP_P2SP_ks0pimpi0_pimpi0_0"] = D2PP_P2SP(); //
361 spinfactor["D2PP_P2SP_ks0pippi0_ks0pip_0"] = D2PP_P2SP();
362 spinfactor["D2PP_P2SP_ks0pippi0_ks0pi0_0"] = D2PP_P2SP();
363 spinfactor["D2PP_P2SP_ks0pippi0_pippi0_0"] = D2PP_P2SP(); //
364 spinfactor["D2PP_P2SP_ks0pippim_ks0pip_0"] = D2PP_P2SP();
365 spinfactor["D2PP_P2SP_ks0pippim_ks0pim_0"] = D2PP_P2SP();
366 spinfactor["D2PP_P2SP_ks0pippim_pippim_0"] = D2PP_P2SP(); //
367
368 //---- D2PP_P2VP
369 spinfactor["D2PP_P2VP_pippimpi0_pimpi0_1"] =
370 D2PP_P2VP( pi0, pim, pip, ks0, 1 ); // sequence of pi0 pi-
371 spinfactor["D2PP_P2VP_pippimpi0_pippi0_1"] = D2PP_P2VP( pip, pi0, pim, ks0, 1 );
372 spinfactor["D2PP_P2VP_pippimpi0_pippim_1"] = D2PP_P2VP( pip, pim, pi0, ks0, 1 ); //
373 spinfactor["D2PP_P2VP_ks0pimpi0_ks0pim_1"] = D2PP_P2VP( ks0, pim, pi0, pip, 1 );
374 spinfactor["D2PP_P2VP_ks0pimpi0_ks0pi0_1"] = D2PP_P2VP( ks0, pi0, pim, pip, 1 );
375 spinfactor["D2PP_P2VP_ks0pimpi0_pimpi0_1"] = D2PP_P2VP( pim, pi0, ks0, pip, 1 ); //
376 spinfactor["D2PP_P2VP_ks0pippi0_ks0pip_1"] = D2PP_P2VP( ks0, pip, pi0, pim, 1 );
377 spinfactor["D2PP_P2VP_ks0pippi0_ks0pi0_1"] = D2PP_P2VP( ks0, pi0, pip, pim, 1 );
378 spinfactor["D2PP_P2VP_ks0pippi0_pippi0_1"] = D2PP_P2VP( pip, pi0, ks0, pim, 1 ); //
379 spinfactor["D2PP_P2VP_ks0pippim_ks0pip_1"] = D2PP_P2VP( ks0, pip, pim, pi0, 1 );
380 spinfactor["D2PP_P2VP_ks0pippim_ks0pim_1"] = D2PP_P2VP( ks0, pim, pip, pi0, 1 );
381 spinfactor["D2PP_P2VP_ks0pippim_pippim_1"] = D2PP_P2VP( pip, pim, ks0, pi0, 1 ); //
382
383 //---- D2VP_V2VP
384 spinfactor["D2VP_V2VP_pippimpi0_pimpi0_1"] =
385 D2VP_V2VP( pi0, pim, pip, ks0, 1 ); // sequence of pi0 pi-
386 spinfactor["D2VP_V2VP_pippimpi0_pippi0_1"] = D2VP_V2VP( pip, pi0, pim, ks0, 1 );
387 spinfactor["D2VP_V2VP_pippimpi0_pippim_1"] = D2VP_V2VP( pip, pim, pi0, ks0, 1 ); //
388 spinfactor["D2VP_V2VP_ks0pimpi0_ks0pim_1"] = D2VP_V2VP( ks0, pim, pi0, pip, 1 );
389 spinfactor["D2VP_V2VP_ks0pimpi0_ks0pi0_1"] = D2VP_V2VP( ks0, pi0, pim, pip, 1 );
390 spinfactor["D2VP_V2VP_ks0pimpi0_pimpi0_1"] = D2VP_V2VP( pim, pi0, ks0, pip, 1 ); //
391 spinfactor["D2VP_V2VP_ks0pippi0_ks0pip_1"] = D2VP_V2VP( ks0, pip, pi0, pim, 1 );
392 spinfactor["D2VP_V2VP_ks0pippi0_ks0pi0_1"] = D2VP_V2VP( ks0, pi0, pip, pim, 1 );
393 spinfactor["D2VP_V2VP_ks0pippi0_pippi0_1"] = D2VP_V2VP( pip, pi0, ks0, pim, 1 ); //
394 spinfactor["D2VP_V2VP_ks0pippim_ks0pip_1"] = D2VP_V2VP( ks0, pip, pim, pi0, 1 );
395 spinfactor["D2VP_V2VP_ks0pippim_ks0pim_1"] = D2VP_V2VP( ks0, pim, pip, pi0, 1 );
396 spinfactor["D2VP_V2VP_ks0pippim_pippim_1"] = D2VP_V2VP( pip, pim, ks0, pi0, 1 ); //
397
398 //---- D2AP_A2SP
399 spinfactor["D2AP_A2SP_pippimpi0_pimpi0_1"] =
400 D2AP_A2SP( pi0, pim, pip, ks0, 1 ); // sequence of pi0 pi-
401 spinfactor["D2AP_A2SP_pippimpi0_pippi0_1"] = D2AP_A2SP( pip, pi0, pim, ks0, 1 );
402 spinfactor["D2AP_A2SP_pippimpi0_pippim_1"] = D2AP_A2SP( pip, pim, pi0, ks0, 1 ); //
403 spinfactor["D2AP_A2SP_ks0pimpi0_ks0pim_1"] = D2AP_A2SP( ks0, pim, pi0, pip, 1 );
404 spinfactor["D2AP_A2SP_ks0pimpi0_ks0pi0_1"] = D2AP_A2SP( ks0, pi0, pim, pip, 1 );
405 spinfactor["D2AP_A2SP_ks0pimpi0_pimpi0_1"] = D2AP_A2SP( pim, pi0, ks0, pip, 1 ); //
406 spinfactor["D2AP_A2SP_ks0pippi0_ks0pip_1"] = D2AP_A2SP( ks0, pip, pi0, pim, 1 );
407 spinfactor["D2AP_A2SP_ks0pippi0_ks0pi0_1"] = D2AP_A2SP( ks0, pi0, pip, pim, 1 );
408 spinfactor["D2AP_A2SP_ks0pippi0_pippi0_1"] = D2AP_A2SP( pip, pi0, ks0, pim, 1 ); //
409 spinfactor["D2AP_A2SP_ks0pippim_ks0pip_1"] = D2AP_A2SP( ks0, pip, pim, pi0, 1 );
410 spinfactor["D2AP_A2SP_ks0pippim_ks0pim_1"] = D2AP_A2SP( ks0, pim, pip, pi0, 1 );
411 spinfactor["D2AP_A2SP_ks0pippim_pippim_1"] = D2AP_A2SP( pip, pim, ks0, pi0, 1 ); //
412
413 //---- D2AP_A2VP, [S,D]
414 spinfactor["D2AP_A2VP_pippimpi0_pimpi0_0"] =
415 D2AP_A2VP( pi0, pim, pip, ks0, 0 ); // sequence of pi0 pi-
416 spinfactor["D2AP_A2VP_pippimpi0_pippi0_0"] = D2AP_A2VP( pip, pi0, pim, ks0, 0 );
417 spinfactor["D2AP_A2VP_pippimpi0_pippim_0"] = D2AP_A2VP( pip, pim, pi0, ks0, 0 ); //
418 spinfactor["D2AP_A2VP_ks0pimpi0_ks0pim_0"] = D2AP_A2VP( ks0, pim, pi0, pip, 0 );
419 spinfactor["D2AP_A2VP_ks0pimpi0_ks0pi0_0"] = D2AP_A2VP( ks0, pi0, pim, pip, 0 );
420 spinfactor["D2AP_A2VP_ks0pimpi0_pimpi0_0"] = D2AP_A2VP( pim, pi0, ks0, pip, 0 ); //
421 spinfactor["D2AP_A2VP_ks0pippi0_ks0pip_0"] = D2AP_A2VP( ks0, pip, pi0, pim, 0 );
422 spinfactor["D2AP_A2VP_ks0pippi0_ks0pi0_0"] = D2AP_A2VP( ks0, pi0, pip, pim, 0 );
423 spinfactor["D2AP_A2VP_ks0pippi0_pippi0_0"] = D2AP_A2VP( pip, pi0, ks0, pim, 0 ); //
424 spinfactor["D2AP_A2VP_ks0pippim_ks0pip_0"] = D2AP_A2VP( ks0, pip, pim, pi0, 0 );
425 spinfactor["D2AP_A2VP_ks0pippim_ks0pim_0"] = D2AP_A2VP( ks0, pim, pip, pi0, 0 );
426 spinfactor["D2AP_A2VP_ks0pippim_pippim_0"] = D2AP_A2VP( pip, pim, ks0, pi0, 0 ); //
427 spinfactor["D2AP_A2VP_pippimpi0_pimpi0_2"] =
428 D2AP_A2VP( pi0, pim, pip, ks0, 2 ); // sequence of pi0 pi-
429 spinfactor["D2AP_A2VP_pippimpi0_pippi0_2"] = D2AP_A2VP( pip, pi0, pim, ks0, 2 );
430 spinfactor["D2AP_A2VP_pippimpi0_pippim_2"] = D2AP_A2VP( pip, pim, pi0, ks0, 2 ); //
431 spinfactor["D2AP_A2VP_ks0pimpi0_ks0pim_2"] = D2AP_A2VP( ks0, pim, pi0, pip, 2 );
432 spinfactor["D2AP_A2VP_ks0pimpi0_ks0pi0_2"] = D2AP_A2VP( ks0, pi0, pim, pip, 2 );
433 spinfactor["D2AP_A2VP_ks0pimpi0_pimpi0_2"] = D2AP_A2VP( pim, pi0, ks0, pip, 2 ); //
434 spinfactor["D2AP_A2VP_ks0pippi0_ks0pip_2"] = D2AP_A2VP( ks0, pip, pi0, pim, 2 );
435 spinfactor["D2AP_A2VP_ks0pippi0_ks0pi0_2"] = D2AP_A2VP( ks0, pi0, pip, pim, 2 );
436 spinfactor["D2AP_A2VP_ks0pippi0_pippi0_2"] = D2AP_A2VP( pip, pi0, ks0, pim, 2 ); //
437 spinfactor["D2AP_A2VP_ks0pippim_ks0pip_2"] = D2AP_A2VP( ks0, pip, pim, pi0, 2 );
438 spinfactor["D2AP_A2VP_ks0pippim_ks0pim_2"] = D2AP_A2VP( ks0, pim, pip, pi0, 2 );
439 spinfactor["D2AP_A2VP_ks0pippim_pippim_2"] = D2AP_A2VP( pip, pim, ks0, pi0, 2 ); //
440
441 //---- D2AP_A2TP
442 spinfactor["D2AP_A2TP_pippimpi0_pimpi0_1"] =
443 D2AP_A2TP( pi0, pim, pip, ks0, 1 ); // sequence of pi0 pi-
444 spinfactor["D2AP_A2TP_pippimpi0_pippi0_1"] = D2AP_A2TP( pip, pi0, pim, ks0, 1 );
445 spinfactor["D2AP_A2TP_pippimpi0_pippim_1"] = D2AP_A2TP( pip, pim, pi0, ks0, 1 ); //
446 spinfactor["D2AP_A2TP_ks0pimpi0_ks0pim_1"] = D2AP_A2TP( ks0, pim, pi0, pip, 1 );
447 spinfactor["D2AP_A2TP_ks0pimpi0_ks0pi0_1"] = D2AP_A2TP( ks0, pi0, pim, pip, 1 );
448 spinfactor["D2AP_A2TP_ks0pimpi0_pimpi0_1"] = D2AP_A2TP( pim, pi0, ks0, pip, 1 ); //
449 spinfactor["D2AP_A2TP_ks0pippi0_ks0pip_1"] = D2AP_A2TP( ks0, pip, pi0, pim, 1 );
450 spinfactor["D2AP_A2TP_ks0pippi0_ks0pi0_1"] = D2AP_A2TP( ks0, pi0, pip, pim, 1 );
451 spinfactor["D2AP_A2TP_ks0pippi0_pippi0_1"] = D2AP_A2TP( pip, pi0, ks0, pim, 1 ); //
452 spinfactor["D2AP_A2TP_ks0pippim_ks0pip_1"] = D2AP_A2TP( ks0, pip, pim, pi0, 1 );
453 spinfactor["D2AP_A2TP_ks0pippim_ks0pim_1"] = D2AP_A2TP( ks0, pim, pip, pi0, 1 );
454 spinfactor["D2AP_A2TP_ks0pippim_pippim_1"] = D2AP_A2TP( pip, pim, ks0, pi0, 1 ); //
455
456 //---- D2TP_T2VP
457 spinfactor["D2TP_T2VP_pippimpi0_pimpi0_2"] =
458 D2TP_T2VP( pi0, pim, pip, ks0, 2 ); // sequence of pi0 pi-
459 spinfactor["D2TP_T2VP_pippimpi0_pippi0_2"] = D2TP_T2VP( pip, pi0, pim, ks0, 2 );
460 spinfactor["D2TP_T2VP_pippimpi0_pippim_2"] = D2TP_T2VP( pip, pim, pi0, ks0, 2 ); //
461 spinfactor["D2TP_T2VP_ks0pimpi0_ks0pim_2"] = D2TP_T2VP( ks0, pim, pi0, pip, 2 );
462 spinfactor["D2TP_T2VP_ks0pimpi0_ks0pi0_2"] = D2TP_T2VP( ks0, pi0, pim, pip, 2 );
463 spinfactor["D2TP_T2VP_ks0pimpi0_pimpi0_2"] = D2TP_T2VP( pim, pi0, ks0, pip, 2 ); //
464 spinfactor["D2TP_T2VP_ks0pippi0_ks0pip_2"] = D2TP_T2VP( ks0, pip, pi0, pim, 2 );
465 spinfactor["D2TP_T2VP_ks0pippi0_ks0pi0_2"] = D2TP_T2VP( ks0, pi0, pip, pim, 2 );
466 spinfactor["D2TP_T2VP_ks0pippi0_pippi0_2"] = D2TP_T2VP( pip, pi0, ks0, pim, 2 ); //
467 spinfactor["D2TP_T2VP_ks0pippim_ks0pip_2"] = D2TP_T2VP( ks0, pip, pim, pi0, 2 );
468 spinfactor["D2TP_T2VP_ks0pippim_ks0pim_2"] = D2TP_T2VP( ks0, pim, pip, pi0, 2 );
469 spinfactor["D2TP_T2VP_ks0pippim_pippim_2"] = D2TP_T2VP( pip, pim, ks0, pi0, 2 ); //
470
471 //---- D2TP_T2TP
472 spinfactor["D2TP_T2TP_pippimpi0_pimpi0_2"] =
473 D2TP_T2TP( pi0, pim, pip, ks0, 2 ); // sequence of pi0 pi-
474 spinfactor["D2TP_T2TP_pippimpi0_pippi0_2"] = D2TP_T2TP( pip, pi0, pim, ks0, 2 );
475 spinfactor["D2TP_T2TP_pippimpi0_pippim_2"] = D2TP_T2TP( pip, pim, pi0, ks0, 2 ); //
476 spinfactor["D2TP_T2TP_ks0pimpi0_ks0pim_2"] = D2TP_T2TP( ks0, pim, pi0, pip, 2 );
477 spinfactor["D2TP_T2TP_ks0pimpi0_ks0pi0_2"] = D2TP_T2TP( ks0, pi0, pim, pip, 2 );
478 spinfactor["D2TP_T2TP_ks0pimpi0_pimpi0_2"] = D2TP_T2TP( pim, pi0, ks0, pip, 2 ); //
479 spinfactor["D2TP_T2TP_ks0pippi0_ks0pip_2"] = D2TP_T2TP( ks0, pip, pi0, pim, 2 );
480 spinfactor["D2TP_T2TP_ks0pippi0_ks0pi0_2"] = D2TP_T2TP( ks0, pi0, pip, pim, 2 );
481 spinfactor["D2TP_T2TP_ks0pippi0_pippi0_2"] = D2TP_T2TP( pip, pi0, ks0, pim, 2 ); //
482 spinfactor["D2TP_T2TP_ks0pippim_ks0pip_2"] = D2TP_T2TP( ks0, pip, pim, pi0, 2 );
483 spinfactor["D2TP_T2TP_ks0pippim_ks0pim_2"] = D2TP_T2TP( ks0, pim, pip, pi0, 2 );
484 spinfactor["D2TP_T2TP_ks0pippim_pippim_2"] = D2TP_T2TP( pip, pim, ks0, pi0, 2 ); //
485
486 //---- D2PTP_PT2SP
487 spinfactor["D2PTP_PT2SP_pippimpi0_pimpi0_2"] =
488 D2PTP_PT2SP( pi0, pim, pip, ks0, 2 ); // sequence of pi0 pi-
489 spinfactor["D2PTP_PT2SP_pippimpi0_pippi0_2"] = D2PTP_PT2SP( pip, pi0, pim, ks0, 2 );
490 spinfactor["D2PTP_PT2SP_pippimpi0_pippim_2"] = D2PTP_PT2SP( pip, pim, pi0, ks0, 2 ); //
491 spinfactor["D2PTP_PT2SP_ks0pimpi0_ks0pim_2"] = D2PTP_PT2SP( ks0, pim, pi0, pip, 2 );
492 spinfactor["D2PTP_PT2SP_ks0pimpi0_ks0pi0_2"] = D2PTP_PT2SP( ks0, pi0, pim, pip, 2 );
493 spinfactor["D2PTP_PT2SP_ks0pimpi0_pimpi0_2"] = D2PTP_PT2SP( pim, pi0, ks0, pip, 2 ); //
494 spinfactor["D2PTP_PT2SP_ks0pippi0_ks0pip_2"] = D2PTP_PT2SP( ks0, pip, pi0, pim, 2 );
495 spinfactor["D2PTP_PT2SP_ks0pippi0_ks0pi0_2"] = D2PTP_PT2SP( ks0, pi0, pip, pim, 2 );
496 spinfactor["D2PTP_PT2SP_ks0pippi0_pippi0_2"] = D2PTP_PT2SP( pip, pi0, ks0, pim, 2 ); //
497 spinfactor["D2PTP_PT2SP_ks0pippim_ks0pip_2"] = D2PTP_PT2SP( ks0, pip, pim, pi0, 2 );
498 spinfactor["D2PTP_PT2SP_ks0pippim_ks0pim_2"] = D2PTP_PT2SP( ks0, pim, pip, pi0, 2 );
499 spinfactor["D2PTP_PT2SP_ks0pippim_pippim_2"] = D2PTP_PT2SP( pip, pim, ks0, pi0, 2 ); //
500
501 //---- D2PTP_PT2VP
502 spinfactor["D2PTP_PT2VP_pippimpi0_pimpi0_2"] =
503 D2PTP_PT2VP( pi0, pim, pip, ks0, 2 ); // sequence of pi0 pi-
504 spinfactor["D2PTP_PT2VP_pippimpi0_pippi0_2"] = D2PTP_PT2VP( pip, pi0, pim, ks0, 2 );
505 spinfactor["D2PTP_PT2VP_pippimpi0_pippim_2"] = D2PTP_PT2VP( pip, pim, pi0, ks0, 2 ); //
506 spinfactor["D2PTP_PT2VP_ks0pimpi0_ks0pim_2"] = D2PTP_PT2VP( ks0, pim, pi0, pip, 2 );
507 spinfactor["D2PTP_PT2VP_ks0pimpi0_ks0pi0_2"] = D2PTP_PT2VP( ks0, pi0, pim, pip, 2 );
508 spinfactor["D2PTP_PT2VP_ks0pimpi0_pimpi0_2"] = D2PTP_PT2VP( pim, pi0, ks0, pip, 2 ); //
509 spinfactor["D2PTP_PT2VP_ks0pippi0_ks0pip_2"] = D2PTP_PT2VP( ks0, pip, pi0, pim, 2 );
510 spinfactor["D2PTP_PT2VP_ks0pippi0_ks0pi0_2"] = D2PTP_PT2VP( ks0, pi0, pip, pim, 2 );
511 spinfactor["D2PTP_PT2VP_ks0pippi0_pippi0_2"] = D2PTP_PT2VP( pip, pi0, ks0, pim, 2 ); //
512 spinfactor["D2PTP_PT2VP_ks0pippim_ks0pip_2"] = D2PTP_PT2VP( ks0, pip, pim, pi0, 2 );
513 spinfactor["D2PTP_PT2VP_ks0pippim_ks0pim_2"] = D2PTP_PT2VP( ks0, pim, pip, pi0, 2 );
514 spinfactor["D2PTP_PT2VP_ks0pippim_pippim_2"] = D2PTP_PT2VP( pip, pim, ks0, pi0, 2 ); //
515
516 //---- D2PTP_PT2TP
517 spinfactor["D2PTP_PT2TP_pippimpi0_pimpi0_2"] =
518 D2PTP_PT2TP( pi0, pim, pip, ks0, 2 ); // sequence of pi0 pi-
519 spinfactor["D2PTP_PT2TP_pippimpi0_pippi0_2"] = D2PTP_PT2TP( pip, pi0, pim, ks0, 2 );
520 spinfactor["D2PTP_PT2TP_pippimpi0_pippim_2"] = D2PTP_PT2TP( pip, pim, pi0, ks0, 2 ); //
521 spinfactor["D2PTP_PT2TP_ks0pimpi0_ks0pim_2"] = D2PTP_PT2TP( ks0, pim, pi0, pip, 2 );
522 spinfactor["D2PTP_PT2TP_ks0pimpi0_ks0pi0_2"] = D2PTP_PT2TP( ks0, pi0, pim, pip, 2 );
523 spinfactor["D2PTP_PT2TP_ks0pimpi0_pimpi0_2"] = D2PTP_PT2TP( pim, pi0, ks0, pip, 2 ); //
524 spinfactor["D2PTP_PT2TP_ks0pippi0_ks0pip_2"] = D2PTP_PT2TP( ks0, pip, pi0, pim, 2 );
525 spinfactor["D2PTP_PT2TP_ks0pippi0_ks0pi0_2"] = D2PTP_PT2TP( ks0, pi0, pip, pim, 2 );
526 spinfactor["D2PTP_PT2TP_ks0pippi0_pippi0_2"] = D2PTP_PT2TP( pip, pi0, ks0, pim, 2 ); //
527 spinfactor["D2PTP_PT2TP_ks0pippim_ks0pip_2"] = D2PTP_PT2TP( ks0, pip, pim, pi0, 2 );
528 spinfactor["D2PTP_PT2TP_ks0pippim_ks0pim_2"] = D2PTP_PT2TP( ks0, pim, pip, pi0, 2 );
529 spinfactor["D2PTP_PT2TP_ks0pippim_pippim_2"] = D2PTP_PT2TP( pip, pim, ks0, pi0, 2 ); //
530
531 //------------------------------------------------------------//
532 // Quasi Two-Body Decay //
533 //------------------------------------------------------------//
534 //---- D2SS
535 spinfactor["D2SS_ks0pim_pippi0_0"] = D2SS();
536 spinfactor["D2SS_ks0pi0_pippim_0"] = D2SS();
537 spinfactor["D2SS_ks0pip_pimpi0_0"] = D2SS();
538
539 //---- D2VS
540 spinfactor["D2VS_ks0pim_pippi0_1"] = D2VS( ks0, pim, pip, pi0, 1 );
541 spinfactor["D2VS_ks0pi0_pippim_1"] = D2VS( ks0, pi0, pip, pim, 1 );
542 spinfactor["D2VS_ks0pip_pimpi0_1"] = D2VS( ks0, pip, pim, pi0, 1 );
543 spinfactor["D2VS_pimpi0_ks0pip_1"] = D2VS( pim, pi0, ks0, pip, 1 );
544 spinfactor["D2VS_pippim_ks0pi0_1"] = D2VS( pip, pim, ks0, pi0, 1 );
545 spinfactor["D2VS_pippi0_ks0pim_1"] = D2VS( pip, pi0, ks0, pim, 1 );
546
547 //---- D2VV, [S,P,D]
548 spinfactor["D2VV_ks0pim_pippi0_0"] = D2VV( ks0, pim, pip, pi0, 0 );
549 spinfactor["D2VV_ks0pi0_pippim_0"] = D2VV( ks0, pi0, pip, pim, 0 );
550 spinfactor["D2VV_ks0pip_pimpi0_0"] = D2VV( ks0, pip, pim, pi0, 0 );
551 spinfactor["D2VV_ks0pim_pippi0_1"] = D2VV( ks0, pim, pip, pi0, 1 );
552 spinfactor["D2VV_ks0pi0_pippim_1"] = D2VV( ks0, pi0, pip, pim, 1 );
553 spinfactor["D2VV_ks0pip_pimpi0_1"] = D2VV( ks0, pip, pim, pi0, 1 );
554 spinfactor["D2VV_ks0pim_pippi0_2"] = D2VV( ks0, pim, pip, pi0, 2 );
555 spinfactor["D2VV_ks0pi0_pippim_2"] = D2VV( ks0, pi0, pip, pim, 2 );
556 spinfactor["D2VV_ks0pip_pimpi0_2"] = D2VV( ks0, pip, pim, pi0, 2 );
557
558 //---- D2TS
559 spinfactor["D2TS_ks0pim_pippi0_2"] = D2TS( ks0, pim, pip, pi0, 2 );
560 spinfactor["D2TS_ks0pi0_pippim_2"] = D2TS( ks0, pi0, pip, pim, 2 );
561 spinfactor["D2TS_ks0pip_pimpi0_2"] = D2TS( ks0, pip, pim, pi0, 2 );
562 spinfactor["D2TS_pimpi0_ks0pip_2"] = D2TS( pim, pi0, ks0, pip, 2 );
563 spinfactor["D2TS_pippim_ks0pi0_2"] = D2TS( pip, pim, ks0, pi0, 2 );
564 spinfactor["D2TS_pippi0_ks0pim_2"] = D2TS( pip, pi0, ks0, pim, 2 );
565
566 //---- D2TV, [P, D]
567 spinfactor["D2TV_ks0pim_pippi0_1"] = D2TV( ks0, pim, pip, pi0, 1 );
568 spinfactor["D2TV_ks0pi0_pippim_1"] = D2TV( ks0, pi0, pip, pim, 1 );
569 spinfactor["D2TV_ks0pip_pimpi0_1"] = D2TV( ks0, pip, pim, pi0, 1 );
570 spinfactor["D2TV_pimpi0_ks0pip_1"] = D2TV( pim, pi0, ks0, pip, 1 );
571 spinfactor["D2TV_pippim_ks0pi0_1"] = D2TV( pip, pim, ks0, pi0, 1 );
572 spinfactor["D2TV_pippi0_ks0pim_1"] = D2TV( pip, pi0, ks0, pim, 1 );
573 spinfactor["D2TV_ks0pim_pippi0_2"] = D2TV( ks0, pim, pip, pi0, 2 );
574 spinfactor["D2TV_ks0pi0_pippim_2"] = D2TV( ks0, pi0, pip, pim, 2 );
575 spinfactor["D2TV_ks0pip_pimpi0_2"] = D2TV( ks0, pip, pim, pi0, 2 );
576 spinfactor["D2TV_pimpi0_ks0pip_2"] = D2TV( pim, pi0, ks0, pip, 2 );
577 spinfactor["D2TV_pippim_ks0pi0_2"] = D2TV( pip, pim, ks0, pi0, 2 );
578 spinfactor["D2TV_pippi0_ks0pim_2"] = D2TV( pip, pi0, ks0, pim, 2 );
579
580 //---- D2TT, [S,P,D]
581 spinfactor["D2TT_ks0pim_pippi0_0"] = D2TT( ks0, pim, pip, pi0, 0 );
582 spinfactor["D2TT_ks0pi0_pippim_0"] = D2TT( ks0, pi0, pip, pim, 0 );
583 spinfactor["D2TT_ks0pip_pimpi0_0"] = D2TT( ks0, pip, pim, pi0, 0 );
584 spinfactor["D2TT_ks0pim_pippi0_1"] = D2TT( ks0, pim, pip, pi0, 1 );
585 spinfactor["D2TT_ks0pi0_pippim_1"] = D2TT( ks0, pi0, pip, pim, 1 );
586 spinfactor["D2TT_ks0pip_pimpi0_1"] = D2TT( ks0, pip, pim, pi0, 1 );
587 spinfactor["D2TT_ks0pim_pippi0_2"] = D2TT( ks0, pim, pip, pi0, 2 );
588 spinfactor["D2TT_ks0pi0_pippim_2"] = D2TT( ks0, pi0, pip, pim, 2 );
589 spinfactor["D2TT_ks0pip_pimpi0_2"] = D2TT( ks0, pip, pim, pi0, 2 );
590}
591
592//============================================================//
593// Propagator //
594//============================================================//
595void EvtD0ToKSpipipi0::createPropagator( EvtVector4R Ks0, EvtVector4R Pip, EvtVector4R Pim,
596 EvtVector4R Pi0 ) {
597
598 double m2_ks0 = Ks0.mass2();
599 double m2_pip = Pip.mass2();
600 double m2_pim = Pim.mass2();
601 double m2_pi0 = Pi0.mass2();
602 double m2_ks0pip = ( Ks0 + Pip ).mass2();
603 double m2_ks0pim = ( Ks0 + Pim ).mass2();
604 double m2_ks0pi0 = ( Ks0 + Pi0 ).mass2();
605 double m2_pippim = ( Pip + Pim ).mass2();
606 double m2_pippi0 = ( Pip + Pi0 ).mass2();
607 double m2_pimpi0 = ( Pim + Pi0 ).mass2();
608 double m2_ks0pippim = ( Ks0 + Pip + Pim ).mass2();
609 double m2_ks0pippi0 = ( Ks0 + Pip + Pi0 ).mass2();
610 double m2_ks0pimpi0 = ( Ks0 + Pim + Pi0 ).mass2();
611 double m2_pippimpi0 = ( Pip + Pim + Pi0 ).mass2();
612
613 std::vector<double> ks0;
614 ks0.clear();
615 ks0.push_back( Ks0.get( 1 ) );
616 ks0.push_back( Ks0.get( 2 ) );
617 ks0.push_back( Ks0.get( 3 ) );
618 ks0.push_back( Ks0.get( 0 ) );
619 std::vector<double> pip;
620 pip.clear();
621 pip.push_back( Pip.get( 1 ) );
622 pip.push_back( Pip.get( 2 ) );
623 pip.push_back( Pip.get( 3 ) );
624 pip.push_back( Pip.get( 0 ) );
625 std::vector<double> pim;
626 pim.clear();
627 pim.push_back( Pim.get( 1 ) );
628 pim.push_back( Pim.get( 2 ) );
629 pim.push_back( Pim.get( 3 ) );
630 pim.push_back( Pim.get( 0 ) );
631 std::vector<double> pi0;
632 pi0.clear();
633 pi0.push_back( Pi0.get( 1 ) );
634 pi0.push_back( Pi0.get( 2 ) );
635 pi0.push_back( Pi0.get( 3 ) );
636 pi0.push_back( Pi0.get( 0 ) );
637
638 propagator["kstarm_892"] =
639 create_RBW_propagator( "kstarm_892", m2_ks0pim, m2_ks0, m2_pim, 1 );
640 propagator["kstar0_892"] =
641 create_RBW_propagator( "kstar0_892", m2_ks0pi0, m2_ks0, m2_pi0, 1 );
642 propagator["kstarp_892"] =
643 create_RBW_propagator( "kstarp_892", m2_ks0pip, m2_ks0, m2_pip, 1 );
644 propagator["k1m_1270"] = create_BW_propagator( "k1m_1270", m2_ks0pimpi0 );
645 propagator["k10_1270"] = create_BW_propagator( "k10_1270", m2_ks0pippim );
646 propagator["k1m_1400"] = create_BW_propagator( "k1m_1400", m2_ks0pimpi0 );
647 propagator["k10_1400"] = create_BW_propagator( "k10_1400", m2_ks0pippim );
648 propagator["kstarm2_1410"] = create_BW_propagator( "kstarm_1410", m2_ks0pim );
649 propagator["kstar02_1410"] = create_BW_propagator( "kstar0_1410", m2_ks0pi0 );
650 propagator["kstarm3_1410"] = create_BW_propagator( "kstarm_1410", m2_ks0pimpi0 );
651 propagator["kstar03_1410"] = create_BW_propagator( "kstar0_1410", m2_ks0pippim );
652 propagator["k0starm_1430"] =
653 create_KPiSLASS_propagator( "k0starm_1430", m2_ks0pim, m2_ks0, m2_pim );
654 propagator["k0star0_1430"] =
655 create_KPiSLASS_propagator( "k0star0_1430", m2_ks0pi0, m2_ks0, m2_pi0 );
656 // propagator["k0starm_1430" ] = create_BW_propagator("k0starm_1430", m2_ks0pim);
657 // propagator["k0star0_1430" ] = create_BW_propagator("k0star0_1430", m2_ks0pi0);
658 propagator["k2starm2_1430"] = create_BW_propagator( "k2starm_1430", m2_ks0pim );
659 propagator["k2star02_1430"] = create_BW_propagator( "k2star0_1430", m2_ks0pi0 );
660 propagator["k2starm3_1430"] = create_BW_propagator( "k2starm_1430", m2_ks0pimpi0 );
661 propagator["k2star03_1430"] = create_BW_propagator( "k2star0_1430", m2_ks0pippim );
662 propagator["km_1460"] = create_BW_propagator( "km_1460", m2_ks0pimpi0 );
663 propagator["k0_1460"] = create_BW_propagator( "k0_1460", m2_ks0pippim );
664 propagator["k2m_1580"] = create_BW_propagator( "k2m_1580", m2_ks0pimpi0 );
665 propagator["k20_1580"] = create_BW_propagator( "k20_1580", m2_ks0pippim );
666 propagator["km_1630"] = create_BW_propagator( "km_1630", m2_ks0pimpi0 );
667 propagator["k0_1630"] = create_BW_propagator( "k0_1630", m2_ks0pippim );
668 propagator["k1m_1650"] = create_BW_propagator( "k1m_1650", m2_ks0pimpi0 );
669 propagator["k10_1650"] = create_BW_propagator( "k10_1650", m2_ks0pippim );
670 propagator["kstarm2_1680"] = create_BW_propagator( "kstarm_1680", m2_ks0pim );
671 propagator["kstar02_1680"] = create_BW_propagator( "kstar0_1680", m2_ks0pi0 );
672 propagator["kstarm3_1680"] = create_BW_propagator( "kstarm_1680", m2_ks0pimpi0 );
673 propagator["kstar03_1680"] = create_BW_propagator( "kstar0_1680", m2_ks0pippim );
674 propagator["k2m_1770"] = create_BW_propagator( "k2m_1770", m2_ks0pimpi0 );
675 propagator["k20_1770"] = create_BW_propagator( "k20_1770", m2_ks0pippim );
676
677 propagator["sigma_500"] = create_sigma_propagator( "sigma_500", m2_pippim );
678 propagator["rhop_770"] = create_RBW_propagator( "rhop_770", m2_pippi0, m2_pip, m2_pi0, 1 );
679 propagator["rhom_770"] = create_RBW_propagator( "rhom_770", m2_pimpi0, m2_pim, m2_pi0, 1 );
680 propagator["rho0_770"] = create_GS_propagator( "rho0_770", m2_pippim );
681 propagator["omega2_782"] = create_BW_propagator( "omega_782", m2_pippim );
682 propagator["f0_980"] =
683 create_Flatte2_propagator( "f0_980", m2_pippim, mpion, mpion, mkaon, mkaon );
684 propagator["f0_1370"] = create_BW_propagator( "f0_1370", m2_pippim );
685 propagator["f2_1270"] = create_BW_propagator( "f2_1270", m2_pippim );
686
687 // // The resolution should be considered in the fit, but when calculating the physical
688 // values, we could only use the pure itude
689 // //propagator["phi_1020" ] = new GPUPropagatorBreitWigner("phi_1021", m2_pippimpi0,
690 // mean, sigma);
691 propagator["phi_1020"] = create_BW_propagator( "phi_1020", m2_pippimpi0 );
692 // //propagator["omega3_782" ] = new GPUPropagatorBreitWigner("omega_782", m2_pippimpi0,
693 // mean, sigma);
694 propagator["omega3_782"] = create_BW_propagator( "omega_782", m2_pippimpi0 );
695 // //propagator["eta_547" ] = new GPUPropagatorBreitWigner("eta_547", m2_pippimpi0,
696 // mean, sigma);
697 propagator["eta_547"] = create_BW_propagator( "eta_547", m2_pippimpi0 );
698
699 // propagator["a10_1260_rhom_770_0"] = create_RBW_propagator("a10_1260", m2_pippimpi0,
700 // m2_pimpi0, m2_pip, 0); propagator["a10_1260_rhop_770_0"] =
701 // create_RBW_propagator("a10_1260", m2_pippimpi0, m2_pippi0, m2_pim, 0);
702 // propagator["a10_1260_rhom_770_2"] = create_RBW_propagator("a10_1260", m2_pippimpi0,
703 // m2_pimpi0, m2_pip, 2); propagator["a10_1260_rhop_770_2"] =
704 // create_RBW_propagator("a10_1260", m2_pippimpi0, m2_pippi0, m2_pim, 2);
705 //
706 propagator["kstarm_phsp"] = create_BW_propagator( "kstarm_phsp", m2_ks0pim );
707 propagator["kstar0_phsp"] = create_BW_propagator( "kstar0_phsp", m2_ks0pi0 );
708 propagator["kstarp_phsp"] = create_BW_propagator( "kstarp_phsp", m2_ks0pip );
709 propagator["rhop_phsp"] = create_BW_propagator( "rhop_phsp", m2_pippi0 );
710 propagator["rhom_phsp"] = create_BW_propagator( "rhom_phsp", m2_pimpi0 );
711 propagator["rho0_phsp"] = create_BW_propagator( "rho0_phsp", m2_pippim );
712 propagator["k1m_phsp"] = create_BW_propagator( "k1m_phsp", m2_ks0pimpi0 );
713 propagator["k10_phsp"] = create_BW_propagator( "k10_phsp", m2_ks0pippim );
714 propagator["k1p_phsp"] = create_BW_propagator( "k1p_phsp", m2_ks0pippi0 );
715 propagator["a10_phsp"] = create_BW_propagator( "a10_phsp", m2_pippimpi0 );
716}
717
718//============================================================//
719// Tensor Contract Function //
720//============================================================//
721// Four-momentum (PX, PY, PZ; E)
722double EvtD0ToKSpipipi0::contract_11_0( vector<double> pa, vector<double> pb ) {
723 assert( pa.size() == pb.size() && pa.size() == 4 );
724
725 double temp = pa[3] * pb[3] - pa[0] * pb[0] - pa[1] * pb[1] - pa[2] * pb[2];
726 return temp;
727}
728
729vector<double> EvtD0ToKSpipipi0::contract_21_1( vector<double> pa, vector<double> pb ) {
730 assert( pa.size() == 16 && pb.size() == 4 );
731
732 vector<double> temp;
733 temp.clear();
734 for ( int i = 0; i < 4; i++ )
735 {
736 double sum = 0;
737 for ( int j = 0; j < 4; j++ )
738 {
739 int idx = i * 4 + j;
740 sum += pa[idx] * pb[j] * g[4 * j + j];
741 }
742 temp.push_back( sum );
743 }
744 return temp;
745}
746
747double EvtD0ToKSpipipi0::contract_22_0( vector<double> pa, vector<double> pb ) {
748 assert( pa.size() == 16 && pb.size() == 16 );
749
750 double temp = 0;
751 for ( int i = 0; i < 4; i++ )
752 {
753 for ( int j = 0; j < 4; j++ )
754 {
755 int idx = i * 4 + j;
756 temp += pa[idx] * pb[idx] * g[4 * i + i] * g[4 * j + j];
757 }
758 }
759 return temp;
760}
761
762vector<double> EvtD0ToKSpipipi0::contract_31_2( vector<double> pa, vector<double> pb ) {
763 assert( pa.size() == 64 && pb.size() == 4 );
764
765 vector<double> temp;
766 temp.clear();
767 for ( int i = 0; i < 16; i++ )
768 {
769 double sum = 0;
770 for ( int j = 0; j < 4; j++ )
771 {
772 int idx = i * 4 + j;
773 sum += pa[idx] * pb[j] * g[4 * j + j];
774 }
775 temp.push_back( sum );
776 }
777 return temp;
778}
779
780vector<double> EvtD0ToKSpipipi0::contract_41_3( vector<double> pa, vector<double> pb ) {
781 assert( pa.size() == 256 && pb.size() == 4 );
782
783 vector<double> temp;
784 temp.clear();
785 for ( int i = 0; i < 64; i++ )
786 {
787 double sum = 0;
788 for ( int j = 0; j < 4; j++ )
789 {
790 int idx = i * 4 + j;
791 sum += pa[idx] * pb[j] * g[4 * j + j];
792 }
793 temp.push_back( sum );
794 }
795 return temp;
796}
797
798vector<double> EvtD0ToKSpipipi0::contract_42_2( vector<double> pa, vector<double> pb ) {
799 assert( pa.size() == 256 && pb.size() == 16 );
800
801 vector<double> temp;
802 temp.clear();
803 for ( int i = 0; i < 16; i++ )
804 {
805 double sum = 0;
806 for ( int j = 0; j < 4; j++ )
807 {
808 for ( int k = 0; k < 4; k++ )
809 {
810 int idxa = i * 16 + j * 4 + k;
811 int idxb = j * 4 + k;
812 sum += pa[idxa] * pb[idxb] * g[4 * j + j] * g[4 * k + k];
813 }
814 }
815 temp.push_back( sum );
816 }
817 return temp;
818}
819
820vector<double> EvtD0ToKSpipipi0::contract_22_2( vector<double> pa, vector<double> pb ) {
821 assert( pa.size() == 16 && pb.size() == 16 );
822
823 vector<double> temp;
824 temp.clear();
825 for ( int i = 0; i < 4; i++ )
826 {
827 for ( int j = 0; j < 4; j++ )
828 {
829 double sum = 0;
830 for ( int k = 0; k < 4; k++ )
831 {
832 int idxa = i * 4 + k;
833 int idxb = j * 4 + k;
834 sum += pa[idxa] * pb[idxb] * g[4 * k + k];
835 }
836 temp.push_back( sum );
837 }
838 }
839 return temp;
840}
841
842//------------------------------------------------------------//
843// Spin Factor //
844//------------------------------------------------------------//
845vector<double> EvtD0ToKSpipipi0::Proj( vector<double> pa, int rank ) {
846 if ( pa.size() != 4 )
847 {
848 cout << "Error: pa" << endl;
849 exit( 0 );
850 }
851 if ( rank < 0 )
852 {
853 cout << "Error: L<0 !!!" << endl;
854 exit( 0 );
855 }
856
857 double msa = contract_11_0( pa, pa );
858
859 // Projection Operator 2-rank
860 vector<double> proj_uv;
861 proj_uv.clear();
862 for ( int i = 0; i < 4; i++ )
863 {
864 for ( int j = 0; j < 4; j++ )
865 {
866 int idx = i * 4 + j;
867 double temp = -g[idx] + pa[i] * pa[j] / msa;
868 proj_uv.push_back( temp );
869 }
870 }
871
872 vector<double> proj;
873 proj.clear();
874 // Orbital Tensors
875 if ( rank == 0 )
876 {
877 vector<double> t;
878 t.clear();
879 t.push_back( 1.0 );
880 proj = t;
881 }
882 else if ( rank == 1 ) { proj = proj_uv; }
883 else if ( rank == 2 )
884 {
885 vector<double> proj_uvmn;
886 proj_uvmn.clear();
887 for ( int i = 0; i < 4; i++ )
888 {
889 for ( int j = 0; j < 4; j++ )
890 {
891 for ( int k = 0; k < 4; k++ )
892 {
893 for ( int l = 0; l < 4; l++ )
894 {
895
896 int idx1_1 = 4 * i + k;
897 int idx1_2 = 4 * i + l;
898 int idx1_3 = 4 * i + j;
899
900 int idx2_1 = 4 * j + l;
901 int idx2_2 = 4 * j + k;
902 int idx2_3 = 4 * k + l;
903
904 double temp = ( 1.0 / 2.0 ) * ( proj_uv[idx1_1] * proj_uv[idx2_1] +
905 proj_uv[idx1_2] * proj_uv[idx2_2] ) -
906 ( 1.0 / 3.0 ) * proj_uv[idx1_3] * proj_uv[idx2_3];
907 proj_uvmn.push_back( temp );
908 }
909 }
910 }
911 }
912 proj = proj_uvmn;
913 }
914 else
915 {
916 cout << "rank>2: please add it by yourself!!!" << endl;
917 exit( 0 );
918 }
919
920 return proj;
921}
922
923// Orbital Tensor
924vector<double> EvtD0ToKSpipipi0::OrbitalTensors( vector<double> pa, vector<double> pb,
925 vector<double> pc, double r, int rank ) {
926 assert( pa.size() == 4 && pb.size() == 4 && pc.size() == 4 );
927 if ( rank < 0 )
928 {
929 cout << "Error: L<0 !!!" << endl;
930 exit( 0 );
931 }
932
933 // relative momentum
934 vector<double> mr;
935 mr.clear();
936
937 for ( int i = 0; i < 4; i++ )
938 {
939 double temp = pb[i] - pc[i];
940 mr.push_back( temp );
941 }
942
943 // "Masses" of particles
944 double msa = contract_11_0( pa, pa );
945 double msb = contract_11_0( pb, pb );
946 double msc = contract_11_0( pc, pc );
947
948 // Q^2_abc
949 double top = msa + msb - msc;
950 double Q2abc = top * top / ( 4.0 * msa ) - msb;
951
952 // Barrier factors
953 // double Q_0 = 0.197321f/r; //mRadius with unit of fm;
954 double Q_0 = 1.0 / r; // mRadius with unit of GeV^-1;
955 double Q_02 = Q_0 * Q_0;
956 double Q_04 = Q_02 * Q_02;
957 // double Q_06 = Q_04*Q_02;
958 // double Q_08 = Q_04*Q_04;
959
960 double Q4abc = Q2abc * Q2abc;
961 // double Q6abc = Q4abc*Q2abc;
962 // double Q8abc = Q4abc*Q4abc;
963
964 double mB1 = sqrt( 2.0 / ( Q2abc + Q_02 ) );
965 double mB2 = sqrt( 13.0 / ( Q4abc + ( 3.0 * Q_02 ) * Q2abc + 9.0 * Q_04 ) );
966 // mB3 = &sqrt(277.0f/(Q6abc + (6.0f*Q_02)*Q4abc + (45.0f*Q_04)*Q2abc + 225.0f*Q_06));
967 // mB4 = &sqrt(12746.0f/(Q8abc + (10.0f*Q_02)*Q6abc + (135.0f*Q_04)*Q4abc +
968 // (1575.0f*Q_06)*Q2abc + 11025.0f*Q_08));
969
970 // Projection Operator 2-rank
971 vector<double> proj_uv;
972 proj_uv.clear();
973 for ( int i = 0; i < 4; i++ )
974 {
975 for ( int j = 0; j < 4; j++ )
976 {
977 int idx = i * 4 + j;
978 double temp = -g[idx] + pa[i] * pa[j] / msa;
979 proj_uv.push_back( temp );
980 }
981 }
982
983 vector<double> Bt;
984 Bt.clear();
985 // Orbital Tensors
986 if ( rank == 0 )
987 {
988 vector<double> t;
989 t.clear();
990 t.push_back( 1.0 );
991 Bt = t;
992 }
993 else if ( rank < 3 )
994 {
995 vector<double> t_u;
996 t_u.clear();
997 vector<double> Bt_u;
998 Bt_u.clear();
999 for ( int i = 0; i < 4; i++ )
1000 {
1001 double temp = 0;
1002 for ( int j = 0; j < 4; j++ )
1003 {
1004 int idx = i * 4 + j;
1005 temp += -proj_uv[idx] * mr[j] * g[j * 4 + j];
1006 }
1007 t_u.push_back( temp );
1008 Bt_u.push_back( temp * mB1 );
1009 }
1010 if ( rank == 1 ) Bt = Bt_u;
1011
1012 double t_u2 = contract_11_0( t_u, t_u );
1013
1014 vector<double> Bt_uv;
1015 Bt_uv.clear();
1016 for ( int i = 0; i < 4; i++ )
1017 {
1018 for ( int j = 0; j < 4; j++ )
1019 {
1020 int idx = 4 * i + j;
1021 double temp = t_u[i] * t_u[j] + ( 1.0 / 3.0 ) * proj_uv[idx] * t_u2;
1022 Bt_uv.push_back( temp * mB2 );
1023 }
1024 }
1025 if ( rank == 2 ) Bt = Bt_uv;
1026 }
1027 else
1028 {
1029 cout << "rank>2: please add it by yorself!!!" << endl;
1030 exit( 0 );
1031 }
1032
1033 return Bt;
1034}
1035
1036// The sequence is 0+, 0-, 1-, 1+, 2+, 2-
1037//------------------------------------------------------------//
1038// Cascade Decay //
1039//------------------------------------------------------------//
1040//---- D -> P P1 [S], P -> S P2 [S]
1041double EvtD0ToKSpipipi0::D2PP_P2SP() {
1042
1043 double SF_D2PP_P2SP = 1.0;
1044
1045 return SF_D2PP_P2SP;
1046}
1047
1048//----D -> P P1 [S], P -> V P2 [P], V -> P3 P4 [P]
1049double EvtD0ToKSpipipi0::D2PP_P2VP( vector<double> p1, vector<double> p2, vector<double> p3,
1050 vector<double> p4, int l ) {
1051 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1052 assert( l >= 0 );
1053
1054 vector<double> V = p1 + p2;
1055 vector<double> P2 = p3;
1056 vector<double> P = V + P2;
1057 vector<double> P1 = p4;
1058
1059 vector<double> P_orbitals_Spin1OrbitalTensor = OrbitalTensors( P, V, P2, rRes, 1 );
1060 vector<double> V_orbitals_Spin1OrbitalTensor = OrbitalTensors( V, p1, p2, rRes, 1 );
1061
1062 double SF_D2PP_P2VP = 0.;
1063 SF_D2PP_P2VP = contract_11_0( P_orbitals_Spin1OrbitalTensor, V_orbitals_Spin1OrbitalTensor );
1064 return SF_D2PP_P2VP;
1065}
1066
1067//----D -> V1 P1 [P], V1 -> V2 P2 [P], V2 -> P3 P4 [P]
1068double EvtD0ToKSpipipi0::D2VP_V2VP( vector<double> p1, vector<double> p2, vector<double> p3,
1069 vector<double> p4, int l ) {
1070 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1071 assert( l >= 0 );
1072
1073 vector<double> V2 = p1 + p2;
1074 vector<double> P2 = p3;
1075 vector<double> V1 = V2 + P2;
1076 vector<double> P1 = p4;
1077 vector<double> D = V1 + P1;
1078
1079 vector<double> D_orbitals_Spin1OrbitalTensor = OrbitalTensors( D, V1, P1, rD, 1 );
1080 vector<double> V1_orbitals_Spin1OrbitalTensor = OrbitalTensors( V1, V2, P2, rRes, 1 );
1081 vector<double> V2_orbitals_Spin1OrbitalTensor = OrbitalTensors( V2, p1, p2, rRes, 1 );
1082 vector<double> Proj1_V1 = Proj( V1, 1 );
1083
1084 double SF_D2VP_V2VP = 0.;
1085 SF_D2VP_V2VP = contract_11_0( contract_21_1( contract_31_2( contract_41_3( epsilon, V1 ),
1086 V2_orbitals_Spin1OrbitalTensor ),
1087 V1_orbitals_Spin1OrbitalTensor ),
1088 contract_21_1( Proj1_V1, D_orbitals_Spin1OrbitalTensor ) );
1089 // SF_D2VP_V2VP = contract_11_0( contract_21_1(contract_31_2(contract_41_3(epsilon, V1),
1090 // V1_orbitals_Spin1OrbitalTensor) | contract_21_1(Proj1_V1, V2_orbitals_Spin1OrbitalTensor))
1091 // | contract_21_1(Proj1_V1, D_orbitals_Spin1OrbitalTensor) );//alternative
1092 return SF_D2VP_V2VP;
1093}
1094
1095//----D -> A P1 [P], A -> S P2 [P], S -> P3 P4 [S]
1096double EvtD0ToKSpipipi0::D2AP_A2SP( vector<double> p1, vector<double> p2, vector<double> p3,
1097 vector<double> p4, int l ) {
1098 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1099 assert( l >= 0 );
1100
1101 vector<double> S = p1 + p2;
1102 vector<double> P2 = p3;
1103 vector<double> A = S + P2;
1104 vector<double> P1 = p4;
1105 vector<double> D = A + P1;
1106
1107 vector<double> D_orbitals_Spin1OrbitalTensor = OrbitalTensors( D, A, P1, rD, 1 );
1108 vector<double> A_orbitals_Spin1OrbitalTensor = OrbitalTensors( A, S, P2, rRes, 1 );
1109
1110 double SF_D2AP_A2SP = 0.;
1111 SF_D2AP_A2SP = contract_11_0( D_orbitals_Spin1OrbitalTensor, A_orbitals_Spin1OrbitalTensor );
1112 return SF_D2AP_A2SP;
1113}
1114
1115//----D -> A P1 [P], A -> V P2 [S, D], V -> P3 P4 [P]
1116double EvtD0ToKSpipipi0::D2AP_A2VP( vector<double> p1, vector<double> p2, vector<double> p3,
1117 vector<double> p4, int l ) {
1118 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1119 assert( l >= 0 );
1120
1121 vector<double> V = p1 + p2;
1122 vector<double> P2 = p3;
1123 vector<double> A = V + P2;
1124 vector<double> P1 = p4;
1125 vector<double> D = A + P1;
1126
1127 vector<double> D_orbitals_Spin1OrbitalTensor = OrbitalTensors( D, A, P1, rD, 1 );
1128 vector<double> A_orbitals_Spin2OrbitalTensor = OrbitalTensors( A, V, P2, rRes, 2 );
1129 vector<double> V_orbitals_Spin1OrbitalTensor = OrbitalTensors( V, p1, p2, rRes, 1 );
1130 vector<double> Proj1_A = Proj( A, 1 );
1131
1132 double SF_D2AP_A2VP = 0.;
1133 if ( l == 0 )
1134 {
1135 SF_D2AP_A2VP = contract_11_0( D_orbitals_Spin1OrbitalTensor,
1136 contract_21_1( Proj1_A, V_orbitals_Spin1OrbitalTensor ) );
1137 }
1138 if ( l == 2 )
1139 {
1140 SF_D2AP_A2VP = contract_11_0(
1141 D_orbitals_Spin1OrbitalTensor,
1142 contract_21_1( A_orbitals_Spin2OrbitalTensor, V_orbitals_Spin1OrbitalTensor ) );
1143 }
1144 return SF_D2AP_A2VP;
1145}
1146
1147//----D -> A P1 [P], A -> T P2 [P], T -> P3 P4 [D]
1148double EvtD0ToKSpipipi0::D2AP_A2TP( vector<double> p1, vector<double> p2, vector<double> p3,
1149 vector<double> p4, int l ) {
1150 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1151 assert( l >= 0 );
1152
1153 vector<double> T = p1 + p2;
1154 vector<double> P2 = p3;
1155 vector<double> A = T + P2;
1156 vector<double> P1 = p4;
1157 vector<double> D = A + P1;
1158
1159 vector<double> D_orbitals_Spin1OrbitalTensor = OrbitalTensors( D, A, P1, rD, 1 );
1160 vector<double> A_orbitals_Spin1OrbitalTensor = OrbitalTensors( A, T, P2, rRes, 1 );
1161 vector<double> T_orbitals_Spin2OrbitalTensor = OrbitalTensors( T, p1, p2, rRes, 2 );
1162 vector<double> Proj2_A = Proj( A, 2 );
1163
1164 double SF_D2AP_A2TP = 0.;
1165 // SF_D2AP_A2TP = &(D_orbitals.Spin1OrbitalTensor() | (T_orbitals.Spin2OrbitalTensor() |
1166 // A_orbitals.Spin1OrbitalTensor()));
1167 SF_D2AP_A2TP =
1168 contract_11_0( D_orbitals_Spin1OrbitalTensor,
1169 contract_21_1( contract_42_2( Proj2_A, T_orbitals_Spin2OrbitalTensor ),
1170 A_orbitals_Spin1OrbitalTensor ) );
1171 return SF_D2AP_A2TP;
1172}
1173
1174//----D -> T P1 [D], T -> V P2 [D], V -> P3 P4 [P]
1175double EvtD0ToKSpipipi0::D2TP_T2VP( vector<double> p1, vector<double> p2, vector<double> p3,
1176 vector<double> p4, int l ) {
1177 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1178 assert( l >= 0 );
1179
1180 vector<double> V = p1 + p2;
1181 vector<double> P2 = p3;
1182 vector<double> T = V + P2;
1183 vector<double> P1 = p4;
1184 vector<double> D = T + P1;
1185
1186 vector<double> D_orbitals_Spin2OrbitalTensor = OrbitalTensors( D, T, P1, rD, 2 );
1187 vector<double> T_orbitals_Spin2OrbitalTensor = OrbitalTensors( T, V, P2, rRes, 2 );
1188 vector<double> V_orbitals_Spin1OrbitalTensor = OrbitalTensors( V, p1, p2, rRes, 1 );
1189 vector<double> Proj1_T = Proj( T, 1 );
1190 vector<double> Proj2_T = Proj( T, 2 );
1191
1192 double SF_D2TP_T2VP = 0.;
1193 SF_D2TP_T2VP = contract_22_0(
1194 contract_42_2( Proj2_T, D_orbitals_Spin2OrbitalTensor ),
1195 contract_22_2( contract_31_2( contract_41_3( epsilon, T ),
1196 contract_21_1( Proj1_T, V_orbitals_Spin1OrbitalTensor ) ),
1197 T_orbitals_Spin2OrbitalTensor ) );
1198 return SF_D2TP_T2VP;
1199}
1200
1201//----D -> T1 P1 [D], T1 -> T2 P2 [P], T2 -> P3 P4 [D]
1202double EvtD0ToKSpipipi0::D2TP_T2TP( vector<double> p1, vector<double> p2, vector<double> p3,
1203 vector<double> p4, int l ) {
1204 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1205 assert( l >= 0 );
1206
1207 vector<double> T2 = p1 + p2;
1208 vector<double> P2 = p3;
1209 vector<double> T1 = T2 + P2;
1210 vector<double> P1 = p4;
1211 vector<double> D = T1 + P1;
1212
1213 vector<double> D_orbitals_Spin2OrbitalTensor = OrbitalTensors( D, T1, P1, rD, 2 );
1214 vector<double> T1_orbitals_Spin1OrbitalTensor = OrbitalTensors( T1, T2, P2, rRes, 1 );
1215 vector<double> T2_orbitals_Spin2OrbitalTensor = OrbitalTensors( T2, p1, p2, rRes, 2 );
1216 vector<double> Proj2_T1 = Proj( T1, 2 );
1217
1218 double SF_D2TP_T2TP = 0.;
1219 SF_D2TP_T2TP = contract_22_0( contract_22_2( contract_31_2( contract_41_3( epsilon, T1 ),
1220 T1_orbitals_Spin1OrbitalTensor ),
1221 T2_orbitals_Spin2OrbitalTensor ),
1222 contract_42_2( Proj2_T1, D_orbitals_Spin2OrbitalTensor ) );
1223 // SF_D2TP_T2TP = &( (((epsilon | T1) | T1_orbitals.Spin1OrbitalTensor()) || (Proj2_T1 |
1224 // D_orbitals.Spin2OrbitalTensor())) | T2_orbitals.Spin2OrbitalTensor() );//alternative
1225 return SF_D2TP_T2TP;
1226}
1227
1228//----D -> PT P1 [D], PT -> S P2 [D], S -> P3 P4 [S]
1229double EvtD0ToKSpipipi0::D2PTP_PT2SP( vector<double> p1, vector<double> p2, vector<double> p3,
1230 vector<double> p4, int l ) {
1231 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1232 assert( l >= 0 );
1233
1234 vector<double> S = p1 + p2;
1235 vector<double> P2 = p3;
1236 vector<double> PT = S + P2;
1237 vector<double> P1 = p4;
1238 vector<double> D = PT + P1;
1239
1240 vector<double> D_orbitals_Spin2OrbitalTensor = OrbitalTensors( D, PT, P1, rD, 2 );
1241 vector<double> PT_orbitals_Spin2OrbitalTensor = OrbitalTensors( PT, S, P2, rRes, 2 );
1242
1243 double SF_D2PTP_PT2SP = 0.;
1244 SF_D2PTP_PT2SP =
1245 contract_22_0( D_orbitals_Spin2OrbitalTensor, PT_orbitals_Spin2OrbitalTensor );
1246 return SF_D2PTP_PT2SP;
1247}
1248
1249//----D -> PT P1 [D], PT -> V P2 [P], V -> P3 P4 [P]
1250double EvtD0ToKSpipipi0::D2PTP_PT2VP( vector<double> p1, vector<double> p2, vector<double> p3,
1251 vector<double> p4, int l ) {
1252 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1253 assert( l >= 0 );
1254
1255 vector<double> V = p1 + p2;
1256 vector<double> P2 = p3;
1257 vector<double> PT = V + P2;
1258 vector<double> P1 = p4;
1259 vector<double> D = PT + P1;
1260
1261 vector<double> D_orbitals_Spin2OrbitalTensor = OrbitalTensors( D, PT, P1, rD, 2 );
1262 vector<double> PT_orbitals_Spin1OrbitalTensor = OrbitalTensors( PT, V, P2, rRes, 1 );
1263 vector<double> V_orbitals_Spin1OrbitalTensor = OrbitalTensors( V, p1, p2, rRes, 1 );
1264 vector<double> Proj2_PT = Proj( PT, 2 );
1265
1266 double SF_D2PTP_PT2VP = 0.;
1267 SF_D2PTP_PT2VP =
1268 contract_22_0( D_orbitals_Spin2OrbitalTensor,
1269 contract_31_2( contract_41_3( Proj2_PT, V_orbitals_Spin1OrbitalTensor ),
1270 PT_orbitals_Spin1OrbitalTensor ) );
1271 return SF_D2PTP_PT2VP;
1272}
1273
1274//----D -> PT P1 [D], PT -> T P2 [S]
1275double EvtD0ToKSpipipi0::D2PTP_PT2TP( vector<double> p1, vector<double> p2, vector<double> p3,
1276 vector<double> p4, int l ) {
1277 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1278 assert( l >= 0 );
1279
1280 vector<double> T = p1 + p2;
1281 vector<double> P2 = p3;
1282 vector<double> PT = T + P2;
1283 vector<double> P1 = p4;
1284 vector<double> D = PT + P1;
1285
1286 vector<double> D_orbitals_Spin2OrbitalTensor = OrbitalTensors( D, PT, P1, rD, 2 );
1287 vector<double> T_orbitals_Spin2OrbitalTensor = OrbitalTensors( T, p1, p2, rRes, 2 );
1288 vector<double> Proj2_PT = Proj( PT, 2 );
1289
1290 double SF_D2PTP_PT2TP = 0.;
1291 SF_D2PTP_PT2TP = contract_22_0( D_orbitals_Spin2OrbitalTensor,
1292 contract_42_2( Proj2_PT, T_orbitals_Spin2OrbitalTensor ) );
1293 return SF_D2PTP_PT2TP;
1294}
1295
1296//------------------------------------------------------------//
1297// Quasi Two-Body Decay //
1298//------------------------------------------------------------//
1299//----D -> S1 S2 [S]
1300double EvtD0ToKSpipipi0::D2SS() {
1301
1302 double SF_D2SS = 1.0;
1303 return SF_D2SS;
1304}
1305
1306//----D->V S [P]
1307double EvtD0ToKSpipipi0::D2VS( vector<double> p1, vector<double> p2, vector<double> p3,
1308 vector<double> p4, int l ) {
1309 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1310 assert( l >= 0 );
1311
1312 vector<double> V = p1 + p2;
1313 vector<double> S = p3 + p4;
1314 vector<double> D = V + S;
1315
1316 vector<double> D_orbitals_Spin1OrbitalTensor = OrbitalTensors( D, V, S, rD, 1 );
1317 vector<double> V_orbitals_Spin1OrbitalTensor = OrbitalTensors( V, p1, p2, rRes, 1 );
1318
1319 double SF_D2VS = 0.;
1320 SF_D2VS = contract_11_0( D_orbitals_Spin1OrbitalTensor, V_orbitals_Spin1OrbitalTensor );
1321 return SF_D2VS;
1322}
1323
1324//----D -> V1 V2 [S, P, D], V1 -> P1 P2 [P], V2 -> P3 P4 [P]
1325double EvtD0ToKSpipipi0::D2VV( vector<double> p1, vector<double> p2, vector<double> p3,
1326 vector<double> p4, int l ) {
1327 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1328 assert( l >= 0 );
1329
1330 vector<double> V1 = p1 + p2;
1331 vector<double> V2 = p3 + p4;
1332 vector<double> D = V1 + V2;
1333
1334 vector<double> D_orbitals_Spin1OrbitalTensor = OrbitalTensors( D, V1, V2, rD, 1 );
1335 vector<double> D_orbitals_Spin2OrbitalTensor = OrbitalTensors( D, V1, V2, rD, 2 );
1336 vector<double> V1_orbitals_Spin1OrbitalTensor = OrbitalTensors( V1, p1, p2, rRes, 1 );
1337 vector<double> V2_orbitals_Spin1OrbitalTensor = OrbitalTensors( V2, p3, p4, rRes, 1 );
1338
1339 double SF_D2VV = 0;
1340 if ( l == 0 )
1341 {
1342 SF_D2VV = contract_11_0( V1_orbitals_Spin1OrbitalTensor, V2_orbitals_Spin1OrbitalTensor );
1343 }
1344 if ( l == 1 )
1345 {
1346 SF_D2VV = contract_11_0( contract_21_1( contract_31_2( contract_41_3( epsilon, D ),
1347 V2_orbitals_Spin1OrbitalTensor ),
1348 V1_orbitals_Spin1OrbitalTensor ),
1349 D_orbitals_Spin1OrbitalTensor );
1350 }
1351 if ( l == 2 )
1352 {
1353 SF_D2VV = contract_11_0(
1354 contract_21_1( D_orbitals_Spin2OrbitalTensor, V2_orbitals_Spin1OrbitalTensor ),
1355 V1_orbitals_Spin1OrbitalTensor );
1356 }
1357
1358 return SF_D2VV;
1359}
1360
1361//----D -> T S [D], T -> P1 P2 [D], S -> P3 P4 [S]
1362double EvtD0ToKSpipipi0::D2TS( vector<double> p1, vector<double> p2, vector<double> p3,
1363 vector<double> p4, int l ) {
1364 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1365 assert( l >= 0 );
1366
1367 vector<double> T = p1 + p2;
1368 vector<double> S = p3 + p4;
1369 vector<double> D = T + S;
1370
1371 vector<double> D_orbitals_Spin2OrbitalTensor = OrbitalTensors( D, T, S, rD, 2 );
1372 vector<double> T_orbitals_Spin2OrbitalTensor = OrbitalTensors( T, p1, p2, rRes, 2 );
1373
1374 double SF_D2TS = 0.;
1375 SF_D2TS = contract_22_0( D_orbitals_Spin2OrbitalTensor, T_orbitals_Spin2OrbitalTensor );
1376 return SF_D2TS;
1377}
1378
1379//----D -> T V [P, D], T -> P1 P2 [D], V -> P3 P4 [P]
1380double EvtD0ToKSpipipi0::D2TV( vector<double> p1, vector<double> p2, vector<double> p3,
1381 vector<double> p4, int l ) {
1382 assert( p1.size() == 4 && p2.size() == 4 && p3.size() == 4 && p4.size() == 4 );
1383 assert( l >= 0 );
1384
1385 vector<double> T = p1 + p2;
1386 vector<double> V = p3 + p4;
1387 vector<double> D = T + V;
1388
1389 vector<double> D_orbitals_Spin1OrbitalTensor = OrbitalTensors( D, T, V, rD, 1 );
1390 vector<double> D_orbitals_Spin2OrbitalTensor = OrbitalTensors( D, T, V, rD, 2 );
1391 vector<double> T_orbitals_Spin2OrbitalTensor = OrbitalTensors( T, p1, p2, rRes, 2 );
1392 vector<double> V_orbitals_Spin1OrbitalTensor = OrbitalTensors( V, p3, p4, rRes, 1 );
1393
1394 double SF_D2TV = 0.;
1395 if ( l == 1 )
1396 {
1397 SF_D2TV = contract_11_0(
1398 D_orbitals_Spin1OrbitalTensor,
1399 contract_21_1( T_orbitals_Spin2OrbitalTensor, V_orbitals_Spin1OrbitalTensor ) );
1400 }
1401 if ( l == 2 )
1402 {
1403 SF_D2TV = contract_22_0(
1404 contract_31_2( contract_41_3( epsilon, D ), V_orbitals_Spin1OrbitalTensor ),
1405 contract_22_2( D_orbitals_Spin2OrbitalTensor, T_orbitals_Spin2OrbitalTensor ) );
1406 }
1407 return SF_D2TV;
1408}
1409
1410//----D -> T1 T2 [S, P, D], T1 -> P1 P2 [D], T2 -> P3 P4 [D]
1411// ok
1412double EvtD0ToKSpipipi0::D2TT( vector<double> p1, vector<double> p2, vector<double> p3,
1413 vector<double> p4, int l ) {
1414
1415 vector<double> T1 = p1 + p2;
1416 vector<double> T2 = p3 + p4;
1417 vector<double> D = T1 + T2;
1418
1419 vector<double> D_orbitals_Spin1OrbitalTensor = OrbitalTensors( D, T1, T2, rD, 1 );
1420 vector<double> D_orbitals_Spin2OrbitalTensor = OrbitalTensors( D, T1, T2, rD, 2 );
1421 vector<double> T1_orbitals_Spin2OrbitalTensor = OrbitalTensors( T1, p1, p2, rRes, 2 );
1422 vector<double> T2_orbitals_Spin2OrbitalTensor = OrbitalTensors( T2, p3, p4, rRes, 2 );
1423
1424 double SF_D2TT = 0.;
1425 if ( l == 0 )
1426 {
1427 SF_D2TT = contract_22_0( T1_orbitals_Spin2OrbitalTensor, T2_orbitals_Spin2OrbitalTensor );
1428 }
1429 if ( l == 1 )
1430 {
1431 SF_D2TT = contract_22_0(
1432 contract_31_2( contract_41_3( epsilon, D ), D_orbitals_Spin1OrbitalTensor ),
1433 contract_22_2( T1_orbitals_Spin2OrbitalTensor, T2_orbitals_Spin2OrbitalTensor ) );
1434 }
1435 if ( l == 2 )
1436 {
1437 SF_D2TT = contract_22_0(
1438 D_orbitals_Spin2OrbitalTensor,
1439 contract_22_2( T1_orbitals_Spin2OrbitalTensor, T2_orbitals_Spin2OrbitalTensor ) );
1440 }
1441 return SF_D2TT;
1442}
1443
1444//------------------------------------------------------------//
1445// Propagators //
1446//------------------------------------------------------------//
1447double EvtD0ToKSpipipi0::fundecaymomentum( double mr2, double m1_2, double m2_2 ) {
1448 double mr = sqrt( mr2 );
1449 double poly = mr2 * mr2 + m1_2 * m1_2 + m2_2 * m2_2 - 2.0 * m1_2 * mr2 - 2.0 * m2_2 * mr2 -
1450 2.0 * m1_2 * m2_2;
1451 double ret = sqrt( poly ) / ( 2.0 * mr );
1452 if ( poly < 0 ) ret = 0.0;
1453 return ret;
1454}
1455
1456double EvtD0ToKSpipipi0::fundecaymomentum2( double mr2, double m1_2, double m2_2 ) {
1457 double poly = mr2 * mr2 + m1_2 * m1_2 + m2_2 * m2_2 - 2.0 * m1_2 * mr2 - 2.0 * m2_2 * mr2 -
1458 2.0 * m1_2 * m2_2;
1459 double ret = poly / ( 4.0 * mr2 );
1460 if ( poly < 0 ) ret = 0.0;
1461 return ret;
1462}
1463
1464double EvtD0ToKSpipipi0::wid( double mass, double sa, double sb, double sc, double r, int l ) {
1465
1466 double widm = 1.0;
1467 double sa0 = mass * mass;
1468 double m = sqrt( sa );
1469 double q = fundecaymomentum2( sa, sb, sc );
1470 double q0 = fundecaymomentum2( sa0, sb, sc );
1471 double z = q * r * r;
1472 double z0 = q0 * r * r;
1473 double F = 0.0;
1474 if ( l == 0 ) F = 1.0;
1475 if ( l == 1 ) F = sqrt( ( 1.0 + z0 ) / ( 1.0 + z ) );
1476 if ( l == 2 ) F = sqrt( ( 9.0 + 3.0 * z0 + z0 * z0 ) / ( 9.0 + 3.0 * z + z * z ) );
1477 if ( l == 3 )
1478 F = sqrt( ( 225.0 + 45.0 * z0 + 6.0 * z0 * z0 + z0 * z0 * z0 ) /
1479 ( 225.0 + 45.0 * z + 6.0 * z * z + z * z * z ) );
1480 if ( l == 4 )
1481 F = sqrt(
1482 ( 11025.0 + 1575.0 * z0 + 135.0 * z0 * z0 + 10.0 * z0 * z0 * z0 + z0 * z0 * z0 * z0 ) /
1483 ( 11025.0 + 1575.0 * z + 135.0 * z * z + 10.0 * z * z * z + z * z * z * z ) );
1484
1485 double t = sqrt( q / q0 );
1486 for ( uint i = 0; i < ( 2 * l + 1 ); i++ ) { widm *= t; }
1487 widm *= ( mass / m * F * F );
1488 return widm;
1489}
1490
1491EvtComplex EvtD0ToKSpipipi0::BW( double mx2, double mr, double wr ) {
1492
1493 double mr2 = mr * mr;
1494 double denom_real = mr2 - mx2;
1495 double denom_imag = mr * wr;
1496 double denom = denom_real * denom_real + denom_imag * denom_imag;
1497
1498 double output_x = 0., output_y = 0.;
1499 if ( wr < 0 ) { output_x = 1 / sqrt( denom ); }
1500 else if ( wr < 10 )
1501 {
1502 output_x = denom_real / denom;
1503 output_y = denom_imag / denom;
1504 }
1505 else
1506 { /* phase space */
1507 output_x = 1.;
1508 output_y = 0.;
1509 }
1510
1511 EvtComplex output( output_x, output_y );
1512
1513 return output;
1514}
1515
1516EvtComplex EvtD0ToKSpipipi0::RBW( double mx2, double mr, double wr, double m1_2, double m2_2,
1517 double r, int l ) {
1518
1519 double mr2 = mr * mr;
1520 double denom_real = mr2 - mx2;
1521 double denom_imag = mr * wr * wid( mr, mx2, m1_2, m2_2, r, l ); // real-i*imag;
1522
1523 double denom = denom_real * denom_real + denom_imag * denom_imag;
1524 double output_x = denom_real / denom;
1525 double output_y = denom_imag / denom;
1526
1527 EvtComplex output( output_x, output_y );
1528
1529 return output;
1530}
1531
1532/* KPi S-wave LASS Model */
1533EvtComplex EvtD0ToKSpipipi0::LASS( double mx2, double m1_2, double m2_2 ) {
1534 double sa = mx2;
1535 double sb = m1_2;
1536 double sc = m2_2;
1537 double m1430 = 1.463;
1538 double w1430 = 0.233;
1539 double sa0 = m1430 * m1430;
1540 double q0 = ( sa0 + sb - sc ) * ( sa0 + sb - sc ) / ( 4.0 * sa0 ) - sb;
1541 if ( q0 < 0 ) q0 = 1e-16;
1542 double qs = ( sa + sb - sc ) * ( sa + sb - sc ) / ( 4.0 * sa ) - sb;
1543 double q = sqrt( qs );
1544 double width = w1430 * q * m1430 / sqrt( sa * q0 );
1545 double temp_R = atan( m1430 * width / ( sa0 - sa ) );
1546 if ( temp_R < 0 ) temp_R += math_pi;
1547 double deltaR = -5.31 + temp_R;
1548 double temp_F = atan( 2.0 * 1.07 * q / ( 2.0 + ( -1.8 ) * 1.07 * qs ) );
1549 if ( temp_F < 0 ) temp_F += math_pi;
1550 double deltaF = 2.33 + temp_F;
1551 double output_x =
1552 0.8 * sin( deltaF ) * cos( deltaF ) + sin( deltaR ) * cos( deltaR + 2.0 * deltaF );
1553 double output_y =
1554 0.8 * sin( deltaF ) * sin( deltaF ) + sin( deltaR ) * sin( deltaR + 2.0 * deltaF );
1555
1556 EvtComplex output( output_x, output_y );
1557 return output;
1558}
1559
1560/* GS propagator */
1561const double mass_Pion = 0.13957;
1562double EvtD0ToKSpipipi0::h( double m, double q ) {
1563 double h = 2.0 / math_pi * q / m * log( ( m + 2.0 * q ) / ( 2.0 * mass_Pion ) );
1564 return h;
1565}
1566
1567double EvtD0ToKSpipipi0::dh( double m0, double q0 ) {
1568 double dh = h( m0, q0 ) * ( 1.0 / ( 8.0 * q0 * q0 ) - 1.0 / ( 2.0 * m0 * m0 ) ) +
1569 1.0 / ( 2.0 * math_pi * m0 * m0 );
1570 return dh;
1571}
1572
1573double EvtD0ToKSpipipi0::f( double m0, double sx, double q0, double q ) {
1574 double m = sqrt( sx );
1575 double f =
1576 m0 * m0 / ( q0 * q0 * q0 ) *
1577 ( q * q * ( h( m, q ) - h( m0, q0 ) ) + ( m0 * m0 - sx ) * q0 * q0 * dh( m0, q0 ) );
1578 return f;
1579}
1580
1581double EvtD0ToKSpipipi0::d( double m0, double q0 ) {
1582 double d = 3.0 / math_pi * mass_Pion * mass_Pion / ( q0 * q0 ) *
1583 log( ( m0 + 2.0 * q0 ) / ( 2.0 * mass_Pion ) ) +
1584 m0 / ( 2.0 * math_pi * q0 ) -
1585 ( mass_Pion * mass_Pion * m0 ) / ( math_pi * q0 * q0 * q0 );
1586 return d;
1587}
1588
1589EvtComplex EvtD0ToKSpipipi0::GS( double mx2, double mr, double wr, double m1_2, double m2_2,
1590 double r, int l ) {
1591
1592 double mr2 = mr * mr;
1593 double q = fundecaymomentum( mx2, m1_2, m2_2 );
1594 double q0 = fundecaymomentum( mr2, m1_2, m2_2 );
1595 double numer = 1.0 + d( mr, q0 ) * wr / mr;
1596 double denom_real = mr2 - mx2 + wr * f( mr, mx2, q0, q );
1597 double denom_imag = mr * wr * wid( mr, mx2, m1_2, m2_2, r, l ); // real-i*imag;
1598
1599 double denom = denom_real * denom_real + denom_imag * denom_imag;
1600 double output_x = denom_real * numer / denom;
1601 double output_y = denom_imag * numer / denom;
1602
1603 EvtComplex output( output_x, output_y );
1604 return output;
1605}
1606
1607EvtComplex EvtD0ToKSpipipi0::create_RBW_propagator( string name, double mx2, double m1_2,
1608 double m2_2, int l ) {
1609 double mr = resonance_par[name + "_mass"];
1610 double wr = resonance_par[name + "_width"];
1611 return RBW( mx2, mr, wr, m1_2, m2_2, rRes, l );
1612}
1613
1614EvtComplex EvtD0ToKSpipipi0::create_BW_propagator( string name, double mx2 ) {
1615 double mr = resonance_par[name + "_mass"];
1616 double wr = resonance_par[name + "_width"];
1617 return BW( mx2, mr, wr );
1618}
1619
1620EvtComplex EvtD0ToKSpipipi0::create_GS_propagator( string name, double mx2 ) { // for rho0
1621 double mr = resonance_par[name + "_mass"];
1622 double wr = resonance_par[name + "_width"];
1623 return GS( mx2, mr, wr, mpip * mpip, mpim * mpim, rRes, 1 );
1624}
1625
1626EvtComplex EvtD0ToKSpipipi0::create_KPiSLASS_propagator( string name, double mx2, double m1_2,
1627 double m2_2 ) {
1628 return LASS( mx2, m1_2, m2_2 );
1629}
1630
1631/* sigma propagator for sigma(500)*/
1632double EvtD0ToKSpipipi0::rho4pi( double s ) {
1633 // double mpi=0.13957018;
1634 double mpi = 0.13957;
1635 double tmp = 1.0 - 16.0 * mpi * mpi / s;
1636 double ret = 0.;
1637 if ( tmp <= 0 ) ret = 0.0;
1638 else ret = sqrt( tmp ) / ( 1.0 + exp( 9.8f - 3.5 * s ) );
1639
1640 return ret;
1641}
1642
1643double EvtD0ToKSpipipi0::rho2pi( double s ) {
1644 // double mpi=0.13957018;
1645 double mpi = 0.13957;
1646 double ret = 0.;
1647 if ( s <= 4.0 * mpi * mpi ) ret = 0.0;
1648 else ret = sqrt( 1.0 - 4.0 * mpi * mpi / s );
1649
1650 return ret;
1651}
1652
1653EvtComplex EvtD0ToKSpipipi0::sigma( double mx2, double mr, double gf ) {
1654
1655 // double mr = 0.9264f;
1656 double mr2 = mr * mr;
1657 double diff = mr2 - mx2; /* m*m-s */
1658 // double mpi=0.13957018;
1659 double mpi = 0.13957;
1660 double g1 = ( 0.5843 + 1.6663 * mx2 ) * ( mx2 - mpi * mpi / 2.0 ) /
1661 ( mr2 - mpi * mpi / 2.0 ) * exp( -( mx2 - mr2 ) / 1.082 );
1662 // double g2=0.0024f;
1663 double g2 = gf;
1664 double w1 = g1 * rho2pi( mx2 ) / rho2pi( mr2 );
1665 double w2 = g2 * rho4pi( mx2 ) / rho4pi( mr2 );
1666 double ws = mr * ( w1 + w2 );
1667 double denom = diff * diff + ws * ws;
1668
1669 double output_x = diff / denom;
1670 double output_y = ws / denom;
1671
1672 EvtComplex output( output_x, output_y );
1673
1674 return output;
1675}
1676
1677EvtComplex EvtD0ToKSpipipi0::create_sigma_propagator( string name, double mx2 ) {
1678
1679 double mr = resonance_par[name + "_mass"];
1680 double wr = resonance_par[name + "_width"];
1681
1682 double gf = resonance_par[name + "_gf"];
1683
1684 return sigma( mx2, mr, gf );
1685}
1686
1687/* Flatte propagator for f0(980) and a0(980)*/
1688EvtComplex EvtD0ToKSpipipi0::irho( double mr2, double m1_2, double m2_2 ) {
1689 double mr = sqrt( mr2 );
1690 double poly = mr2 * mr2 + m1_2 * m1_2 + m2_2 * m2_2 - 2 * m1_2 * mr2 - 2 * m2_2 * mr2 -
1691 2 * m1_2 * m2_2;
1692
1693 double output_x = 0., output_y = 0.;
1694 if ( poly >= 0 )
1695 {
1696 output_x = 2.0 * sqrt( poly ) / ( 2.0 * mr2 );
1697 output_y = 0.;
1698 }
1699 if ( poly < 0 )
1700 {
1701 output_x = 0.;
1702 output_y = 2.0 * sqrt( -1.0 * poly ) / ( 2.0 * mr2 );
1703 }
1704
1705 EvtComplex output( output_x, output_y );
1706
1707 return output;
1708}
1709
1710EvtComplex EvtD0ToKSpipipi0::Flatte2( double mx2, double mr2, double g1, double m1a,
1711 double m1b, double g2, double m2a, double m2b ) {
1712
1713 double diff = mr2 - mx2; /* m*m-s */
1714 g2 = g2 * g1;
1715 EvtComplex ps1 = irho( mx2, m1a * m1a, m1b * m1b );
1716 EvtComplex ps2 = irho( mx2, m2a * m2a, m2b * m2b );
1717 EvtComplex iws = g1 * ps1 + g2 * ps2; /*mass dependent width */
1718
1719 diff += imag( iws );
1720 double ws = real( iws );
1721 double denom = diff * diff + ws * ws; /* common denominator*/
1722
1723 double output_x = diff / denom; /* real part*/
1724 double output_y = ws / denom; /* imaginary part*/
1725
1726 EvtComplex output( output_x, output_y );
1727
1728 return output;
1729}
1730
1731EvtComplex EvtD0ToKSpipipi0::create_Flatte2_propagator( string name, double mx2, double m1a,
1732 double m1b, double m2a, double m2b ) {
1733 double mr = resonance_par[name + "_mass"];
1734 double wr = resonance_par[name + "_width"];
1735 double mr2 = mr * mr;
1736
1737 double g1 = resonance_par[name + "_g1"];
1738 double g2 = resonance_par[name + "_g2"];
1739
1740 return Flatte2( mx2, mr2, g1, m1a, m1b, g2, m2a, m2b );
1741}
1742
1743void EvtD0ToKSpipipi0::readInputCoeff() {
1744 coefficient["phasespace"] = 0;
1745 coefficient["<kstarm_892|rhop_770>_0_mag"] = 25;
1746 coefficient["<kstarm_892|rhop_770>_0_phase"] = 0;
1747 coefficient["<kstarm_892|rhop_770>_1_mag"] = 11.2665;
1748 coefficient["<kstarm_892|rhop_770>_1_phase"] = 0.548994;
1749 coefficient["<kstarm_892|rhop_770>_2_mag"] = 6.73954;
1750 coefficient["<kstarm_892|rhop_770>_2_phase"] = 1.96845;
1751 coefficient["<kstar0_892|rho0_770>_0_mag"] = 9.62259;
1752 coefficient["<kstar0_892|rho0_770>_0_phase"] = 4.39676;
1753 coefficient["<kstar0_892|rho0_770>_1_mag"] = 3.49377;
1754 coefficient["<kstar0_892|rho0_770>_1_phase"] = -2.56438;
1755 coefficient["<kstar0_892|rho0_770>_2_mag"] = 3.30796;
1756 coefficient["<kstar0_892|rho0_770>_2_phase"] = 5.17417;
1757 coefficient["<omega3_782|rho_770>_1_mag"] = 14.1209;
1758 coefficient["<omega3_782|rho_770>_1_phase"] = 2.06989;
1759 coefficient["<phi_1020|rho_770>_1_mag"] = 0.835627;
1760 coefficient["<phi_1020|rho_770>_1_phase"] = -2.92636;
1761 coefficient["<eta_547|rhom_phsp>_1_mag"] = 81.3879;
1762 coefficient["<eta_547|rhom_phsp>_1_phase"] = -1.51303;
1763 coefficient["<eta_547|rhom_phsp>_0_mag"] = 72.5434;
1764 coefficient["<eta_547|rhom_phsp>_0_phase"] = 3.89337;
1765 coefficient["<k1m_1270|rhom_770>_0_mag"] = 48.1147;
1766 coefficient["<k1m_1270|rhom_770>_0_phase"] = -0.0202615;
1767 coefficient["<k1m_1270|kstar_892>_0_mag"] = 5.73667;
1768 coefficient["<k1m_1270|kstar_892>_0_phase"] = 3.92441;
1769 coefficient["<k1m_1270|k0star_1430>_1_mag"] = 72.8988;
1770 coefficient["<k1m_1270|k0star_1430>_1_phase"] = -1.12332;
1771 coefficient["<k10_1270|rho0_770>_0_mag"] = 14.6425;
1772 coefficient["<k10_1270|rho0_770>_0_phase"] = 0.882522;
1773 coefficient["<k10_1270|kstarm_892>_0_mag"] = 4.40263;
1774 coefficient["<k10_1270|kstarm_892>_0_phase"] = 0.299025;
1775 coefficient["<k10_1270|k0starm_1430>_1_mag"] = 57.4271;
1776 coefficient["<k10_1270|k0starm_1430>_1_phase"] = 1.25342;
1777 coefficient["<k1m_1400|kstar_892>_0_mag"] = 7.41837;
1778 coefficient["<k1m_1400|kstar_892>_0_phase"] = 2.60505;
1779 coefficient["<k10_1400|kstarm_892>_0_mag"] = -38.728;
1780 coefficient["<k10_1400|kstarm_892>_0_phase"] = 2.64084;
1781 coefficient["<kstar02_1410|rho0_770>_1_mag"] = 32.6332;
1782 coefficient["<kstar02_1410|rho0_770>_1_phase"] = 0.266067;
1783 coefficient["<kstarm2_1410|rhop_770>_0_mag"] = 84.8711;
1784 coefficient["<kstarm2_1410|rhop_770>_0_phase"] = -1.22924;
1785 coefficient["<kstarm2_1410|rhop_770>_1_mag"] = 30.6599;
1786 coefficient["<kstarm2_1410|rhop_770>_1_phase"] = -3.02164;
1787 coefficient["<kstarm2_1410|rhop_770>_2_mag"] = 23.5133;
1788 coefficient["<kstarm2_1410|rhop_770>_2_phase"] = 4.90416;
1789 coefficient["<kstarm3_1410|kstar_892>_1_mag"] = -5.00895;
1790 coefficient["<kstarm3_1410|kstar_892>_1_phase"] = 0.75297;
1791 coefficient["<kstarm3_1410|rhom_770>_1_mag"] = 21.0146;
1792 coefficient["<kstarm3_1410|rhom_770>_1_phase"] = 2.11827;
1793 coefficient["<km_1460|kstar_892>_1_mag"] = -18.5953;
1794 coefficient["<km_1460|kstar_892>_1_phase"] = 1.72652;
1795 coefficient["<km_1460|k0star_1430>_0_mag"] = 500;
1796 coefficient["<km_1460|k0star_1430>_0_phase"] = 3.22602;
1797 coefficient["<k0_1460|kstarm_892>_1_mag"] = 14.7543;
1798 coefficient["<k0_1460|kstarm_892>_1_phase"] = 0.510956;
1799 coefficient["<k1m_1650|kstar_892>_0_mag"] = 11.7635;
1800 coefficient["<k1m_1650|kstar_892>_0_phase"] = 2.59972;
1801 coefficient["<k10_1650|kstarm_892>_0_mag"] = 27.9014;
1802 coefficient["<k10_1650|kstarm_892>_0_phase"] = -1.07364;
1803 coefficient["<k1m_1650|rhom_770>_0_mag"] = 23.2323;
1804 coefficient["<k1m_1650|rhom_770>_0_phase"] = -1.33583;
1805 coefficient["<kstar03_1680|kstarm_892>_1_mag"] = 23.4824;
1806 coefficient["<kstar03_1680|kstarm_892>_1_phase"] = 0.139691;
1807 coefficient["<kstarm3_1680|kstar_892>_1_mag"] = 11.8553;
1808 coefficient["<kstarm3_1680|kstar_892>_1_phase"] = -3.47852;
1809 coefficient["<k0star0_1430|rho0_phsp>_0_mag"] = 500;
1810 coefficient["<k0star0_1430|rho0_phsp>_0_phase"] = -5.1306;
1811 coefficient["<k1p_phsp|rhop_770>_1_mag"] = 197.192;
1812 coefficient["<k1p_phsp|rhop_770>_1_phase"] = 0.35067;
1813 coefficient["<k1m_phsp|rhom_770>_1_mag"] = 151.864;
1814 coefficient["<k1m_phsp|rhom_770>_1_phase"] = -5.11224;
1815 coefficient["<k1m_phsp|rhom_770>_2_mag"] = 21.6383;
1816 coefficient["<k1m_phsp|rhom_770>_2_phase"] = -0.843621;
1817 coefficient["<k0starm_1430|rhop_phsp>_0_mag"] = -500;
1818 coefficient["<k0starm_1430|rhop_phsp>_0_phase"] = 4.00277;
1819 coefficient["<k0_1460|sigma_500>_0_mag"] = 500;
1820 coefficient["<k0_1460|sigma_500>_0_phase"] = -5.70929;
1821 coefficient["<k10_1270|f0_980>_1_mag"] = 106.293;
1822 coefficient["<k10_1270|f0_980>_1_phase"] = 0.947063;
1823 coefficient["<k10_1400|f0_1370>_1_mag"] = 184.717;
1824 coefficient["<k10_1400|f0_1370>_1_phase"] = -0.0197254;
1825 coefficient["<rho0_770|k0star0_1430>_1_mag"] = 159.521;
1826 coefficient["<rho0_770|k0star0_1430>_1_phase"] = 0.164038;
1827 coefficient["<k2starm3_1430|kstar_892>_2_mag"] = 0.400134;
1828 coefficient["<k2starm3_1430|kstar_892>_2_phase"] = 0.533053;
1829 coefficient["<k2star02_1430|f2_1270>_1_mag"] = 45.2478;
1830 coefficient["<k2star02_1430|f2_1270>_1_phase"] = 1.81556;
1831 coefficient["<kstar0_892|f0_980>_1_mag"] = 24.2878;
1832 coefficient["<kstar0_892|f0_980>_1_phase"] = 1.93217;
1833 coefficient["<kstarm_phsp|rhop_phsp>_0_mag"] = 500;
1834 coefficient["<kstarm_phsp|rhop_phsp>_0_phase"] = 1.48706;
1835 coefficient["<rho0_770|kstar0_phsp>_1_mag"] = 58.8657;
1836 coefficient["<rho0_770|kstar0_phsp>_1_phase"] = 3.54134;
1837}
1838
1839double EvtD0ToKSpipipi0::amps() {
1840
1841 initPar();
1842
1843 // (E; PX, PY, PZ)
1844 EvtVector4R ks0 = _pd[0];
1845 EvtVector4R pip = _pd[1];
1846 EvtVector4R pim = _pd[2];
1847 EvtVector4R pi0 = _pd[3];
1848
1849 createPropagator( ks0, pip, pim, pi0 );
1850 createSpinfactor( ks0, pip, pim, pi0 );
1851
1852 //----------------------------------------------------------------------//
1853 // Add the partial waves here //
1854 //----------------------------------------------------------------------//
1855
1856 double f_isoConj = -1. * sqrt( 2.0 );
1857 double one = 1.0, minus_one = -1.0;
1858 addPartialWave( spinfactor["D2VV_ks0pim_pippi0_0"], propagator["kstarm_892"],
1859 propagator["rhop_770"], coefficient["<kstarm_892|rhop_770>_0_mag"],
1860 coefficient["<kstarm_892|rhop_770>_0_phase"] );
1861 addPartialWave( spinfactor["D2VV_ks0pim_pippi0_1"], propagator["kstarm_892"],
1862 propagator["rhop_770"], coefficient["<kstarm_892|rhop_770>_1_mag"],
1863 coefficient["<kstarm_892|rhop_770>_1_phase"] );
1864 addPartialWave( spinfactor["D2VV_ks0pim_pippi0_2"], propagator["kstarm_892"],
1865 propagator["rhop_770"], coefficient["<kstarm_892|rhop_770>_2_mag"],
1866 coefficient["<kstarm_892|rhop_770>_2_phase"] );
1867 addPartialWave( spinfactor["D2VV_ks0pi0_pippim_0"], propagator["kstar0_892"],
1868 propagator["rho0_770"], coefficient["<kstar0_892|rho0_770>_0_mag"],
1869 coefficient["<kstar0_892|rho0_770>_0_phase"] );
1870 addPartialWave( spinfactor["D2VV_ks0pi0_pippim_1"], propagator["kstar0_892"],
1871 propagator["rho0_770"], coefficient["<kstar0_892|rho0_770>_1_mag"],
1872 coefficient["<kstar0_892|rho0_770>_1_phase"] );
1873 addPartialWave( spinfactor["D2VV_ks0pi0_pippim_2"], propagator["kstar0_892"],
1874 propagator["rho0_770"], coefficient["<kstar0_892|rho0_770>_2_mag"],
1875 coefficient["<kstar0_892|rho0_770>_2_phase"] );
1876 EvtComplex omega_ks0_propagator =
1877 propagator["rhop_770"] * spinfactor["D2VP_V2VP_pippimpi0_pippi0_1"] +
1878 minus_one * propagator["rho0_770"] * spinfactor["D2VP_V2VP_pippimpi0_pippim_1"] +
1879 propagator["rhom_770"] * spinfactor["D2VP_V2VP_pippimpi0_pimpi0_1"];
1880 addPartialWave( omega_ks0_propagator, propagator["omega3_782"],
1881 coefficient["<omega3_782|rho_770>_1_mag"],
1882 coefficient["<omega3_782|rho_770>_1_phase"] );
1883 EvtComplex phi_ks0_propagator =
1884 propagator["rhop_770"] * spinfactor["D2VP_V2VP_pippimpi0_pippi0_1"] +
1885 minus_one * propagator["rho0_770"] * spinfactor["D2VP_V2VP_pippimpi0_pippim_1"] +
1886 propagator["rhom_770"] * spinfactor["D2VP_V2VP_pippimpi0_pimpi0_1"];
1887 addPartialWave( phi_ks0_propagator, propagator["phi_1020"],
1888 coefficient["<phi_1020|rho_770>_1_mag"],
1889 coefficient["<phi_1020|rho_770>_1_phase"] );
1890 addPartialWave( spinfactor["D2PP_P2VP_pippimpi0_pimpi0_1"], propagator["eta_547"],
1891 propagator["rhom_phsp"], coefficient["<eta_547|rhom_phsp>_1_mag"],
1892 coefficient["<eta_547|rhom_phsp>_1_phase"] );
1893 addPartialWave( spinfactor["D2PP_P2SP_pippimpi0_pimpi0_0"], propagator["eta_547"],
1894 propagator["rhom_phsp"], coefficient["<eta_547|rhom_phsp>_0_mag"],
1895 coefficient["<eta_547|rhom_phsp>_0_phase"] );
1896 addPartialWave( spinfactor["D2AP_A2VP_ks0pimpi0_pimpi0_0"], propagator["k1m_1270"],
1897 propagator["rhom_770"], coefficient["<k1m_1270|rhom_770>_0_mag"],
1898 coefficient["<k1m_1270|rhom_770>_0_phase"] );
1899 addPartialWave( spinfactor["D2AP_A2VP_ks0pimpi0_ks0pim_0"], propagator["k1m_1270"],
1900 propagator["kstarm_892"], spinfactor["D2AP_A2VP_ks0pimpi0_ks0pi0_0"],
1901 propagator["k1m_1270"], propagator["kstar0_892"], f_isoConj,
1902 coefficient["<k1m_1270|kstar_892>_0_mag"],
1903 coefficient["<k1m_1270|kstar_892>_0_phase"] );
1904 addPartialWave( spinfactor["D2AP_A2SP_ks0pimpi0_ks0pim_1"], propagator["k1m_1270"],
1905 propagator["k0starm_1430"], spinfactor["D2AP_A2SP_ks0pimpi0_ks0pi0_1"],
1906 propagator["k1m_1270"], propagator["k0star0_1430"], f_isoConj,
1907 coefficient["<k1m_1270|k0star_1430>_1_mag"],
1908 coefficient["<k1m_1270|k0star_1430>_1_phase"] );
1909 addPartialWave( spinfactor["D2AP_A2VP_ks0pippim_pippim_0"], propagator["k10_1270"],
1910 propagator["rho0_770"], coefficient["<k10_1270|rho0_770>_0_mag"],
1911 coefficient["<k10_1270|rho0_770>_0_phase"] );
1912 addPartialWave( spinfactor["D2AP_A2VP_ks0pippim_ks0pim_0"], propagator["k10_1270"],
1913 propagator["kstarm_892"], coefficient["<k10_1270|kstarm_892>_0_mag"],
1914 coefficient["<k10_1270|kstarm_892>_0_phase"] );
1915 addPartialWave( spinfactor["D2AP_A2SP_ks0pippim_ks0pim_1"], propagator["k10_1270"],
1916 propagator["k0starm_1430"], coefficient["<k10_1270|k0starm_1430>_1_mag"],
1917 coefficient["<k10_1270|k0starm_1430>_1_phase"] );
1918 addPartialWave( spinfactor["D2AP_A2VP_ks0pimpi0_ks0pim_0"], propagator["k1m_1400"],
1919 propagator["kstarm_892"], spinfactor["D2AP_A2VP_ks0pimpi0_ks0pi0_0"],
1920 propagator["k1m_1400"], propagator["kstar0_892"], f_isoConj,
1921 coefficient["<k1m_1400|kstar_892>_0_mag"],
1922 coefficient["<k1m_1400|kstar_892>_0_phase"] );
1923 addPartialWave( spinfactor["D2AP_A2VP_ks0pippim_ks0pim_0"], propagator["k10_1400"],
1924 propagator["kstarm_892"], coefficient["<k10_1400|kstarm_892>_0_mag"],
1925 coefficient["<k10_1400|kstarm_892>_0_phase"] );
1926 addPartialWave( spinfactor["D2VV_ks0pi0_pippim_1"], propagator["kstar02_1410"],
1927 propagator["rho0_770"], coefficient["<kstar02_1410|rho0_770>_1_mag"],
1928 coefficient["<kstar02_1410|rho0_770>_1_phase"] );
1929 addPartialWave( spinfactor["D2VV_ks0pim_pippi0_0"], propagator["kstarm2_1410"],
1930 propagator["rhop_770"], coefficient["<kstarm2_1410|rhop_770>_0_mag"],
1931 coefficient["<kstarm2_1410|rhop_770>_0_phase"] );
1932 addPartialWave( spinfactor["D2VV_ks0pim_pippi0_1"], propagator["kstarm2_1410"],
1933 propagator["rhop_770"], coefficient["<kstarm2_1410|rhop_770>_1_mag"],
1934 coefficient["<kstarm2_1410|rhop_770>_1_phase"] );
1935 addPartialWave( spinfactor["D2VV_ks0pim_pippi0_2"], propagator["kstarm2_1410"],
1936 propagator["rhop_770"], coefficient["<kstarm2_1410|rhop_770>_2_mag"],
1937 coefficient["<kstarm2_1410|rhop_770>_2_phase"] );
1938 addPartialWave( spinfactor["D2VP_V2VP_ks0pimpi0_ks0pim_1"], propagator["kstarm3_1410"],
1939 propagator["kstarm_892"], spinfactor["D2VP_V2VP_ks0pimpi0_ks0pi0_1"],
1940 propagator["kstarm3_1410"], propagator["kstar0_892"], f_isoConj,
1941 coefficient["<kstarm3_1410|kstar_892>_1_mag"],
1942 coefficient["<kstarm3_1410|kstar_892>_1_phase"] );
1943 addPartialWave( spinfactor["D2VP_V2VP_ks0pimpi0_pimpi0_1"], propagator["kstarm3_1410"],
1944 propagator["rhom_770"], coefficient["<kstarm3_1410|rhom_770>_1_mag"],
1945 coefficient["<kstarm3_1410|rhom_770>_1_phase"] );
1946 addPartialWave( spinfactor["D2PP_P2VP_ks0pimpi0_ks0pim_1"], propagator["km_1460"],
1947 propagator["kstarm_892"], spinfactor["D2PP_P2VP_ks0pimpi0_ks0pi0_1"],
1948 propagator["km_1460"], propagator["kstar0_892"], f_isoConj,
1949 coefficient["<km_1460|kstar_892>_1_mag"],
1950 coefficient["<km_1460|kstar_892>_1_phase"] );
1951 addPartialWave( spinfactor["D2PP_P2SP_ks0pimpi0_ks0pim_0"], propagator["km_1460"],
1952 propagator["k0starm_1430"], spinfactor["D2PP_P2SP_ks0pimpi0_ks0pi0_0"],
1953 propagator["km_1460"], propagator["k0star0_1430"], f_isoConj,
1954 coefficient["<km_1460|k0star_1430>_0_mag"],
1955 coefficient["<km_1460|k0star_1430>_0_phase"] );
1956 addPartialWave( spinfactor["D2PP_P2VP_ks0pippim_ks0pim_1"], propagator["k0_1460"],
1957 propagator["kstarm_892"], coefficient["<k0_1460|kstarm_892>_1_mag"],
1958 coefficient["<k0_1460|kstarm_892>_1_phase"] );
1959 addPartialWave( spinfactor["D2AP_A2VP_ks0pimpi0_ks0pim_0"], propagator["k1m_1650"],
1960 propagator["kstarm_892"], spinfactor["D2AP_A2VP_ks0pimpi0_ks0pi0_0"],
1961 propagator["k1m_1650"], propagator["kstar0_892"], f_isoConj,
1962 coefficient["<k1m_1650|kstar_892>_0_mag"],
1963 coefficient["<k1m_1650|kstar_892>_0_phase"] );
1964 addPartialWave( spinfactor["D2AP_A2VP_ks0pippim_ks0pim_0"], propagator["k10_1650"],
1965 propagator["kstarm_892"], coefficient["<k10_1650|kstarm_892>_0_mag"],
1966 coefficient["<k10_1650|kstarm_892>_0_phase"] );
1967 addPartialWave( spinfactor["D2AP_A2VP_ks0pimpi0_pimpi0_0"], propagator["k1m_1650"],
1968 propagator["rhom_770"], coefficient["<k1m_1650|rhom_770>_0_mag"],
1969 coefficient["<k1m_1650|rhom_770>_0_phase"] );
1970 addPartialWave( spinfactor["D2VP_V2VP_ks0pippim_ks0pim_1"], propagator["kstar03_1680"],
1971 propagator["kstarm_892"], coefficient["<kstar03_1680|kstarm_892>_1_mag"],
1972 coefficient["<kstar03_1680|kstarm_892>_1_phase"] );
1973 addPartialWave( spinfactor["D2VP_V2VP_ks0pimpi0_ks0pim_1"], propagator["kstarm3_1680"],
1974 propagator["kstarm_892"], spinfactor["D2VP_V2VP_ks0pimpi0_ks0pi0_1"],
1975 propagator["kstarm3_1680"], propagator["kstar0_892"], f_isoConj,
1976 coefficient["<kstarm3_1680|kstar_892>_1_mag"],
1977 coefficient["<kstarm3_1680|kstar_892>_1_phase"] );
1978 addPartialWave( spinfactor["D2SS_ks0pi0_pippim_0"], propagator["k0star0_1430"],
1979 propagator["rho0_phsp"], coefficient["<k0star0_1430|rho0_phsp>_0_mag"],
1980 coefficient["<k0star0_1430|rho0_phsp>_0_phase"] );
1981 addPartialWave( spinfactor["D2PP_P2VP_ks0pippi0_pippi0_1"], propagator["k1p_phsp"],
1982 propagator["rhop_770"], coefficient["<k1p_phsp|rhop_770>_1_mag"],
1983 coefficient["<k1p_phsp|rhop_770>_1_phase"] );
1984 addPartialWave( spinfactor["D2PP_P2VP_ks0pimpi0_pimpi0_1"], propagator["k1m_phsp"],
1985 propagator["rhom_770"], coefficient["<k1m_phsp|rhom_770>_1_mag"],
1986 coefficient["<k1m_phsp|rhom_770>_1_phase"] );
1987 addPartialWave( spinfactor["D2AP_A2VP_ks0pimpi0_pimpi0_2"], propagator["k1m_phsp"],
1988 propagator["rhom_770"], coefficient["<k1m_phsp|rhom_770>_2_mag"],
1989 coefficient["<k1m_phsp|rhom_770>_2_phase"] );
1990 addPartialWave( spinfactor["D2SS_ks0pim_pippi0_0"], propagator["k0starm_1430"],
1991 propagator["rhop_phsp"], coefficient["<k0starm_1430|rhop_phsp>_0_mag"],
1992 coefficient["<k0starm_1430|rhop_phsp>_0_phase"] );
1993 addPartialWave( spinfactor["D2PP_P2SP_ks0pippim_pippim_0"], propagator["k0_1460"],
1994 propagator["sigma_500"], coefficient["<k0_1460|sigma_500>_0_mag"],
1995 coefficient["<k0_1460|sigma_500>_0_phase"] );
1996 addPartialWave( spinfactor["D2AP_A2SP_ks0pippim_pippim_1"], propagator["k10_1270"],
1997 propagator["f0_980"], coefficient["<k10_1270|f0_980>_1_mag"],
1998 coefficient["<k10_1270|f0_980>_1_phase"] );
1999 addPartialWave( spinfactor["D2AP_A2SP_ks0pippim_pippim_1"], propagator["k10_1400"],
2000 propagator["f0_1370"], coefficient["<k10_1400|f0_1370>_1_mag"],
2001 coefficient["<k10_1400|f0_1370>_1_phase"] );
2002 addPartialWave( spinfactor["D2VS_pippim_ks0pi0_1"], propagator["rho0_770"],
2003 propagator["k0star0_1430"], coefficient["<rho0_770|k0star0_1430>_1_mag"],
2004 coefficient["<rho0_770|k0star0_1430>_1_phase"] );
2005 addPartialWave( spinfactor["D2TP_T2VP_ks0pimpi0_ks0pim_2"], propagator["k2starm3_1430"],
2006 propagator["kstarm_892"], spinfactor["D2TP_T2VP_ks0pimpi0_ks0pi0_2"],
2007 propagator["k2starm3_1430"], propagator["kstar0_892"], f_isoConj,
2008 coefficient["<k2starm3_1430|kstar_892>_2_mag"],
2009 coefficient["<k2starm3_1430|kstar_892>_2_phase"] );
2010 addPartialWave( spinfactor["D2TT_ks0pi0_pippim_1"], propagator["k2star02_1430"],
2011 propagator["f2_1270"], coefficient["<k2star02_1430|f2_1270>_1_mag"],
2012 coefficient["<k2star02_1430|f2_1270>_1_phase"] );
2013 addPartialWave( spinfactor["D2VS_ks0pi0_pippim_1"], propagator["kstar0_892"],
2014 propagator["f0_980"], coefficient["<kstar0_892|f0_980>_1_mag"],
2015 coefficient["<kstar0_892|f0_980>_1_phase"] );
2016 addPartialWave( spinfactor["D2SS_ks0pim_pippi0_0"], propagator["kstarm_phsp"],
2017 propagator["rhop_phsp"], coefficient["<kstarm_phsp|rhop_phsp>_0_mag"],
2018 coefficient["<kstarm_phsp|rhop_phsp>_0_phase"] );
2019 addPartialWave( spinfactor["D2VS_pippim_ks0pi0_1"], propagator["rho0_770"],
2020 propagator["kstar0_phsp"], coefficient["<rho0_770|kstar0_phsp>_1_mag"],
2021 coefficient["<rho0_770|kstar0_phsp>_1_phase"] );
2022
2023 double Amps = abs2( totalAmp );
2024 return Amps;
2025}
2026
2027double EvtD0ToKSpipipi0::AmplitudeSquare() {
2028 _pd[0] = GetDaugMomCM( 0 );
2029 _pd[1] = GetDaugMomCM( 1 );
2030 _pd[2] = GetDaugMomCM( 2 );
2031 _pd[3] = GetDaugMomCM( 3 );
2032 double prob = amps();
2033 return ( prob <= 0 ) ? 1e-20 : prob;
2034}
p2
double p2[4]
Definition BabayagaNLO.cxx:41
p1
double p1[4]
Definition BabayagaNLO.cxx:41
P
double P(RecMdcKalTrack *trk)
Definition CalibEventSelect.cxx:84
mass
double mass
Definition CosmicGenerator.cxx:128
g1
TF1 * g1
Definition DataBase/tau_mode.c:58
imag
double imag(const EvtComplex &c)
Definition EvtComplex.hh:207
abs2
double abs2(const EvtComplex &c)
Definition EvtComplex.hh:193
exp
EvtComplex exp(const EvtComplex &c)
Definition EvtComplex.hh:209
operator+
std::vector< double > operator+(std::vector< double > &lhs, std::vector< double > &rhs)
Definition EvtD0ToKSpipipi0.cc:48
toString
string toString(const T &t)
Definition EvtD0ToKSpipipi0.cc:42
mass_Pion
const double mass_Pion
Definition EvtD0ToKSpipipi0.cc:1561
mpi
double mpi
Definition EvtPhokharaDef.hh:17
output
*******INTEGER m_nBinMax INTEGER m_NdiMax !No of bins in histogram for cell exploration division $ !Last vertex $ !Last active cell $ !Last cell in buffer $ !No of sampling when dividing cell $ !No of function total $ !Flag for random ceel for $ !Flag for type of for WtMax $ !Flag which decides whether vertices are included in the sampling $ entire domain is hyp !Maximum effective eevents per saves r n generator level $ !Flag for chat level in output
Definition FoamA.h:89
s
XmlRpcServer s
Definition HelloServer.cpp:11
sin
double sin(const BesAngle a)
Definition InstallArea/x86_64-el9-gcc13-dbg/include/MdcGeom/BesAngle.h:185
cos
double cos(const BesAngle a)
Definition InstallArea/x86_64-el9-gcc13-dbg/include/MdcGeom/BesAngle.h:187
one
const DifNumber one
q
****INTEGER imax DOUBLE PRECISION m_pi *DOUBLE PRECISION m_amfin DOUBLE PRECISION m_Chfin DOUBLE PRECISION m_Xenph DOUBLE PRECISION m_sinw2 DOUBLE PRECISION m_GFermi DOUBLE PRECISION m_MfinMin DOUBLE PRECISION m_ta2 INTEGER m_out INTEGER m_KeyFSR INTEGER m_KeyQCD *COMMON c_Semalib $ !copy of input $ !CMS energy $ !beam mass $ !final mass $ !beam charge $ !final charge $ !smallest final mass $ !Z mass $ !Z width $ !EW mixing angle $ !Gmu Fermi $ alphaQED at q
Definition KKsem.h:33
EvtD0ToKSpipipi0::~EvtD0ToKSpipipi0
virtual ~EvtD0ToKSpipipi0()
Definition EvtD0ToKSpipipi0.cc:58
EvtD0ToKSpipipi0::clone
EvtDecayBase * clone()
Definition EvtD0ToKSpipipi0.cc:62
EvtD0ToKSpipipi0::getName
void getName(std::string &name)
Definition EvtD0ToKSpipipi0.cc:60
EvtD0ToKSpipipi0::decay
void decay(EvtParticle *p)
Definition EvtD0ToKSpipipi0.cc:81
EvtDecayBase::checkSpinParent
void checkSpinParent(EvtSpinType::spintype sp)
Definition EvtDecayBase.cc:479
EvtDecayBase::setProbMax
void setProbMax(double prbmx)
Definition EvtDecayBase.cc:261
EvtDecayBase::checkNDaug
void checkNDaug(int d1, int d2=-1)
Definition EvtDecayBase.cc:465
EvtDecayBase::getDaugs
EvtId * getDaugs()
Definition EvtDecayBase.hh:68
EvtDecayBase::checkNArg
void checkNArg(int a1, int a2=-1, int a3=-1, int a4=-1)
Definition EvtDecayBase.cc:449
EvtDecayProb::setProb
void setProb(double prob)
Definition EvtDecayProb.hh:40
EvtParticle::getP4Lab
EvtVector4R getP4Lab()
Definition EvtParticle.cc:694
EvtParticle::getP4
const EvtVector4R & getP4() const
Definition EvtParticle.cc:108
EvtParticle::getDaug
EvtParticle * getDaug(int i)
Definition EvtParticle.cc:80
EvtParticle::initializePhaseSpace
double initializePhaseSpace(int numdaughter, EvtId *daughters, double poleSize=-1., int whichTwo1=0, int whichTwo2=1)
Definition EvtParticle.cc:1042
EvtVector4R::get
double get(int i) const
Definition EvtVector4R.hh:171
EvtVector4R::mass2
double mass2() const
Definition EvtVector4R.hh:110
t
int t()
Definition t.c:1