10#include "std_make_vector.h"
11#include "gslpp_function_adapter.h"
12#include <gsl/gsl_sf_zeta.h>
13#include <boost/bind/bind.hpp>
15#include <TFitResult.h>
16#include <gsl/gsl_sf_gegenbauer.h>
17#include <gsl/gsl_sf_expint.h>
19using namespace boost::placeholders;
22: mySM(SM_i), ig(ROOT::Math::Integration::kADAPTIVESINGULAR, ROOT::Math::Integration::kGAUSS51), wf()
33 checkcache_int_tau =
false;
34 checkcache_int_mu =
false;
35 checkcache_int_el =
false;
37 double max_double = std::numeric_limits<double>::max();
39 fplusz0_cache = max_double;
40 rho1to2_cache = max_double;
41 N_0_cache = max_double;
42 alpha_0_cache = max_double;
43 alpha_p_cache = max_double;
44 beta_0_cache = max_double;
45 beta_p_cache = max_double;
46 gamma_0_cache = max_double;
47 gamma_p_cache = max_double;
48 af0_1_cache = max_double;
49 af0_2_cache = max_double;
50 afplus_0_cache = max_double;
51 afplus_1_cache = max_double;
52 afplus_2_cache = max_double;
55 af0_3_cache = max_double;
56 afplus_3_cache = max_double;
59 afplus1_cache = max_double;
60 afplus2_cache = max_double;
61 afplus3_cache = max_double;
62 afplus4_cache = max_double;
63 afplus5_cache = max_double;
64 afplus6_cache = max_double;
65 afplus7_cache = max_double;
66 afplus8_cache = max_double;
67 afplus9_cache = max_double;
68 afplus10_cache = max_double;
69 afzero1_cache = max_double;
70 afzero2_cache = max_double;
71 afzero3_cache = max_double;
72 afzero4_cache = max_double;
73 afzero5_cache = max_double;
74 afzero6_cache = max_double;
75 afzero7_cache = max_double;
76 afzero8_cache = max_double;
77 afzero9_cache = max_double;
78 afzero10_cache = max_double;
79 bfplus1_cache = max_double;
80 bfplus2_cache = max_double;
81 bfplus3_cache = max_double;
82 bfplus4_cache = max_double;
83 bfplus5_cache = max_double;
84 bfplus6_cache = max_double;
85 bfplus7_cache = max_double;
86 bfplus8_cache = max_double;
87 bfplus9_cache = max_double;
88 bfplus10_cache = max_double;
89 bfzero1_cache = max_double;
90 bfzero2_cache = max_double;
91 bfzero3_cache = max_double;
92 bfzero4_cache = max_double;
93 bfzero5_cache = max_double;
94 bfzero6_cache = max_double;
95 bfzero7_cache = max_double;
96 bfzero8_cache = max_double;
97 bfzero9_cache = max_double;
98 bfzero10_cache = max_double;
100 CS_cache = max_double;
101 CSp_cache = max_double;
102 CP_cache = max_double;
103 CPp_cache = max_double;
104 CV_cache = max_double;
105 CVp_cache = max_double;
106 CA_cache = max_double;
107 CAp_cache = max_double;
108 CT_cache = max_double;
109 CTp_cache = max_double;
111 ig.SetRelTolerance(1.E-6);
112 ig.SetAbsTolerance(1.E-6);
126 if (
CLNflag +
BGLflag +
DMflag !=
true)
throw std::runtime_error(
"MPlnu: Set only one among CLNflag, BGLflag, DMflag to true");
131 <<
"fplusz0" <<
"rho1to2"
132 <<
"N_0" <<
"alpha_0" <<
"alpha_p" <<
"beta_0" <<
"beta_p" <<
"gamma_0" <<
"gamma_p";
137 <<
"af0_1" <<
"af0_2" <<
"afplus_0" <<
"afplus_1" <<
"afplus_2"
139 <<
"af0_3" <<
"afplus_3"
141 <<
"mBc1m_1" <<
"mBc1m_2" <<
"mBc1m_3" <<
"mBc1m_4"
142 <<
"mBc0p_1" <<
"mBc0p_2" <<
"chitildeT" <<
"chiL" <<
"n_I";
147 <<
"afplus1" <<
"afplus2" <<
"afplus3" <<
"afplus4" <<
"afplus5"
148 <<
"afplus6" <<
"afplus7" <<
"afplus8" <<
"afplus9" <<
"afplus10"
149 <<
"afzero1" <<
"afzero2" <<
"afzero3" <<
"afzero4" <<
"afzero5"
150 <<
"afzero6" <<
"afzero7" <<
"afzero8" <<
"afzero9" <<
"afzero10"
151 <<
"bfplus1" <<
"bfplus2" <<
"bfplus3" <<
"bfplus4" <<
"bfplus5"
152 <<
"bfplus6" <<
"bfplus7" <<
"bfplus8" <<
"bfplus9" <<
"bfplus10"
153 <<
"bfzero1" <<
"bfzero2" <<
"bfzero3" <<
"bfzero4" <<
"bfzero5"
154 <<
"bfzero6" <<
"bfzero7" <<
"bfzero8" <<
"bfzero9" <<
"bfzero10";
157 std::stringstream out;
159 throw std::runtime_error(
"MPlnu: vector " + out.str() +
" not implemented");
167void MPlnu::updateParameters()
183 amplsq_factor = 1./(64.*M_PI*M_PI*M_PI*
MM*
MM*
MM);
188 gslpp::complex norm = 4./sqrt(2.);
190 CV = (*(allcoeff_bclnu[
LO]))(0)/norm*(1.+ale_mub/M_PI*log(
mySM.
getMz()/
mu_b))/2.;
192 CVp = (*(allcoeff_bclnu[
LO]))(1)/norm/2.;
194 CS = (*(allcoeff_bclnu[
LO]))(2)/norm/2.;
195 CSp = (*(allcoeff_bclnu[
LO]))(3)/norm/2.;
200 CT = (*(allcoeff_bclnu[
LO]))(4)/norm/2.;
319 std::stringstream out;
321 throw std::runtime_error(
"MPlnu: vector " + out.str() +
" not implemented");
324 if ((fplusz0 != fplusz0_cache) || (rho1to2 != rho1to2_cache)
325 || (N_0 != N_0_cache)
326 || (alpha_0 != alpha_0_cache) || (alpha_p != alpha_p_cache)
327 || (beta_0 != beta_0_cache) || (beta_p != beta_p_cache)
328 || (gamma_0 != gamma_0_cache) || (gamma_p != gamma_p_cache)
329 || (af0_1 != af0_1_cache) || (af0_2 != af0_2_cache)
330 || (afplus_0 != afplus_0_cache) || (afplus_1 != afplus_1_cache)
331 || (afplus_2 != afplus_2_cache)
333 || (af0_3 != af0_3_cache) || (afplus_3 != afplus_3_cache)
335 || (afplus1 != afplus1_cache) || (afplus2 != afplus2_cache)
336 || (afplus3 != afplus3_cache) || (afplus4 != afplus4_cache)
337 || (afplus5 != afplus5_cache) || (afplus6 != afplus6_cache)
338 || (afplus7 != afplus7_cache) || (afplus8 != afplus8_cache)
339 || (afplus9 != afplus9_cache) || (afplus10 != afplus10_cache)
340 || (afzero1 != afzero1_cache) || (afzero2 != afzero2_cache)
341 || (afzero3 != afzero3_cache) || (afzero4 != afzero4_cache)
342 || (afzero5 != afzero5_cache) || (afzero6 != afzero6_cache)
343 || (afzero7 != afzero7_cache) || (afzero8 != afzero8_cache)
344 || (afzero9 != afzero9_cache) || (afzero10 != afzero10_cache)
345 || (bfplus1 != bfplus1_cache) || (bfplus2 != bfplus2_cache)
346 || (bfplus3 != bfplus3_cache) || (bfplus4 != bfplus4_cache)
347 || (bfplus5 != bfplus5_cache) || (bfplus6 != bfplus6_cache)
348 || (bfplus7 != bfplus7_cache) || (bfplus8 != bfplus8_cache)
349 || (bfplus9 != bfplus9_cache) || (bfplus10 != bfplus10_cache)
350 || (bfzero1 != bfzero1_cache) || (bfzero2 != bfzero2_cache)
351 || (bfzero3 != bfzero3_cache) || (bfzero4 != bfzero4_cache)
352 || (bfzero5 != bfzero5_cache) || (bfzero6 != bfzero6_cache)
353 || (bfzero7 != bfzero7_cache) || (bfzero8 != bfzero8_cache)
354 || (bfzero9 != bfzero9_cache) || (bfzero10 != bfzero10_cache)
355 || (CS != CS_cache) || (CSp != CSp_cache)
356 || (CP != CP_cache) || (CPp != CPp_cache)
357 || (CV != CV_cache) || (CVp != CVp_cache)
358 || (CA != CA_cache) || (CAp != CAp_cache)
359 || (CT != CT_cache) || (CTp != CTp_cache)) {
360 checkcache_int_tau =
false;
361 checkcache_int_mu =
false;
362 checkcache_int_el =
false;
365 if (!checkcache_int_tau || !checkcache_int_mu || !checkcache_int_el) {
367 cached_intJ1_tau = integrateJ(1, q2min, q2max);
368 cached_intJ3_tau = integrateJ(3, q2min, q2max);
369 cached_intJ2_tau = 0.;
371 checkcache_int_tau =
true;
374 cached_intJ1_mu = integrateJ(1, q2min, q2max);
375 cached_intJ3_mu = integrateJ(3, q2min, q2max);
376 cached_intJ2_mu = 0.;
378 checkcache_int_mu =
true;
381 cached_intJ1_el = integrateJ(1, q2min, q2max);
382 cached_intJ3_el = integrateJ(3, q2min, q2max);
383 cached_intJ2_el = 0.;
385 checkcache_int_el =
true;
389 fplusz0_cache = fplusz0;
390 rho1to2_cache = rho1to2;
392 alpha_0_cache = alpha_0;
393 alpha_p_cache = alpha_p;
394 beta_0_cache = beta_0;
395 beta_p_cache = beta_p;
396 gamma_0_cache = gamma_0;
397 gamma_p_cache = gamma_p;
403 afplus_0_cache = afplus_0;
404 afplus_1_cache = afplus_1;
405 afplus_2_cache = afplus_2;
408 afplus_3_cache = afplus_3;
411 z0 = (sqrt(
w0+1.)-sqrt(2.))/(sqrt(
w0+1.)+sqrt(2.));
414 af0_0 -= (af0_1 *
z0 + af0_2 *
z0 *
z0 + af0_3 *
z0 *
z0 *
z0) / phi_f0(
z0) / ((
z0 - z0p_1) / (1. -
z0 * z0p_1)*(
z0 - z0p_2) / (1. -
z0 * z0p_2));
416 af0_0 -= (af0_1 *
z0 + af0_2 *
z0 *
z0) / phi_f0(
z0) / ((
z0 - z0p_1) / (1. -
z0 * z0p_1)*(
z0 - z0p_2) / (1. -
z0 * z0p_2));
418 af0_0 *= phi_f0(
z0)*((
z0 - z0p_1) / (1. -
z0 * z0p_1)*(
z0 - z0p_2) / (1. -
z0 * z0p_2));
421 afplus1_cache = afplus1;
422 afplus2_cache = afplus2;
423 afplus3_cache = afplus3;
424 afplus4_cache = afplus4;
425 afplus5_cache = afplus5;
426 afplus6_cache = afplus6;
427 afplus7_cache = afplus7;
428 afplus8_cache = afplus8;
429 afplus9_cache = afplus9;
430 afplus10_cache = afplus10;
431 afzero1_cache = afzero1;
432 afzero2_cache = afzero2;
433 afzero3_cache = afzero3;
434 afzero4_cache = afzero4;
435 afzero5_cache = afzero5;
436 afzero6_cache = afzero6;
437 afzero7_cache = afzero7;
438 afzero8_cache = afzero8;
439 afzero9_cache = afzero9;
440 afzero10_cache = afzero10;
441 bfplus1_cache = bfplus1;
442 bfplus2_cache = bfplus2;
443 bfplus3_cache = bfplus3;
444 bfplus4_cache = bfplus4;
445 bfplus5_cache = bfplus5;
446 bfplus6_cache = bfplus6;
447 bfplus7_cache = bfplus7;
448 bfplus8_cache = bfplus8;
449 bfplus9_cache = bfplus9;
450 bfplus10_cache = bfplus10;
451 bfzero1_cache = bfzero1;
452 bfzero2_cache = bfzero2;
453 bfzero3_cache = bfzero3;
454 bfzero4_cache = bfzero4;
455 bfzero5_cache = bfzero5;
456 bfzero6_cache = bfzero6;
457 bfzero7_cache = bfzero7;
458 bfzero8_cache = bfzero8;
459 bfzero9_cache = bfzero9;
460 bfzero10_cache = bfzero10;
484double MPlnu::lambda_half(
double a,
double b,
double c)
486 return sqrt(a*a+b*b+c*c-2.*(a*b+a*c+b*c));
492double MPlnu::phi_fplus(
double z)
495 double prefac = 8. * (
MP /
MM)*(
MP /
MM) /
MM * sqrt(8. * nI / (3. * M_PI * chitildeT));
496 double num = (1. + z)*(1. + z) * sqrt(1. - z);
497 double den = (1. +
MP /
MM)*(1. - z) + 2. * sqrt(
MP /
MM)*(1. + z);
498 double den5 = den * den * den * den*den;
499 return prefac * num / den5;
502double MPlnu::phi_f0(
double z)
505 double prefac = (
MP /
MM)*(1. - (
MP /
MM)*(
MP /
MM)) * sqrt(8. * nI / (M_PI * chiL));
506 double num = (1. - z * z) * sqrt((1. - z));
507 double den = (1. +
MP /
MM)*(1. - z) + 2. * sqrt(
MP /
MM)*(1. + z);
508 double den4 = den * den * den*den;
509 return prefac * num / den4;
512double MPlnu::fplus(
double q2)
514 double w =
w0 - q2 / (2. *
MM *
MP);
515 double z = (sqrt(w + 1.) - M_SQRT2) / (sqrt(w + 1.) + M_SQRT2);
517 return fplusz0 * N_0 * (1. - alpha_p*8. * rho1to2 * z + beta_p*(51. * rho1to2 - 10.) * z * z - gamma_p*(252. * rho1to2 - 84.) * z * z * z);
520 double P_fplus = (z - z1m_1) / (1. - z * z1m_1)*(z - z1m_2) / (1. - z * z1m_2)*(z - z1m_3) / (1. - z * z1m_3);
522 return (afplus_0 + afplus_1 * z + afplus_2 * z * z + afplus_3 * z * z * z) / phi_fplus(z) / P_fplus;
524 return (afplus_0 + afplus_1 * z + afplus_2 * z * z) / phi_fplus(z) / P_fplus;
528 double fp_formfac = 0.;
529 if (w<1.06) {fp_formfac = afplus1 + bfplus1 * w;}
530 else if (w<1.12) {fp_formfac = afplus2 + bfplus2 * w;}
531 else if (w<1.18) {fp_formfac = afplus3 + bfplus3 * w;}
532 else if (w<1.24) {fp_formfac = afplus4 + bfplus4 * w;}
533 else if (w<1.30) {fp_formfac = afplus5 + bfplus5 * w;}
534 else if (w<1.36) {fp_formfac = afplus6 + bfplus6 * w;}
535 else if (w<1.42) {fp_formfac = afplus7 + bfplus7 * w;}
536 else if (w<1.48) {fp_formfac = afplus8 + bfplus8 * w;}
537 else if (w<1.54) {fp_formfac = afplus9 + bfplus9 * w;}
538 else {fp_formfac = afplus10 + bfplus10 * w;}
544double MPlnu::f0(
double q2)
546 double w =
w0 - q2 / (2. *
MM *
MP);
547 double z = (sqrt(w + 1.) - M_SQRT2) / (sqrt(w + 1.) + M_SQRT2);
549 double prefac0 = 2. * (
MP /
MM) / (1. +
MP /
MM)/ (1. +
MP /
MM)*(1. +
w0);
550 double norm = prefac0 * (1. - alpha_0*0.0068 * (
w0 - 1.) + beta_0*0.0017 * (
w0 - 1.)*(
w0 - 1.) - gamma_0*0.0013 * (
w0 - 1.)*(
w0 - 1.)*(
w0 - 1.));
551 double prefac = fplus(q2)* prefac0/(1. +
w0)*(1. + w);
553 return prefac/norm * (1. - alpha_0*0.0068 * (w - 1.) + beta_0*0.0017 * (w - 1.)*(w - 1.) - gamma_0*0.0013 * (w - 1.)*(w - 1.)*(w - 1.));
556 double P_f0 = (z - z0p_1) / (1. - z * z0p_1)*(z - z0p_2) / (1. - z * z0p_2);
558 return (
af0_0 + af0_1 * z + af0_2 * z * z + af0_3 * z * z * z) / phi_f0(z) / P_f0;
560 return (
af0_0 + af0_1 * z + af0_2 * z * z) / phi_f0(z) / P_f0;
564 double f0_formfac = 0.;
565 if (w<1.06) {f0_formfac = afzero1 + bfzero1 * w;}
566 else if (w<1.12) {f0_formfac = afzero2 + bfzero2 * w;}
567 else if (w<1.18) {f0_formfac = afzero3 + bfzero3 * w;}
568 else if (w<1.24) {f0_formfac = afzero4 + bfzero4 * w;}
569 else if (w<1.30) {f0_formfac = afzero5 + bfzero5 * w;}
570 else if (w<1.36) {f0_formfac = afzero6 + bfzero6 * w;}
571 else if (w<1.42) {f0_formfac = afzero7 + bfzero7 * w;}
572 else if (w<1.48) {f0_formfac = afzero8 + bfzero8 * w;}
573 else if (w<1.54) {f0_formfac = afzero9 + bfzero9 * w;}
574 else {f0_formfac = afzero10 + bfzero10 * w;}
580double MPlnu::fT(
double q2)
590 return lambda_half(
MM*
MM,
MP*
MP, q2) / 2. / sqrt(q2)*((CV + CVp) * fplus(q2) + 2. *
Mb / (
MM +
MP)*(C7 + C7p) * fT(q2));
593gslpp::complex MPlnu::HA(
double q2)
595 return lambda_half(
MM*
MM,
MP*
MP, q2) / 2. / sqrt(q2)*(CA + CAp) * fplus(q2);
598gslpp::complex MPlnu::HP(
double q2)
600 return (
MM *
MM -
MP *
MP) / 2. * ((CP + CPp) / (
Mb -
Mc)+(
Mlep +
Mnu) / q2 * (CA + CAp)) * f0(q2);
603gslpp::complex MPlnu::HS(
double q2)
605 return (
MM *
MM -
MP *
MP) / 2. * ((CS + CSp) / (
Mb -
Mc)+(
Mlep -
Mnu) / q2 * (CV + CVp)) * f0(q2);
608gslpp::complex MPlnu::HT(
double q2)
610 return -gslpp::complex::i() * lambda_half(
MM*
MM,
MP*
MP, q2) / sqrt(2.) / (
MM +
MP)*(CT - CTp) * fT(q2);
613gslpp::complex MPlnu::HTt(
double q2)
615 return -gslpp::complex::i() * lambda_half(
MM*
MM,
MP*
MP, q2) / 2. / (
MM +
MP)*(CT + CTp) * fT(q2);
621gslpp::complex MPlnu::G0(
double q2)
623 double lambda_MM = lambda_half(
MM*
MM,
MP*
MP, q2);
625 double lambda_lep2 = lambda_lep*lambda_lep;
626 double Elep = sqrt(
Mlep *
Mlep + lambda_lep2 / (4. * q2));
627 double Enu = sqrt(
Mnu *
Mnu + lambda_lep2 / (4. * q2));
628 double Gprefactor = lambda_MM * lambda_lep / q2;
630 return Gprefactor * ((4. * (Elep * Enu +
Mlep *
Mnu) + lambda_lep2 / (3. * q2)) *
HV(q2).abs2()+
631 (4. * (Elep * Enu -
Mlep *
Mnu) + lambda_lep2 / (3. * q2)) * HA(q2).abs2()+
632 (4. * (Elep * Enu -
Mlep *
Mnu) + lambda_lep2 / q2) * HS(q2).abs2()+
633 (4. * (Elep * Enu +
Mlep *
Mnu) + lambda_lep2 / q2) * HP(q2).abs2()+
634 (8. * (Elep * Enu -
Mlep *
Mnu) - lambda_lep2 / (12. * q2)) * HT(q2).abs2()+
635 (16. * (Elep * Enu +
Mlep *
Mnu) - lambda_lep2 / (12. * q2)) * HTt(q2).abs2() +
636 8. * M_SQRT2 *(Enu *
Mlep - Elep *
Mnu)*(HA(q2) * HT(q2).conjugate()).imag() +
637 16. * (Enu *
Mlep + Elep *
Mnu)*(
HV(q2) * HTt(q2).conjugate()).imag());
640gslpp::complex MPlnu::G1(
double q2)
642 double lambda_MM = lambda_half(
MM*
MM,
MP*
MP, q2);
644 double Gprefactor = lambda_MM * lambda_lep / q2;
646 return -Gprefactor * 4. * lambda_lep * ((HA(q2)*(
Mlep -
Mnu) * HP(q2).conjugate()
647 +
HV(q2)*(
Mlep +
Mnu) * HS(q2).conjugate()).real() / sqrt(q2)
648 -(sqrt(2.) * HT(q2) * HP(q2).conjugate() + 2. * HTt(q2) * HS(q2).conjugate()).imag());
651gslpp::complex MPlnu::G2(
double q2)
653 double lambda_MM = lambda_half(
MM*
MM,
MP*
MP, q2);
655 double lambda_lep2 = lambda_lep*lambda_lep;
656 double Gprefactor = lambda_MM * lambda_lep / q2;
658 return -Gprefactor * (4. * lambda_lep2 / (3. * q2))*(HA(q2).abs2() +
HV(q2).abs2() - 2. * HT(q2).abs2() - 4. * HTt(q2).abs2());
664double MPlnu::J1(
double q2)
667 return amplsq_factor * (G0(q2) - G2(q2) / 2).real();
670double MPlnu::J2(
double q2)
673 return amplsq_factor * G1(q2).real();
676double MPlnu::J3(
double q2)
679 return amplsq_factor * 3. / 2. * G2(q2).real();
682double MPlnu::dGammadq2(
double q2)
684 if ((q2 < q2min) or (q2 > (
MM -
MP)*(
MM -
MP)))
return 0.;
685 double sqlambdaB = lambda_half(q2,
MM*
MM,
MP*
MP);
686 double prefac = (CV-CA).abs2()*
MM*sqlambdaB/192./M_PI/M_PI/M_PI;
689 double TotAmp2 = coeff_fp*fplus(q2)*fplus(q2)+coeff_f0*f0(q2)*f0(q2);
697double MPlnu::integrateJ(
int i,
double q2_min,
double q2_max)
702 if (
lep ==
StandardModel::TAU)
if (checkcache_int_tau && (q2_min == q2min) && (q2_max == q2max))
return cached_intJ1_tau;
703 if (
lep ==
StandardModel::MU)
if (checkcache_int_mu && (q2_min == q2min) && (q2_max == q2max))
return cached_intJ1_mu;
705 wf=ROOT::Math::Functor1D(&(*
this),&MPlnu::J1);
707 J_res = ig.Integral(q2_min, q2_max);
708 if (ig.Status() != 0)
return std::numeric_limits<double>::quiet_NaN();
712 if (
lep ==
StandardModel::TAU)
if (checkcache_int_tau && (q2_min == q2min) && (q2_max == q2max))
return cached_intJ2_tau;
713 if (
lep ==
StandardModel::MU)
if (checkcache_int_mu && (q2_min == q2min) && (q2_max == q2max))
return cached_intJ2_mu;
715 wf=ROOT::Math::Functor1D(&(*
this),&MPlnu::J2);
717 J_res = ig.Integral(q2_min, q2_max);
718 if (ig.Status() != 0)
return std::numeric_limits<double>::quiet_NaN();
722 if (
lep ==
StandardModel::TAU)
if (checkcache_int_tau && (q2_min == q2min) && (q2_max == q2max))
return cached_intJ3_tau;
723 if (
lep ==
StandardModel::MU)
if (checkcache_int_mu && (q2_min == q2min) && (q2_max == q2max))
return cached_intJ3_mu;
725 wf=ROOT::Math::Functor1D(&(*
this),&MPlnu::J3);
727 J_res = ig.Integral(q2_min, q2_max);
728 if (ig.Status() != 0)
return std::numeric_limits<double>::quiet_NaN();
731 wf=ROOT::Math::Functor1D(&(*
this),&MPlnu::dGammadq2);
733 J_res = ig.Integral(q2_min, q2_max);
734 if (ig.Status() != 0)
return std::numeric_limits<double>::quiet_NaN();
738 std::stringstream out;
740 throw std::runtime_error(
"MPlnu::integrateJ: index " + out.str() +
" not implemented");
745double MPlnu::dGammadw(
double q2)
749 return 2. * (J1(q2) + J3(q2) / 3.);
757 double q2_min = (2. *
MM *
MP)*(
w0 - w_max);
758 double q2_max = (2. *
MM *
MP)*(
w0 - w_min);
760 double intJ1 = integrateJ(1, q2_min, q2_max);
761 double intJ3 = integrateJ(3, q2_min, q2_max);
763 return 2. * (intJ1 + intJ3 / 3.);
774 return afplus_0 * afplus_0 + afplus_1 * afplus_1 + afplus_2 * afplus_2 + afplus_3 * afplus_3;
776 return afplus_0 * afplus_0 + afplus_1 * afplus_1 + afplus_2 * afplus_2;
785 return af0_0 *
af0_0 + af0_1 * af0_1 + af0_2 * af0_2 + af0_3 * af0_3;
787 return af0_0 *
af0_0 + af0_1 * af0_1 + af0_2*af0_2;
796 return 1707.54 * afplus_0 * afplus_0 + 1299.57 * afplus_0 * afplus_1 + 442.82 * afplus_1 * afplus_1 - 356.01 * afplus_0 * afplus_2
797 - 101.62 * afplus_1 * afplus_2 + 34.947 * afplus_2 * afplus_2 - 206.767 * afplus_0 * afplus_3 - 127.668 * afplus_1 * afplus_3
798 + 33.234 * afplus_2 * afplus_3 + 16.475 * afplus_3 * afplus_3;
800 return 442.82 * afplus_0 * afplus_0 - 101.619 * afplus_0 * afplus_1 + 34.947* afplus_1 * afplus_1 - 127.668 * afplus_0 * afplus_2 + 33.234 * afplus_1 * afplus_2 + 16.4754 * afplus_2 * afplus_2;
void setUpdateFlag(QCD::meson meson_i, QCD::meson meson_j, QCD::lepton lep_i, bool updated_i) const
sets the update flag for the initial and final state dependent object for .
bool getUpdateFlag(QCD::meson meson_i, QCD::meson meson_j, QCD::lepton lep_i) const
gets the update flag for the initial and final state dependent object for .
gslpp::vector< gslpp::complex > ** ComputeCoeffdiujlknu(int i, int j, int k, double mu) const
Computes the Wilson coefficient for the Hamiltonian transitions in the JMS basis ordered as CnueduVL...
virtual ~MPlnu()
Destructor.
double get_fT(double q2)
return fT form factor at
std::vector< std::string > mplnuParameters
double getDeltaGammaDeltaw(double w_min, double w_max)
The integral of from to .
std::vector< std::string > initializeMPlnuParameters()
double get_fplus(double q2)
return fplus form factor at
MPlnu(const StandardModel &SM_i, QCD::meson meson_i, QCD::meson pseudoscalar_i, QCD::lepton lep_i)
Constructor.
double get_unitarity_1min_BGL()
Weak Unitarity constraint for BGL parameters related to 1- resonances.
double get_strong_unitarity_BGL()
Strong Unitarity constraint for BGL parameters using HQET.
double get_unitarity_0plus_BGL()
Weak Unitarity constraint for BGL parameters related to 0+ resonances.
const StandardModel & mySM
double get_f0(double q2)
return f0 form factor at
double computeWidth() const
A method to compute the width of the meson from its lifetime.
std::string getModelName() const
A method to fetch the name of the model.
const double & getMass() const
A get method to access the particle mass.
meson
An enum type for mesons.
const double getOptionalParameter(std::string name) const
A method to get parameters that are specific to only one set of observables.
const Meson & getMesons(const QCD::meson m) const
A get method to access a meson as an object of the type Meson.
lepton
An enum type for leptons.
const Particle & getQuarks(const QCD::quark q) const
A get method to access a quark as an object of the type Particle.
void initializeMeson(QCD::meson meson_i) const
A method to initialize a meson.
A model class for the Standard Model.
virtual const double getCCC1() const
A virtual implementation for the RealWeakEFTCC class.
const Particle & getLeptons(const QCD::lepton p) const
A get method to retrieve the member object of a lepton.
virtual const double getCCC3() const
A virtual implementation for the RealWeakEFTCC class.
const double getMz() const
A get method to access the mass of the boson .
const double Ale(double mu, orders order, bool Nf_thr=true) const
The running electromagnetic coupling in the scheme.
virtual const double getCCC4() const
A virtual implementation for the RealWeakEFTCC class.
const Flavour & getFlavour() const
virtual const double getCCC2() const
A virtual implementation for the RealWeakEFTCC class.
virtual const double getCCC5() const
A virtual implementation for the RealWeakEFTCC class.