41 model_name=
"DTopipi0pi0";
113 int GG[4][4] = { {1,0,0,0}, {0,-1,0,0}, {0,0,-1,0}, {0,0,0,-1} };
114 for (
int i=0; i<4; i++) {
115 for (
int j=0; j<4; j++) {
132 double P1[4], P2[4], P3[4];
133 P1[0] = D1.
get(0); P1[1] = D1.
get(1); P1[2] = D1.
get(2); P1[3] = D1.
get(3);
134 P2[0] = D2.
get(0); P2[1] = D2.
get(1); P2[2] = D2.
get(2); P2[3] = D2.
get(3);
135 P3[0] = D3.
get(0); P3[1] = D3.
get(1); P3[2] = D3.
get(2); P3[3] = D3.
get(3);
139 int spin[4]={2,1,1,0};
141 double r0[4]={3,3,3,3};
142 double r1[4]={5,5,5,5};
144 calEva(P1, P2, P3, mass, width, rho, phi, g0, spin, modetype, nstates, value,r0,r1);
151void EvtDTopipi0pi0::Com_Multi(
double a1[2],
double a2[2],
double res[2])
153 res[0] = a1[0]*a2[0]-a1[1]*a2[1];
154 res[1] = a1[1]*a2[0]+a1[0]*a2[1];
156void EvtDTopipi0pi0::Com_Divide(
double a1[2],
double a2[2],
double res[2])
158 double tmp = a2[0]*a2[0]+a2[1]*a2[1];
159 res[0] = (a1[0]*a2[0]+a1[1]*a2[1])/tmp;
160 res[1] = (a1[1]*a2[0]-a1[0]*a2[1])/tmp;
163double EvtDTopipi0pi0::CalRho4pi(double_t
s)
166 return sqrt((
s-16.*mass_Pion*mass_Pion)/
s);
169 double_t s0 = 1.2274+0.00370909/(
s*
s) - (0.111203)/(
s) - 6.39017*
s +16.8358*
s*
s - 21.8845*
s*
s*
s + 11.3153*
s*
s*
s*
s;
170 double_t gam = s0*sqrt(1.0-(16.0*mass_Pion*mass_Pion));
177double EvtDTopipi0pi0::SCADot(
double a1[4],
double a2[4])
179 double _cal = a1[0]*a2[0]-a1[1]*a2[1]-a1[2]*a2[2]-a1[3]*a2[3];
182double EvtDTopipi0pi0::barrier(
int l,
double sa,
double sb,
double sc,
double r,
double mass)
184 double q = (sa+sb-sc)*(sa+sb-sc)/(4*sa)-sb;
190 double q0 = (sa0+sb-sc)*(sa0+sb-sc)/(4*sa0)-sb;
191 if(q0 < 0) q0 = 1e-16;
195 if(l == 1) F = sqrt((1+z0)/(1+z));
196 if(l == 2) F = sqrt((9+3*z0+z0*z0)/(9+3*z+z*z));
200void EvtDTopipi0pi0::calt1(
double daug1[4],
double daug2[4],
double t1[4])
204 for(
int i=0; i<4; i++) {
205 pa[i] = daug1[i] + daug2[i];
206 qa[i] = daug1[i] - daug2[i];
211 for(
int i=0; i<4; i++) {
212 t1[i] = qa[i] - tmp*pa[i];
215void EvtDTopipi0pi0::calt2(
double daug1[4],
double daug2[4],
double t2[4][4])
219 calt1(daug1,daug2,t1);
220 r = SCADot(t1,t1)/3.0;
221 for(
int i=0; i<4; i++) {
222 pa[i] = daug1[i] + daug2[i];
225 for(
int i=0; i<4; i++) {
226 for(
int j=0; j<4; j++) {
227 t2[i][j] = t1[i]*t1[j] - r*(G[i][j]-pa[i]*pa[j]/p);
232void EvtDTopipi0pi0::propagator(
double mass2,
double mass,
double width,
double sx,
double prop[2])
239 Com_Divide(a,
b,prop);
241double EvtDTopipi0pi0::wid(
double mass2,
double mass,
double sa,
double sb,
double sc,
double r2,
int l)
246 double tmp1 = sa+tmp;
247 double q = 0.25*tmp1*tmp1/sa-sb;
249 double tmp2 = mass2+tmp;
250 double q0 = 0.25*tmp2*tmp2/mass2-sb;
255 if(l == 0) {widm = sqrt(
t)*
mass/m;}
256 else if(l == 1) {widm =
t*sqrt(
t)*
mass/m*(1+z0)/(1+z);}
257 else if(l == 2) {widm =
t*
t*sqrt(
t)*
mass/m*(9+3*z0+z0*z0)/(9+3*z+z*z);}
260double EvtDTopipi0pi0::widl1(
double mass2,
double mass,
double sa,
double sb,
double sc,
double r2)
265 double tmp1 = sa+tmp;
266 double q = 0.25*tmp1*tmp1/sa-sb;
268 double tmp2 = mass2+tmp;
269 double q0 = 0.25*tmp2*tmp2/mass2-sb;
273 double F = (1+z0)/(1+z);
275 widm =
t*sqrt(
t)*
mass/m*F;
278void EvtDTopipi0pi0::propagatorRBW(
double mass2,
double mass,
double width,
double sa,
double sb,
double sc,
double r2,
int l,
double prop[2])
285 b[1] = -
mass*width*wid(mass2,mass,sa,sb,sc,r2,l);
286 Com_Divide(a,
b,prop);
289void EvtDTopipi0pi0::propagatorFlatte(
double mass,
double width,
double sa,
double sb,
double sc,
int r,
double prop[2]){
291 double rhoab[2], rhoKsK[2];
292 q = 0.25*(sa+sb-sc)*(sa+sb-sc)/sa-sb;
293 if(r == 0) qKsK = 0.25*sa - 0.49368*0.49368;
294 if(r == 1) qKsK = 0.25*(sa + 3.899750596e-03)*(sa + 3.899750596e-03)/sa - 0.497614*0.497614;
296 rhoab[0] = 2*sqrt(
q/sa);
301 rhoab[1] = 2*sqrt(-
q/sa);
304 rhoKsK[0] = 2*sqrt(qKsK/sa);
309 rhoKsK[1] = 2*sqrt(-qKsK/sa);
314 b[0] =
mass*
mass - sa + 0.341*rhoab[1] + 0.892*0.341*rhoKsK[1];
315 b[1] = - (0.341*rhoab[0] + 0.892*0.341*rhoKsK[0]);
316 Com_Divide(a,
b,prop);
319void EvtDTopipi0pi0::rhoab(
double sa,
double sb,
double sc,
double res[2]) {
320 double tmp = sa+sb-sc;
321 double q = 0.25*tmp*tmp/sa-sb;
323 res[0]=2.0*sqrt(
q/sa);
327 res[1]=2.0*sqrt(-
q/sa);
330void EvtDTopipi0pi0::rho4Pi(
double sa,
double res[2]) {
331 double temp = 1.0-0.3116765584/sa;
333 res[0]=sqrt(temp)/(1.0+
exp(9.8-3.5*sa));
337 res[1]=sqrt(-temp)/(1.0+
exp(9.8-3.5*sa));
342void EvtDTopipi0pi0::propagatorsigma500(
double sa,
double sb,
double sc,
double prop[2]) {
343 double f = 0.5843+1.6663*sa;
344 const double M = 0.9264;
345 const double mass2 = 0.85821696;
346 const double mpi2d2 = 0.00973989245;
347 double g1 =
f*(sa-mpi2d2)/(mass2-mpi2d2)*
exp((mass2-sa)/1.082);
348 double rho1s[2], rho1M[2], rho2s[2], rho2M[2], rho1[2], rho2[2];
349 rhoab(sa,sb,sc,rho1s);
350 rhoab(mass2,sb,sc,rho1M);
353 Com_Divide(rho1s,rho1M,rho1);
354 Com_Divide(rho2s,rho2M,rho2);
358 b[0] = mass2-sa+M*(
g1*rho1[1]+0.0024*rho2[1]);
359 b[1] = -M*(
g1*rho1[0]+0.0024*rho2[0]);
360 Com_Divide(a,
b,prop);
362void EvtDTopipi0pi0::Flatte_rhoab(
double sa,
double sb,
double rho[2])
364 double q = 1.0-(4*sb/sa);
376void EvtDTopipi0pi0::propagator980(
double mass,
double sx,
double *sb,
double prop[2])
379 double gpipi1 = 2.0/3.0*0.165;
380 double gpipi2 = 1.0/3.0*0.165;
381 double gKK1 = 0.5*0.69465;
382 double gKK2 = 0.5*0.69465;
384 double unit[2]={1.0};
387 Flatte_rhoab(sx,sb[0],rho1);
389 Flatte_rhoab(sx,sb[1],rho2);
391 Flatte_rhoab(sx,sb[2],rho3);
393 Flatte_rhoab(sx,sb[3],rho4);
395 double tmp1[2]={gpipi1,0};
397 double tmp2[2]={gpipi2,0};
399 double tmp3[2]={gKK1,0};
401 double tmp4[2]={gKK2,0};
404 Com_Multi(tmp1,rho1,tmp11);
405 Com_Multi(tmp2,rho2,tmp22);
406 Com_Multi(tmp3,rho3,tmp33);
407 Com_Multi(tmp4,rho4,tmp44);
409 double tmp5[2]={tmp11[0]+tmp22[0]+tmp33[0]+tmp44[0],tmp11[1]+tmp22[1]+tmp33[1]+tmp44[1]};
411 Com_Multi(tmp5, ci,tmp51);
412 double tmp6[2]={
mass*
mass-sx-tmp51[0], -1.0*tmp51[1]};
414 Com_Divide(
unit,tmp6, prop);
419void EvtDTopipi0pi0::propagatorGS(
double mass2,
double mass,
double width,
double sa,
double sb,
double sc,
double r2,
double prop[2])
423 double tmp1 = sa+tmp;
424 double q2 = 0.25*tmp1*tmp1/sa-sb;
427 double tmp2 = mass2+tmp;
428 double q02 = 0.25*tmp2*tmp2/mass2-sb;
429 if(q02<0) q02 = 1e-16;
432 double q0 = sqrt(q02);
435 double tmp3 = log(mass+2*q0)+1.2926305904;
437 double h = GS1*
q/m*(log(m+2*
q)+1.2926305904);
438 double h0 = GS1*q0/
mass*tmp3;
439 double dh = h0*(0.125/q02-0.5/mass2)+GS3/mass2;
440 double d = GS2/q02*tmp3+GS3*
mass/q0-GS4*
mass/q03;
441 double f = mass2/q03*(q2*(h-h0)+(mass2-sa)*q02*dh);
443 a[0] = 1.0+d*width/
mass;
445 b[0] = mass2-sa+width*
f;
446 b[1] = -
mass*width*widl1(mass2,mass,sa,sb,sc,r2);
447 Com_Divide(a,
b,prop);
450double EvtDTopipi0pi0::DDalitz(
double P1[4],
double P2[4],
double P3[4],
int Ang,
double mass,
double rRES0,
double rRES1){
452 double temp_PDF, v_re;
455 double B[2], s1, s2, s3, sR, sD;
456 for(
int i=0; i<4; i++){
457 pR[i] = P1[i] + P2[i];
458 pD[i] = pR[i] + P3[i];
466 for(
int i=0; i!=4; i++){
467 for(
int j=0; j!=4; j++){
469 if(i==0) G[i][j] = 1;
481 B[0] = barrier(1,sR,s1,s2,rRES0,mass);
482 B[1] = barrier(1,sD,sR,s3,rRES1,1.86966);
487 for(
int i=0; i<4; i++){
488 temp_PDF += t1[i]*T1[i]*G[i][i];
492 B[0] = barrier(2,sR,s1,s2,rRES0,mass);
493 B[1] = barrier(2,sD,sR,s3,rRES1,1.86966);
494 double t2[4][4], T2[4][4];
498 for(
int i=0; i<4; i++){
499 for(
int j=0; j<4; j++){
500 temp_PDF += t2[i][j]*T2[j][i]*G[i][i]*G[j][j];
504 v_re = temp_PDF*
B[0]*
B[1];
508TComplex EvtDTopipi0pi0::ResonanceSkm(Double_t&
m2)
510 Double_t g11 = 0.22889, g12 = -0.55377, g13 = 0.00, g14 = -0.39899, g15 = -0.34639;
511 Double_t g21 = 0.94128, g22 = 0.55095, g23 = 0.00, g24 = 0.39065, g25 = 0.31503;
512 Double_t g31 = 0.36856, g32 = 0.23888, g33 = 0.55639, g34 = 0.18340, g35 = 0.18681;
513 Double_t g41 = 0.33650, g42 = 0.40907, g43 = 0.85679, g44 = 0.19906, g45 = -0.00984;
514 Double_t g51 = 0.18171, g52 = -0.17558, g53 = -0.79658, g54 = -0.00355, g55 = 0.22358;
516 Double_t fs11 = 0.23399, fs12 = 0.15044, fs13 =-0.20545, fs14 = 0.32825, fs15 = 0.35412;
517 Double_t m_1 = 0.65100, m_2 = 1.20360, m_3 = 1.55817, m_4 = 1.21000, m_5 = 1.82206;
518 Double_t ss0 = -3.92637;
519 Double_t sp0 = -0.07;
520 double_t sA0 = -0.15;
523 double_t km11 = (g11*g11/(m_1*m_1-
m2)+g21*g21/(m_2*m_2-
m2)+g31*g31/(m_3*m_3-
m2)+g41*g41/(m_4*m_4-
m2)+g51*g51/(m_5*m_5-
m2)+fs11*(1-ss0)/(
m2-ss0))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
524 double_t km12 = (g11*g12/(m_1*m_1-
m2)+g21*g22/(m_2*m_2-
m2)+g31*g32/(m_3*m_3-
m2)+g41*g42/(m_4*m_4-
m2)+g51*g52/(m_5*m_5-
m2)+fs12*(1-ss0)/(
m2-ss0))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
525 double_t km13 = (g11*g13/(m_1*m_1-
m2)+g21*g23/(m_2*m_2-
m2)+g31*g33/(m_3*m_3-
m2)+g41*g43/(m_4*m_4-
m2)+g51*g53/(m_5*m_5-
m2)+fs13*(1-ss0)/(
m2-ss0))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
526 double_t km14 = (g11*g14/(m_1*m_1-
m2)+g21*g24/(m_2*m_2-
m2)+g31*g34/(m_3*m_3-
m2)+g41*g44/(m_4*m_4-
m2)+g51*g54/(m_5*m_5-
m2)+fs14*(1-ss0)/(
m2-ss0))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
527 double_t km15 = (g11*g15/(m_1*m_1-
m2)+g21*g25/(m_2*m_2-
m2)+g31*g35/(m_3*m_3-
m2)+g41*g45/(m_4*m_4-
m2)+g51*g55/(m_5*m_5-
m2)+fs15*(1-ss0)/(
m2-ss0))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
529 double_t km21 = km12, km31 = km13, km41 = km14, km51 = km15;
530 double_t km23 = (g12*g13/(m_1*m_1-
m2)+g22*g23/(m_2*m_2-
m2)+g32*g33/(m_3*m_3-
m2)+g42*g43/(m_4*m_4-
m2)+g52*g53/(m_5*m_5-
m2))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
531 double_t km24 = (g12*g14/(m_1*m_1-
m2)+g22*g24/(m_2*m_2-
m2)+g32*g34/(m_3*m_3-
m2)+g42*g44/(m_4*m_4-
m2)+g52*g54/(m_5*m_5-
m2))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
532 double_t km25 = (g12*g15/(m_1*m_1-
m2)+g22*g25/(m_2*m_2-
m2)+g32*g35/(m_3*m_3-
m2)+g42*g45/(m_4*m_4-
m2)+g52*g55/(m_5*m_5-
m2))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
533 double_t km32 = km23, km42 = km24, km52 = km25;
534 double_t km34 = (g13*g14/(m_1*m_1-
m2)+g23*g24/(m_2*m_2-
m2)+g33*g34/(m_3*m_3-
m2)+g43*g44/(m_4*m_4-
m2)+g53*g54/(m_5*m_5-
m2))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
535 double_t km35 = (g13*g15/(m_1*m_1-
m2)+g23*g25/(m_2*m_2-
m2)+g33*g35/(m_3*m_3-
m2)+g43*g45/(m_4*m_4-
m2)+g53*g55/(m_5*m_5-
m2))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
536 double_t km43 = km34, km53 = km35;
537 double_t km45 = (g14*g15/(m_1*m_1-
m2)+g24*g25/(m_2*m_2-
m2)+g34*g35/(m_3*m_3-
m2)+g44*g45/(m_4*m_4-
m2)+g54*g55/(m_5*m_5-
m2))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
538 double_t km54 = km45;
539 double_t km22 = (g12*g12/(m_1*m_1-
m2)+g22*g22/(m_2*m_2-
m2)+g32*g32/(m_3*m_3-
m2)+g42*g42/(m_4*m_4-
m2)+g52*g52/(m_5*m_5-
m2))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
540 double_t km33 = (g13*g13/(m_1*m_1-
m2)+g23*g23/(m_2*m_2-
m2)+g33*g33/(m_3*m_3-
m2)+g43*g43/(m_4*m_4-
m2)+g53*g53/(m_5*m_5-
m2))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
541 double_t km44 = (g14*g14/(m_1*m_1-
m2)+g24*g24/(m_2*m_2-
m2)+g34*g34/(m_3*m_3-
m2)+g44*g44/(m_4*m_4-
m2)+g54*g54/(m_5*m_5-
m2))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
542 double_t km55 = (g15*g15/(m_1*m_1-
m2)+g25*g25/(m_2*m_2-
m2)+g35*g35/(m_3*m_3-
m2)+g45*g45/(m_4*m_4-
m2)+g55*g55/(m_5*m_5-
m2))*(1-sA0)/(
m2-sA0)*(
m2-sA*mass_Pion*mass_Pion*0.5);
559 TComplex P1 = fp11*(1-sp0)/(
m2-sp0)+beta1*g11/(m_1*m_1-
m2)+beta2*g21/(m_2*m_2-
m2)+beta3*g31/(m_3*m_3-
m2)+beta4*g41/(m_4*m_4-
m2)+beta5*g51/(m_5*m_5-
m2);
560 TComplex P2 = fp12*(1-sp0)/(
m2-sp0)+beta1*g12/(m_1*m_1-
m2)+beta2*g22/(m_2*m_2-
m2)+beta3*g32/(m_3*m_3-
m2)+beta4*g42/(m_4*m_4-
m2)+beta5*g52/(m_5*m_5-
m2);
561 TComplex P3 = fp13*(1-sp0)/(
m2-sp0)+beta1*g13/(m_1*m_1-
m2)+beta2*g23/(m_2*m_2-
m2)+beta3*g33/(m_3*m_3-
m2)+beta4*g43/(m_4*m_4-
m2)+beta5*g53/(m_5*m_5-
m2);
562 TComplex P4 = fp14*(1-sp0)/(
m2-sp0)+beta1*g14/(m_1*m_1-
m2)+beta2*g24/(m_2*m_2-
m2)+beta3*g34/(m_3*m_3-
m2)+beta4*g44/(m_4*m_4-
m2)+beta5*g54/(m_5*m_5-
m2);
563 TComplex P5 = fp15*(1-sp0)/(
m2-sp0)+beta1*g15/(m_1*m_1-
m2)+beta2*g25/(m_2*m_2-
m2)+beta3*g35/(m_3*m_3-
m2)+beta4*g45/(m_4*m_4-
m2)+beta5*g55/(m_5*m_5-
m2);
565 TMatrixD MI(5,5), MA(5,5), MA_invt(5,5), MB(5,5), KM(5,5), RhoA(5,5), RhoB(5,5), MRe(5,5), MIm(5,5);
569 TMatrixDRow(KM,0)(0) = km11;TMatrixDRow(KM,1)(0) = km21;TMatrixDRow(KM,2)(0) = km31;TMatrixDRow(KM,3)(0) = km41;TMatrixDRow(KM,4)(0)= km51;
570 TMatrixDRow(KM,0)(1) = km12;TMatrixDRow(KM,1)(1) = km22;TMatrixDRow(KM,2)(1) = km32;TMatrixDRow(KM,3)(1) = km42;TMatrixDRow(KM,4)(1)= km52;
571 TMatrixDRow(KM,0)(2) = km13;TMatrixDRow(KM,1)(2) = km23;TMatrixDRow(KM,2)(2) = km33;TMatrixDRow(KM,3)(2) = km43;TMatrixDRow(KM,4)(2)= km53;
572 TMatrixDRow(KM,0)(3) = km14;TMatrixDRow(KM,1)(3) = km24;TMatrixDRow(KM,2)(3) = km34;TMatrixDRow(KM,3)(3) = km44;TMatrixDRow(KM,4)(3)= km54;
573 TMatrixDRow(KM,0)(4) = km15;TMatrixDRow(KM,1)(4) = km25;TMatrixDRow(KM,2)(4) = km35;TMatrixDRow(KM,3)(4) = km45;TMatrixDRow(KM,4)(4)= km55;
575 if ( (1.0-4.0*mass_Pion*mass_Pion/
m2) > 0) {
576 TMatrixDRow(RhoA,0)(0) = sqrt(1.0-4.0*mass_Pion*mass_Pion/
m2);
577 TMatrixDRow(RhoB,0)(0) = 0.0;
580 TMatrixDRow(RhoA,0)(0) = 0.0;
581 TMatrixDRow(RhoB,0)(0) = sqrt(4.0*mass_Pion*mass_Pion/
m2-1.0);
584 if ( (1.0-4.0*mass_Kaon*mass_Kaon/
m2) > 0) {
585 TMatrixDRow(RhoA,1)(1) = sqrt(1.0-4.0*mass_Kaon*mass_Kaon/
m2);
586 TMatrixDRow(RhoB,1)(1) = 0.0;
589 TMatrixDRow(RhoA,1)(1) = 0.0;
590 TMatrixDRow(RhoB,1)(1) = sqrt(4.0*mass_Kaon*mass_Kaon/
m2-1.0);
593 TMatrixDRow(RhoA,2)(2) = CalRho4pi(
m2);
594 TMatrixDRow(RhoB,2)(2) = 0.0;
595 if ( (1.0-4.0*mass_Eta*mass_Eta/
m2) > 0) {
596 TMatrixDRow(RhoA,3)(3) = sqrt(1.0-4.0*mass_Eta*mass_Eta/
m2);
597 TMatrixDRow(RhoB,3)(3) = 0.0;
599 TMatrixDRow(RhoA,3)(3) = 0.0;
600 TMatrixDRow(RhoB,3)(3) = sqrt(4.0*mass_Eta*mass_Eta/
m2-1.0);
604 if ( (1.0-(mass_Eta + mass_Etap)*(mass_Eta + mass_Etap)/
m2) > 0) {
605 TMatrixDRow(RhoA,4)(4) = sqrt(1.0-(mass_Eta + mass_Etap)*(mass_Eta + mass_Etap)/
m2);
606 TMatrixDRow(RhoB,4)(4) = 0.0;
608 TMatrixDRow(RhoA,4)(4) = 0.0;
609 TMatrixDRow(RhoB,4)(4) = sqrt((mass_Eta + mass_Etap)*(mass_Eta + mass_Etap)/
m2-1.0);
618 MRe = MA+MB*MA_invt*MB;
620 MIm = MA_invt*MB*MRe;
626 amp = (ciR*TMatrixDRow(MRe,0)(0)-ciM*TMatrixDRow(MIm,0)(0))*P1+
627 (ciR*TMatrixDRow(MRe,0)(1)-ciM*TMatrixDRow(MIm,0)(1))*P2+
628 (ciR*TMatrixDRow(MRe,0)(2)-ciM*TMatrixDRow(MIm,0)(2))*P3+
629 (ciR*TMatrixDRow(MRe,0)(3)-ciM*TMatrixDRow(MIm,0)(3))*P4+
630 (ciR*TMatrixDRow(MRe,0)(4)-ciM*TMatrixDRow(MIm,0)(4))*P5;
635void EvtDTopipi0pi0::calEva(
double* Pic,
double* Pi01,
double* Pi02,
double *mass1,
double *width1,
double *amp,
double *phase,
int* g0,
int* spin,
int* modetype,
int nstates,
double & Result ,
double*r0 ,
double*r1)
637 double P12[4], P23[4], P13[4];
639 double cof[2], amp_PDF[2], PDF[2];
640 double spi01, spi02, scpi;
644 double Realpipis, Imagpipis;
645 for(
int i=0; i<4; i++){
646 P23[i] = Pi01[i] + Pi02[i];
647 P12[i] = Pic[i] + Pi01[i];
648 P13[i] = Pic[i] + Pi02[i];
650 scpi = SCADot(Pic,Pic);
651 spi01 = SCADot(Pi01,Pi01);
652 spi02 = SCADot(Pi02,Pi02);
653 s23 = SCADot(P23,P23);
654 s12 = SCADot(P12,P12);
655 s13 = SCADot(P13,P13);
661 double pro[2], temp_PDF, pro1[2], temp_PDF1, pro2[2], temp_PDF2, amp_tmp[2];
669 for(
int i=0; i<nstates; i++) {
672 mass1sq = mass1[i]*mass1[i];
673 cof[0] = amp[i]*
cos(phase[i]);
674 cof[1] = amp[i]*
sin(phase[i]);
676 if(modetype[i] == 23){
677 temp_PDF = DDalitz(Pi01, Pi02, Pic, spin[i], mass1[i], r0[i],r1[i]);
678 if(g0[i]==2) propagatorsigma500(
s23,spi01,spi02,pro);
679 if(g0[i]==1) propagatorRBW(mass1sq, mass1[i],width1[i],
s23,spi01,spi02,r0[i]*r0[i],spin[i],pro);
680 if(g0[i]==3) propagator980(mass1[i],
s23,spi012,pro);
682 tmpswave = ResonanceSkm(
s23);
683 Realpipis = tmpswave.Re();
684 Imagpipis = tmpswave.Im();
685 amp_tmp[0] = temp_PDF*Realpipis;
686 amp_tmp[1] = temp_PDF*Imagpipis;
692 if(g0[i]!=6) amp_tmp[0] = temp_PDF*pro[0];
693 if(g0[i]!=6) amp_tmp[1] = temp_PDF*pro[1];
695 if(modetype[i] == 12){
696 temp_PDF1 = DDalitz(Pic, Pi01, Pi02, spin[i], mass1[i], r0[i], r1[i]);
697 if(g0[i]==1) propagatorRBW(mass1sq, mass1[i],width1[i],
s12,scpi,spi01,r0[i]*r0[i],spin[i],pro1);
698 if(g0[i]==2) propagatorGS(mass1sq,mass1[i],width1[i],
s12,scpi,spi01,r0[i]*r0[i],pro1);
704 temp_PDF2 = DDalitz(Pic, Pi02, Pi01, spin[i], mass1[i], r0[i], r1[i]);
705 if(g0[i]==1) propagatorRBW(mass1sq, mass1[i],width1[i],s13,scpi,spi02,r0[i]*r0[i],spin[i],pro2);
706 if(g0[i]==2) propagatorGS(mass1sq,mass1[i],width1[i],s13,scpi,spi02,r0[i]*r0[i],pro2);
711 amp_tmp[0] = temp_PDF1*pro1[0] + temp_PDF2*pro2[0];
713 amp_tmp[1] = temp_PDF1*pro1[1] + temp_PDF2*pro2[1];
715 Com_Multi(amp_tmp,cof,amp_PDF);
716 PDF[0] += amp_PDF[0];
717 PDF[1] += amp_PDF[1];
719 double value = PDF[0]*PDF[0] + PDF[1]*PDF[1];
720 if(value <=0) value = 1e-20;
double sin(const BesAngle a)
double cos(const BesAngle a)
TFile f("ana_bhabha660a_dqa_mcPat_zy_old.root")
EvtComplex exp(const EvtComplex &c)
*******INTEGER m_nBinMax INTEGER m_NdiMax !No of bins in histogram for cell exploration division $ !Last vertex $ !Last active cell $ !Last cell in buffer $ !No of sampling when dividing cell $ !No of function total $ !Flag for random ceel for $ !Flag for type of for WtMax $ !Flag which decides whether vertices are included in the sampling $ entire domain is hyp !Maximum effective eevents per saves r n generator level $ !Flag for chat level in !Latex Output unit
****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
virtual ~EvtDTopipi0pi0()
void decay(EvtParticle *p)
void getName(std::string &name)
void checkSpinDaughter(int d1, EvtSpinType::spintype sp)
void checkSpinParent(EvtSpinType::spintype sp)
void setProbMax(double prbmx)
void checkNDaug(int d1, int d2=-1)
void checkNArg(int a1, int a2=-1, int a3=-1, int a4=-1)
void setProb(double prob)
const EvtVector4R & getP4() const
EvtParticle * getDaug(int i)
double initializePhaseSpace(int numdaughter, EvtId *daughters, double poleSize=-1., int whichTwo1=0, int whichTwo2=1)
double double double double * s12
double double double double double * s23