a Code for the Combination of Indirect and Direct Constraints on High Energy Physics Models Logo
StandardModel.h
Go to the documentation of this file.
1/*
2 * Copyright (C) 2012 HEPfit Collaboration
3 *
4 *
5 * For the licensing terms see doc/COPYING.
6 */
7
8#ifndef STANDARDMODEL_H
9#define STANDARDMODEL_H
10
11#include "Flavour.h"
12#include "QCD.h"
13#include "CKM.h"
14#include "PMNS.h"
16#include "Matching.h"
17#include <gsl/gsl_integration.h>
18
19class EWSMcache;
20class EWSMOneLoopEW;
21class EWSMTwoLoopQCD;
22class EWSMTwoLoopEW;
25class EWSMThreeLoopEW;
27class LeptonFlavour;
28/* BEGIN: REMOVE FROM THE PACKAGE */
30/* END: REMOVE FROM THE PACKAGE */
31
32
521class StandardModel : public QCD {
522public:
523
524
525 // Radiative Corrections for the LEP-II observables
526 enum LEP2RCs {
527 Weak = 0,
533 };
534
540 EW1 = 0,
547 };
548
549
550
551// static const int NSMvars = 41; ///< The number of the model parameters in %StandardModel. !!! PMNS INCLUDED
552 static const int NSMvars = 28;
556 static std::string SMvars[NSMvars];
557
558 static const double GeVminus2_to_nb;
559
564 static const double Mw_error;
565
570
574 virtual ~StandardModel();
575
576
578 // Initialization
579
586 virtual bool InitializeModel();
587
588
590 // Model parameters
591
598 virtual bool Init(const std::map<std::string, double>& DPars);
599
607 virtual bool PreUpdate();
608
616 virtual bool Update(const std::map<std::string, double>& DPars);
617
627 virtual bool PostUpdate();
628
629 const int getIterationNo() const
630 {
631 return iterationNo;
632 }
633
641 virtual bool CheckParameters(const std::map<std::string, double>& DPars);
642
643
645 // Flags
646
653 virtual bool setFlag(const std::string name, const bool value);
654
661 virtual bool setFlagStr(const std::string name, const std::string value);
662
667 virtual bool CheckFlags() const;
668
679 {
681 }
682
691 const bool IsFlagNoApproximateGammaZ() const
692 {
694 }
695
697 {
698 this->FlagNoApproximateGammaZ = FlagNoApproximateGammaZ;
699 }
700
706 const std::string getFlagMw() const
707 {
708 return FlagMw;
709 }
710
716 const std::string getFlagRhoZ() const
717 {
718 return FlagRhoZ;
719 }
720
726 const std::string getFlagKappaZ() const
727 {
728 return FlagKappaZ;
729 }
730
743 {
744 this->FlagCacheInStandardModel = FlagCacheInStandardModel;
745 }
746
747
749 // get and set methods for class members
750
756 const Particle& getLeptons(const QCD::lepton p) const
757 {
758 return leptons[p];
759 }
760
765 const double getMz() const
766 {
767 return Mz;
768 }
769
774 const double getMw() const
775 {
776 return Mw_inp;
777 }
778
783 const double getAlsMz() const
784 {
785 return AlsMz;
786 }
787
792 const double getGF() const
793 {
794 return GF;
795 }
796
801 const double getAle() const
802 {
803 return ale;
804 }
805
812 const double getDAle5Mz() const
813 {
814 return dAle5Mz;
815 }
816
821 virtual const double getMHl() const
822 {
823 return mHl;
824 }
825
831 const double getDelMw() const
832 {
833 return delMw;
834 }
835
842 const double getDelSin2th_l() const
843 {
844 return delSin2th_l;
845 }
846
853 const double getDelSin2th_q() const
854 {
855 return delSin2th_q;
856 }
857
864 const double getDelSin2th_b() const
865 {
866 return delSin2th_b;
867 }
868
874 const double getDelGammaZ() const
875 {
876 return delGammaZ;
877 }
878
884 const double getDelSigma0H() const
885 {
886 return delsigma0H;
887 }
888
894 const double getDelR0l() const
895 {
896 return delR0l;
897 }
898
904 const double getDelR0c() const
905 {
906 return delR0c;
907 }
908
914 const double getDelR0b() const
915 {
916 return delR0b;
917 }
918
924 const double getDelGammaWlv() const
925 {
926 return delGammaWlv;
927 }
928
934 const double getDelGammaWqq() const
935 {
936 return delGammaWqq;
937 }
938
943 const gslpp::matrix<gslpp::complex> getVCKM() const
944 {
945 return myCKM.getCKM();
946 }
947
952 const CKM& getCKM() const
953 {
954 return myCKM;
955 }
956
957
958
963 const gslpp::matrix<gslpp::complex> getUPMNS() const
964 {
965 return myPMNS.getPMNS();
966 }
967
973 const gslpp::matrix<gslpp::complex>& getYn() const
974 {
975 return Yn;
976 }
977
983 const double getMuw() const
984 {
985 return muw;
986 }
987
988 virtual const StandardModel& getTrueSM() const
989 {
990 throw std::runtime_error("StandardModel::getTrueSM() must be overridden by the NP extension.");
991 }
992
998 {
999 return SMM.getObj();
1000 }
1001
1007 {
1008 return myEWSMcache;
1009 }
1010
1016 {
1017 return myOneLoopEW;
1018 }
1019
1025 {
1026 return myApproximateFormulae;
1027 }
1028
1029 /* BEGIN: REMOVE FROM THE PACKAGE */
1035 {
1036 return myTwoFermionsLEP2;
1037 }
1038 /* END: REMOVE FROM THE PACKAGE */
1039
1041 {
1042 return myThreeLoopEW;
1043 }
1044
1046 {
1047 return myThreeLoopEW2QCD;
1048 }
1049
1051 {
1052 return myThreeLoopQCD;
1053 }
1054
1056 {
1057 return myTwoLoopEW;
1058 }
1059
1061 {
1062 return myTwoLoopQCD;
1063 }
1064
1065 const Flavour& getFlavour() const
1066 {
1067 return SMFlavour;
1068 }
1069
1071 {
1072 return myLeptonFlavour;
1073 }
1075 // QED coupling
1076
1088 const double ale_OS(const double mu, orders order = FULLNLO) const;
1089
1096 const double Beta_s(int nm, unsigned int nf) const;
1097
1104 const double Beta_e(int nm, unsigned int nf) const;
1105
1106 //make the QCD implementation of Als available here
1107 using QCD::Als;
1108
1120 const double Als(const double mu, const orders order, const bool Nf_thr, const bool qed_flag) const
1121 {
1122 if (qed_flag && order == FULLNNNLO)
1123 return AlsE(mu, order, Nf_thr);
1124
1125 return Als(mu, order, Nf_thr);
1126 }
1127
1137 const double Ale(double mu, orders order, bool Nf_thr = true) const;
1138
1145 const double DeltaAlphaLepton(const double s) const;
1146
1159 const double DeltaAlphaL5q() const;
1160
1167 const double DeltaAlphaTop(const double s) const;
1168
1181 const double DeltaAlpha() const;
1182
1192 virtual const double alphaMz() const;
1193
1200 const double Alstilde5(const double mu) const;
1201
1205 virtual const double getCCC1() const { return 0.; };
1206
1210 virtual const double getCCC2() const { return 0.; };
1211
1215 virtual const double getCCC3() const { return 0.; };
1216
1220 virtual const double getCCC4() const { return 0.; };
1221
1225 virtual const double getCCC5() const { return 0.; };
1226
1227
1228
1230 // Higgs VEV
1231
1240 const double v() const;
1241
1242
1244 // The W-boson mass
1245
1250 const double Mw_tree() const;
1251
1267 const double s02() const;
1268
1281 const double c02() const;
1282
1309 virtual const double Mw() const;
1310
1311
1320 virtual const double Dalpha5hMz() const;
1321
1331 virtual const double cW2(const double Mw_i) const;
1332 virtual const double cW2() const;
1333
1343 virtual const double sW2(const double Mw_i) const;
1344 const double sW2() const;
1345
1352 const double sW2_MSbar_Approx() const;
1353
1360 const double sW2_ND() const;
1361
1362
1363
1385 virtual const double DeltaR() const;
1386
1395 void ComputeDeltaRho(const double Mw_i, double DeltaRho[orders_EW_size]) const;
1396
1405 void ComputeDeltaR_rem(const double Mw_i, double DeltaR_rem[orders_EW_size]) const;
1406
1407
1409 // The W and Z masses in the complex-pole/fixed-width scheme
1410
1441 double Mzbar() const;
1442
1463 const double MwbarFromMw(const double Mw) const;
1464
1488 const double MwFromMwbar(const double Mwbar) const;
1489
1507 virtual const double DeltaRbar() const;
1508
1509
1511 // The W-boson decay width
1512
1521 virtual const double rho_GammaW(const Particle fi, const Particle fj) const;
1522
1549 virtual const double GammaW(const Particle fi, const Particle fj) const;
1550
1555 virtual const double GammaW() const;
1556
1561 virtual const double BrW(const Particle fi, const Particle fj) const;
1562
1567 virtual const double RWlilj(const Particle li, const Particle lj) const;
1568
1569 virtual const double Ruc() const; //AG:added
1570
1575 virtual const double RWc() const;
1576
1578 // EWPO at Z-pole
1579
1598 virtual const double A_f(const Particle f) const;
1599
1605 virtual const double AFB(const Particle f) const;
1606
1626 virtual const double sin2thetaEff(const Particle f) const;
1627
1653 virtual const double GammaZ(const Particle f) const;
1654
1665 virtual const double Gamma_inv() const;
1666
1681 virtual const double Gamma_had() const;
1682
1696 virtual const double Gamma_Z() const;
1697
1702 virtual const double RZlilj(const Particle li, const Particle lj) const;
1703
1717 virtual const double sigma0_had() const;
1718
1733 virtual const double R0_f(const Particle f) const;
1734
1743 virtual const double R_inv() const;
1744
1754 virtual const double N_nu() const;
1755
1756
1758 // Zff effective couplings
1759
1769 virtual const gslpp::complex gV_f(const Particle f) const;
1770
1780 virtual const gslpp::complex gA_f(const Particle f) const;
1781
1797 virtual const gslpp::complex rhoZ_f(const Particle f) const;
1798
1826 virtual const gslpp::complex kappaZ_f(const Particle f) const;
1827
1848 virtual const gslpp::complex deltaRhoZ_f(const Particle f) const;
1849
1874 virtual const gslpp::complex deltaKappaZ_f(const Particle f) const;
1875
1876
1878 // Epsilon parameters for EWPO
1879
1891 virtual const double epsilon1() const;
1892
1913 virtual const double epsilon2() const;
1914
1932 virtual const double epsilon3() const;
1933
1951 virtual const double epsilonb() const;
1952
1953
1955 // EW low-energy observables: Parity violation
1956
1957
1963 virtual const double amuon() const;
1964
1972 virtual const double Qwemoller(const double q2, const double y) const;
1973
1974
1982 virtual const double alrmoller(const double q2, const double y) const;
1983
1984
1991 virtual const double Qwp() const;
1992
1993
2000 virtual const double Qwn() const;
2001
2003 // EW low-energy observables: neutrino-scattering
2004
2011 virtual const double gLnuN2() const;
2012
2013
2020 virtual const double gRnuN2() const;
2021
2022
2029 virtual const double ThetaLnuN() const;
2030
2031
2038 virtual const double ThetaRnuN() const;
2039
2040
2047 virtual const double gVnue() const;
2048
2049
2056 virtual const double gAnue() const;
2057
2058
2060 // Lepton decays
2061
2062 // Muon decays
2063
2069 virtual const double Gamma_muon() const;
2070
2071
2072 // Tau decays
2073
2079 virtual const double Gamma_tau_l_nunu(const Particle l) const;
2080
2081
2082 // Lepton Flavor universality tests in Tau decays
2083
2089 virtual const double TauLFU_gmuge() const;
2090
2091
2097 virtual const double TauLFU_gtaugmu() const;
2098
2099
2105 virtual const double TauLFU_gtauge() const;
2106
2112 virtual const double TauLFU_gtaugmuPi() const;
2113
2119 virtual const double TauLFU_gtaugmuK() const;
2120
2121
2123 // For EWPO caches
2124
2132 static const int NumSMParamsForEWPO = 35;
2133
2148 bool checkSMparamsForEWPO();
2149
2151 // Several Higgs-related quantities used in Higgs coupling analysis
2153
2159 gslpp::complex f_triangle(const double tau) const;
2165 gslpp::complex g_triangle(const double tau) const;
2171 gslpp::complex I_triangle_1(const double tau, const double lambda) const;
2177 gslpp::complex I_triangle_2(const double tau, const double lambda) const;
2185 gslpp::complex AH_f(const double tau) const;
2186
2194 gslpp::complex AH_W(const double tau) const;
2195
2202 gslpp::complex AHZga_f(const double tau, const double lambda) const;
2203
2210 gslpp::complex AHZga_W(const double tau, const double lambda) const;
2211
2213
2219 virtual const double SigmaeeZH(const double sqrt_s, const double Pe, const double Pp) const;
2220
2226 virtual const double SigmaeeHvv(const double sqrt_s, const double Pe, const double Pp) const;
2227
2233 virtual const double SigmaeeHee(const double sqrt_s, const double Pe, const double Pp) const;
2234
2236
2242 virtual const double GammaHtogg() const;
2243
2249 virtual const double GammaHtoZZstar() const;
2250
2256 virtual const double GammaHtoWWstar() const;
2257
2263 virtual const double GammaHtoZga() const;
2264
2270 virtual const double GammaHtogaga() const;
2271
2277 virtual const double GammaHtomumu() const;
2278
2284 virtual const double GammaHtotautau() const;
2285
2291 virtual const double GammaHtocc() const;
2292
2298 virtual const double GammaHtoss() const;
2299
2305 virtual const double GammaHtobb() const;
2306
2312 virtual const double GammaHTot() const;
2313
2315
2321 virtual const double BrHtogg() const;
2322
2328 virtual const double BrHtoZZstar() const;
2329
2335 virtual const double BrHtoWWstar() const;
2336
2342 virtual const double BrHtoZga() const;
2343
2349 virtual const double BrHtogaga() const;
2350
2356 virtual const double BrHtomumu() const;
2357
2363 virtual const double BrHtotautau() const;
2364
2370 virtual const double BrHtocc() const;
2371
2377 virtual const double BrHtoss() const;
2378
2384 virtual const double BrHtobb() const;
2385
2387
2403 const double computeSigmaggH(const double sqrt_s) const
2404 {
2405 if (sqrt_s == 7.0) {
2406 return 16.83; // in pb for Mh=125.1 GeV
2407 } else if (sqrt_s == 8.0) {
2408 return 21.40; // in pb for Mh=125.1 GeV
2409 } else if (sqrt_s == 13.0) {
2410 return 48.61; // in pb for Mh=125.09 GeV
2411 } else if (sqrt_s == 13.6) {
2412 return 52.17; // in pb for Mh=125.09 GeV
2413 } else if (sqrt_s == 14.0) {
2414 return 54.72; // in pb for Mh=125.09 GeV
2415 } else if (sqrt_s == 27.0) {
2416 return 146.65; // in pb for Mh=125.09 GeV
2417 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2418 return 740.3; // in pb for Mh=125. GeV
2419 } else if (sqrt_s == 1.96) {
2420 return 0.9493; // in pb for Mh=125 GeV
2421 } else
2422 throw std::runtime_error("Bad argument in StandardModel::computeSigmaggH()");
2423 }
2424
2435 const double computeSigmaggH_tt(const double sqrt_s) const
2436 {
2437 if (sqrt_s == 7.0) {
2438 return 16.69; // in pb for Mh=125.09 GeV
2439 } else if (sqrt_s == 8.0) {
2440 return 21.20; // in pb for Mh=125.09 GeV
2441 } else if (sqrt_s == 13.0) {
2442 return 47.94; // in pb for Mh=125.09 GeV
2443 } else if (sqrt_s == 13.6) {
2444 return 51.534; // in pb for Mh=125.09 GeV (interpolation between 13 and 14 TeV)
2445 } else if (sqrt_s == 14.0) {
2446 return 53.93; // in pb for Mh=125.09 GeV
2447 } else if (sqrt_s == 27.0) {
2448 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_tt(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2449 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2450 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_tt(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2451 } else
2452 throw std::runtime_error("Bad argument in StandardModel::computeSigmaggH_tt()");
2453 }
2454
2465 const double computeSigmaggH_bb(const double sqrt_s) const
2466 {
2467 if (sqrt_s == 7.0) {
2468 return 0.04; // in pb for Mh=125.09 GeV
2469 } else if (sqrt_s == 8.0) {
2470 return 0.05; // in pb for Mh=125.09 GeV
2471 } else if (sqrt_s == 13.0) {
2472 return 0.10; // in pb for Mh=125.09 GeV
2473 } else if (sqrt_s == 13.6) {
2474 return 0.106; // in pb for Mh=125.09 GeV (interpolation between 13 and 14 TeV)
2475 } else if (sqrt_s == 14.0) {
2476 return 0.11; // in pb for Mh=125.09 GeV
2477 } else if (sqrt_s == 27.0) {
2478 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_bb(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2479 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2480 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_bb(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2481 } else
2482 throw std::runtime_error("Bad argument in StandardModel::computeSigmaggH_bb()");
2483 }
2484
2495 const double computeSigmaggH_tb(const double sqrt_s) const
2496 {
2497 if (sqrt_s == 7.0) {
2498 return -0.66; // in pb for Mh=125.09 GeV
2499 } else if (sqrt_s == 8.0) {
2500 return -0.82; // in pb for Mh=125.09 GeV
2501 } else if (sqrt_s == 13.0) {
2502 return -1.73; // in pb for Mh=125.09 GeV
2503 } else if (sqrt_s == 13.6) {
2504 return -1.844; // in pb for Mh=125.09 GeV (interpolation between 13 and 14 TeV)
2505 } else if (sqrt_s == 14.0) {
2506 return -1.92; // in pb for Mh=125.09 GeV
2507 } else if (sqrt_s == 27.0) {
2508 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_tb(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2509 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2510 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_tb(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2511 } else
2512 throw std::runtime_error("Bad argument in StandardModel::computeSigmaggH_tb()");
2513 }
2514
2528 const double computeSigmaVBF(const double sqrt_s) const
2529 {
2530 if (sqrt_s == 7.0) {
2531 return 1.241; // in pb for Mh=125.09 GeV
2532 } else if (sqrt_s == 8.0) {
2533 return 1.601; // in pb for Mh=125.09 GeV
2534 } else if (sqrt_s == 13.0) {
2535 return 3.766; // in pb for Mh=125.09 GeV
2536 } else if (sqrt_s == 13.6) {
2537 return 4.075; // in pb for Mh=125.09 GeV
2538 } else if (sqrt_s == 14.0) {
2539 return 4.260; // in pb for Mh=125.09 GeV
2540 } else if (sqrt_s == 27.0) {
2541 return 11.838; // in pb for Mh=125.09 GeV
2542 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2543 return 82.0; // in pb for Mh=125. GeV
2544 } else if (sqrt_s == 1.96) {
2545 return 0.0653; // in pb for Mh=125 GeV
2546 } else
2547 throw std::runtime_error("Bad argument in StandardModel::computeSigmaVBF()");
2548 }
2549
2561 const double computeSigmaWF(const double sqrt_s) const
2562 {
2563 if (sqrt_s == 7.0) {
2564 return 0.946; // in pb for Mh=125 GeV
2565 } else if (sqrt_s == 8.0) {
2566 return 1.220; // in pb for Mh=125 GeV
2567 } else if (sqrt_s == 13.0) {
2568 return 2.882; // in pb for Mh=125 GeV
2569 } else if (sqrt_s == 13.6) {
2570 return 3.1088; // in pb for Mh=125 GeV (interpolation between 13 and 14 TeV)
2571 } else if (sqrt_s == 14.0) {
2572 return 3.260; // in pb for Mh=125 GeV
2573 } else if (sqrt_s == 27.0) {
2574 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(14.) * computeSigmaWF(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2575 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2576 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(14.) * computeSigmaWF(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2577 } else if (sqrt_s == 1.96) {
2578 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(7.) * computeSigmaWF(7.); // in the absence of individual cross sections for TeVatron we rescale the LHC ones
2579 } else
2580 throw std::runtime_error("Bad argument in StandardModel::computeSigmaWF()");
2581 }
2582
2594 const double computeSigmaZF(const double sqrt_s) const
2595 {
2596 if (sqrt_s == 7.0) {
2597 return 0.333; // in pb for Mh=125 GeV
2598 } else if (sqrt_s == 8.0) {
2599 return 0.432; // in pb for Mh=125 GeV
2600 } else if (sqrt_s == 13.0) {
2601 return 1.049; // in pb for Mh=125 GeV
2602 } else if (sqrt_s == 13.6) {
2603 return 1.1342; // in pb for Mh=125 GeV (interpolation between 13 and 14 TeV)
2604 } else if (sqrt_s == 14.0) {
2605 return 1.191; // in pb for Mh=125 GeV
2606 } else if (sqrt_s == 27.0) {
2607 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(14.) * computeSigmaZF(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2608 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2609 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(14.) * computeSigmaZF(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2610 } else if (sqrt_s == 1.96) {
2611 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(7.) * computeSigmaZF(7.); // in the absence of individual cross sections for TeVatron we rescale the LHC ones
2612 } else
2613 throw std::runtime_error("Bad argument in StandardModel::computeSigmaZF()");
2614 }
2615
2623 const double computeSigmaZWF(const double sqrt_s) const
2624 {
2625 return 0.;
2626 }
2627
2641 const double computeSigmaWH(const double sqrt_s) const
2642 {
2643 if (sqrt_s == 7.0) {
2644 return 0.577; // in pb for Mh=125.1 GeV
2645 //return 0.5688; // in pb for Mh=125.6 GeV
2646 } else if (sqrt_s == 8.0) {
2647 return 0.7027; // in pb for Mh=125.1 GeV
2648 //return 0.6931; // in pb for Mh=125.6 GeV
2649 } else if (sqrt_s == 13.0) {
2650 return 1.358; // in pb for Mh=125.09 GeV
2651 } else if (sqrt_s == 13.6) {
2652 return 1.453; // in pb for Mh=125.09 GeV
2653 } else if (sqrt_s == 14.0) {
2654 return 1.498; // in pb for Mh=125.09 GeV
2655 } else if (sqrt_s == 27.0) {
2656 return 3.397; // in pb for Mh=125.09 GeV
2657 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2658 return 15.9; // in pb for Mh=125. GeV
2659 } else if (sqrt_s == 1.96) {
2660 return 0.1295; // in pb for Mh=125 GeV
2661 } else
2662 throw std::runtime_error("Bad argument in StandardModel::computeSigmaWH()");
2663 }
2664
2678 const double computeSigmaZH(const double sqrt_s) const
2679 {
2680 if (sqrt_s == 7.0) {
2681 return 0.3341; // in pb for Mh=125.1 GeV
2682 //return 0.3299; // in pb for Mh=125.6 GeV
2683 } else if (sqrt_s == 8.0) {
2684 return 0.4142; // in pb for Mh=125.1 GeV
2685 //return 0.4091; // in pb for Mh=125.6 GeV
2686 } else if (sqrt_s == 13.0) {
2687 return 0.880; // in pb for Mh=125.09 GeV
2688 } else if (sqrt_s == 13.6) {
2689 return 0.9422; // in pb for Mh=125.09 GeV
2690 } else if (sqrt_s == 14.0) {
2691 return 0.981; // in pb for Mh=125.09 GeV
2692 } else if (sqrt_s == 27.0) {
2693 return 2.463; // in pb for Mh=125.09 GeV
2694 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2695 return 11.26; // in pb for Mh=125. GeV
2696 } else if (sqrt_s == 1.96) {
2697 return 0.0785; // in pb for Mh=125 GeV
2698 } else
2699 throw std::runtime_error("Bad argument in StandardModel::computeSigmaZH()");
2700 }
2701
2718 const double computeSigmattH(const double sqrt_s) const
2719 {
2720 if (sqrt_s == 7.0) {
2721 return 0.0861; // in pb for Mh=125.1 GeV
2722 //return 0.0851; // in pb for Mh=125.6 GeV
2723 } else if (sqrt_s == 8.0) {
2724 return 0.129; // in pb for Mh=125.1 GeV
2725 //return 0.1274; // in pb for Mh=125.6 GeV
2726 } else if (sqrt_s == 13.0) {
2727 return 0.5060; // in pb for Mh=125.1 GeV
2728 } else if (sqrt_s == 13.6) {
2729 return 0.5688; // in pb for Mh=125.09 GeV
2730 } else if (sqrt_s == 14.0) {
2731 return 0.6128; // in pb for Mh=125.09 GeV
2732 } else if (sqrt_s == 27.0) {
2733 return 2.86; // in pb for Mh=125.09 GeV
2734 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2735 return 37.9; // in pb for Mh=125. GeV
2736 } else if (sqrt_s == 1.96) {
2737 return 0.0043; // in pb for Mh=125 GeV
2738 } else
2739 throw std::runtime_error("Bad argument in StandardModel::computeSigmattH()");
2740 }
2741
2748 const double computeSigmatHq(const double sqrt_s) const
2749 {
2750 if (sqrt_s == 13.0) {
2751 return 0.07426; // in pb for Mh=125.09 GeV
2752 } else
2753 throw std::runtime_error("Bad argument in StandardModel::computeSigmatHq()");
2754 }
2755
2756
2757
2767 const double computeSigmabbH(const double sqrt_s) const
2768 {
2769 if (sqrt_s == 13.0){
2770 return 0.4863; // in pb for Mh=125.09 GeV
2771 }
2772 else if (sqrt_s == 13.6) {
2773 return 0.566; // in pb for Mh=125.09 GeV (NLO+NNLLpart+ybyt matching)
2774 } else
2775 throw std::runtime_error("Bad argument in StandardModel::computeSigmabbH()");
2776 }
2777
2784 const double computeBrHtogg() const
2785 {
2786 return 8.179e-2; // Mh=125.1 GeV
2787 }
2788
2795 const double computeBrHtoWW() const
2796 {
2797 //return 2.23e-1; // Mh=125.5 GeV
2798 return 2.154e-1; // Mh=125.1 GeV
2799 }
2800
2807 const double computeBrHtoZZ() const
2808 {
2809 return 2.643e-2; // Mh=125.1 GeV
2810 //return 2.79e-2; // Mh=125.6 GeV
2811 }
2812
2819 const double computeBrHtoZga() const
2820 {
2821 return 1.541e-3; // Mh=125.1 GeV
2822 //return 1.59e-3; // Mh=125.6 GeV
2823 }
2824
2831 const double computeBrHtogaga() const
2832 {
2833 return 2.27e-3; // Mh=125.1 GeV
2834 }
2835
2842 const double computeBrHtomumu() const
2843 {
2844 return 2.17e-4; // Mh=125.1 GeV
2845 }
2846
2853 const double computeBrHtotautau() const
2854 {
2855 return 6.256e-2; // Mh=125.1 GeV
2856 //return 6.22e-2; // Mh=125.6 GeV
2857 }
2858
2865 const double computeBrHtocc() const
2866 {
2867 return 2.883e-2; // Mh=125.1 GeV
2868 //return 2.86e-2; // Mh=125.6 GeV
2869 }
2870
2877 const double computeBrHtoss() const
2878 {
2879 return 4.0e-4;
2880 }
2881
2888 const double computeBrHtobb() const
2889 {
2890 return 5.807e-1; // Mh=125.1 GeV
2891 //return 5.67e-1; // Mh=125.6 GeV
2892 }
2893
2894 // H decays into 4-fermions
2895
2901 const double computeBrHto4v() const
2902 {
2903 return 1.06e-3;
2904 }
2905
2912 const double computeBrHtoevmuv() const
2913 {
2914 return 2.539e-03; // Mh=125.1 GeV
2915 }
2916
2923 const double computeBrHtollvv2() const
2924 {
2925 return 1.063e-02; // Mh=125.1 GeV
2926 }
2927
2934 const double computeBrHtollvv3() const
2935 {
2936 return 2.356e-02; // Mh=125.1 GeV
2937 }
2938
2945 const double computeBrHto4l2() const
2946 {
2947 return 1.252e-04; // Mh=125.1 GeV
2948 }
2949
2956 const double computeBrHto4l3() const
2957 {
2958 return 2.771e-04; // Mh=125.1 GeV
2959 }
2960
2967 const double computeBrHto4q() const
2968 {
2969 return 1.098e-01; // Mh=125.1 GeV
2970 }
2971
2978 const double computeBrHto4f() const
2979 {
2980 return 2.406e-01; // Mh=125.1 GeV
2981 }
2982
2989 const double computeGammaHTotal() const
2990 {
2991 return 4.101e-3; // Mh=125.1 GeV
2992 //return 4.15e-3; // Mh=125.6 GeV
2993 }
2994
3000 const double computeGammaHgg_tt() const
3001 {
3002 return 380.8; // in keV for Mh=125 GeV
3003 //return 389.6; // in keV for Mh=126 GeV
3004 }
3005
3011 const double computeGammaHgg_bb() const
3012 {
3013 return 3.96; // in keV for Mh=125 GeV
3014 //return 3.95; // in keV for Mh=126 GeV
3015 }
3016
3022 const double computeGammaHgg_tb() const
3023 {
3024 return -42.1; // in keV for Mh=125 GeV
3025 //return -42.7; // in keV for Mh=126 GeV
3026 }
3027
3033 const double computeGammaHZga_tt() const
3034 {
3035 return 21.74; // in eV for Mh=125 GeV
3036 //return 23.51; // in eV for Mh=126 GeV
3037 }
3038
3044 const double computeGammaHZga_WW() const
3045 {
3046 return 7005.6; // in eV for Mh=125 GeV
3047 //return 7648.4; // in eV for Mh=126 GeV
3048 }
3049
3055 const double computeGammaHZga_tW() const
3056 {
3057 return -780.4; // in eV for Mh=125 GeV
3058 //return -848.1; // in eV for Mh=126 GeV
3059 }
3060
3066 const double computeGammaHgaga_tt() const
3067 {
3068 return 662.84; // in eV for Mh=125 GeV
3069 //return 680.39; // in eV for Mh=126 GeV
3070 }
3071
3077 const double computeGammaHgaga_WW() const
3078 {
3079 return 14731.86; // in eV for Mh=125 GeV
3080 //return 15221.98; // in eV for Mh=126 GeV
3081 }
3082
3088 const double computeGammaHgaga_tW() const
3089 {
3090 return -6249.93; // in eV for Mh=125 GeV
3091 //return -6436.35; // in eV for Mh=126 GeV
3092 }
3093
3098 virtual const double getCBd() const
3099 {
3100 return 1.;
3101 }
3102
3107 virtual const double getCBs() const
3108 {
3109 return 1.;
3110 }
3111
3116 virtual const double getCDMK() const
3117 {
3118 return 1.;
3119 }
3120
3125 virtual const double getCepsK() const
3126 {
3127 return 1.;
3128 }
3129
3134 virtual const double getPhiBs() const
3135 {
3136 return 0.;
3137 }
3138
3143 virtual const double getPhiBd() const
3144 {
3145 return 0.;
3146 }
3147
3149 //Generic e+e- -> ff Inclusive Observables
3150
3151// For f!=e
3152
3153// Helicity amplitudes squared
3154 const double MLR2eeff(const Particle f, const double s) const;
3155 const double MRL2eeff(const Particle f, const double s) const;
3156
3157 const double MLL2eeff(const Particle f, const double s, const double t) const;
3158 const double MRR2eeff(const Particle f, const double s, const double t) const;
3159
3160// Some simple functions for cos \theta integrals
3161 const double tovers2(const double cosmin, const double cosmax) const;
3162 const double uovers2(const double cosmin, const double cosmax) const;
3163
3164// For f=e
3165
3166// Integrals of the SM squared amplitudes x (t/s)^2, (s/t)^2, (u/s)^2 in [t0, t1]
3167 const double intMLR2eeeets2(const double s, const double t0, const double t1) const;
3168 const double intMLRtilde2eeeest2(const double s, const double t0, const double t1) const;
3169 const double intMLL2eeeeus2(const double s, const double t0, const double t1) const;
3170 const double intMRR2eeeeus2(const double s, const double t0, const double t1) const;
3171
3172// Cross sections
3173
3174 const double eeffsigmaEbin(const double pol_e, const double pol_p, const double s, const double cosmin, const double cosmax) const;
3175
3176 virtual const double eeffsigmaE(const double pol_e, const double pol_p, const double s) const
3177 {
3178 double cosmin = -0.90; // As in LEP2
3179 double cosmax = 0.90; // As in LEP2
3180
3181 return eeffsigmaEbin(pol_e, pol_p, s, cosmin, cosmax);
3182 }
3183
3184 virtual const double eeffsigma(const Particle f, const double pol_e, const double pol_p, const double s, const double cosmin, const double cosmax) const;
3185
3186 virtual const double eeffsigmaEtsub(const double pol_e, const double pol_p, const double s) const
3187 {
3188 double cosmin = -0.90; // As in LEP2
3189 double cosmax = 0.90; // As in LEP2
3190
3191 return eeffsigma(leptons[ELECTRON], pol_e, pol_p, s, cosmin, cosmax);
3192 }
3193
3194 virtual const double eeffsigmaMu(const double pol_e, const double pol_p, const double s) const
3195 {
3196 return eeffsigma(leptons[MU], pol_e, pol_p, s, -1.0, 1.0);
3197 }
3198
3199 virtual const double eeffsigmaTau(const double pol_e, const double pol_p, const double s) const
3200 {
3201 return eeffsigma(leptons[TAU], pol_e, pol_p, s, -1.0, 1.0);
3202 }
3203
3204 virtual const double eeffsigmaHadron(const double pol_e, const double pol_p, const double s) const
3205 {
3206 return ( eeffsigma(quarks[UP], pol_e, pol_p, s, -1.0, 1.0) + eeffsigma(quarks[DOWN], pol_e, pol_p, s, -1.0, 1.0)
3207 + eeffsigma(quarks[CHARM], pol_e, pol_p, s, -1.0, 1.0) + eeffsigma(quarks[STRANGE], pol_e, pol_p, s, -1.0, 1.0)
3208 + eeffsigma(quarks[BOTTOM], pol_e, pol_p, s, -1.0, 1.0) );
3209 }
3210
3211 virtual const double eeffsigmaStrange(const double pol_e, const double pol_p, const double s) const
3212 {
3213 return eeffsigma(quarks[STRANGE], pol_e, pol_p, s, -1.0, 1.0);
3214 }
3215
3216 virtual const double eeffsigmaCharm(const double pol_e, const double pol_p, const double s) const
3217 {
3218 return eeffsigma(quarks[CHARM], pol_e, pol_p, s, -1.0, 1.0);
3219 }
3220
3221 virtual const double eeffsigmaBottom(const double pol_e, const double pol_p, const double s) const
3222 {
3223 return eeffsigma(quarks[BOTTOM], pol_e, pol_p, s, -1.0, 1.0);
3224 }
3225
3226// Ratios
3227
3228 virtual const double eeffRelectron(const double pol_e, const double pol_p, const double s) const
3229 {
3230 return ( eeffsigmaHadron(pol_e, pol_p, s) / eeffsigmaE(pol_e, pol_p, s) );
3231 }
3232 virtual const double eeffRelectrontsub(const double pol_e, const double pol_p, const double s) const
3233 {
3234 return ( eeffsigmaHadron(pol_e, pol_p, s) / eeffsigmaEtsub(pol_e, pol_p, s) );
3235 }
3236 virtual const double eeffRmuon(const double pol_e, const double pol_p, const double s) const
3237 {
3238 return ( eeffsigmaHadron(pol_e, pol_p, s) / eeffsigmaMu(pol_e, pol_p, s) );
3239 }
3240 virtual const double eeffRtau(const double pol_e, const double pol_p, const double s) const
3241 {
3242 return ( eeffsigmaHadron(pol_e, pol_p, s) / eeffsigmaTau(pol_e, pol_p, s) );
3243 }
3244 virtual const double eeffRstrange(const double pol_e, const double pol_p, const double s) const
3245 {
3246 return ( eeffsigmaStrange(pol_e, pol_p, s) / eeffsigmaHadron(pol_e, pol_p, s) );
3247 }
3248 virtual const double eeffRcharm(const double pol_e, const double pol_p, const double s) const
3249 {
3250 return ( eeffsigmaCharm(pol_e, pol_p, s) / eeffsigmaHadron(pol_e, pol_p, s) );
3251 }
3252 virtual const double eeffRbottom(const double pol_e, const double pol_p, const double s) const
3253 {
3254 return ( eeffsigmaBottom(pol_e, pol_p, s) / eeffsigmaHadron(pol_e, pol_p, s) );
3255 }
3256
3257// FB asymmetries
3258
3259 virtual const double eeffAFBe(const double pol_e, const double pol_p, const double s) const
3260 {
3261 double cosmin = -0.90; // As in LEP2
3262 double cosmax = 0.90; // As in LEP2
3263
3264 return (( eeffsigmaEbin(pol_e, pol_p, s, 0.0 , cosmax) - eeffsigmaEbin(pol_e, pol_p, s, cosmin, 0.0) ) / eeffsigmaEbin(pol_e, pol_p, s, cosmin, cosmax));
3265 }
3266 virtual const double eeffAFBetsub(const double pol_e, const double pol_p, const double s) const
3267 {
3268 double cosmin = -0.90; // As in LEP2
3269 double cosmax = 0.90; // As in LEP2
3270
3271 return ( ( eeffsigma(leptons[ELECTRON], pol_e, pol_p, s, 0.0 , cosmax) - eeffsigma(leptons[ELECTRON], pol_e, pol_p, s, cosmin, 0.0) ) / eeffsigma(leptons[ELECTRON], pol_e, pol_p, s, cosmin, cosmax) );
3272 }
3273 virtual const double eeffAFBmu(const double pol_e, const double pol_p, const double s) const
3274 {
3275 return ( ( eeffsigma(leptons[MU], pol_e, pol_p, s, 0.0, 1.0) - eeffsigma(leptons[MU], pol_e, pol_p, s, -1.0, 0.0) ) / eeffsigma(leptons[MU], pol_e, pol_p, s, -1.0, 1.0) );
3276 }
3277 virtual const double eeffAFBtau(const double pol_e, const double pol_p, const double s) const
3278 {
3279 return ( ( eeffsigma(leptons[TAU], pol_e, pol_p, s, 0.0, 1.0) - eeffsigma(leptons[TAU], pol_e, pol_p, s, -1.0, 0.0) ) / eeffsigma(leptons[TAU], pol_e, pol_p, s, -1.0, 1.0) );
3280 }
3281 virtual const double eeffAFBstrange(const double pol_e, const double pol_p, const double s) const
3282 {
3283 return ( ( eeffsigma(quarks[STRANGE], pol_e, pol_p, s, 0.0, 1.0) - eeffsigma(quarks[STRANGE], pol_e, pol_p, s, -1.0, 0.0) ) / eeffsigma(quarks[STRANGE], pol_e, pol_p, s, -1.0, 1.0) );
3284 }
3285 virtual const double eeffAFBcharm(const double pol_e, const double pol_p, const double s) const
3286 {
3287 return ( ( eeffsigma(quarks[CHARM], pol_e, pol_p, s, 0.0, 1.0) - eeffsigma(quarks[CHARM], pol_e, pol_p, s, -1.0, 0.0) ) / eeffsigma(quarks[CHARM], pol_e, pol_p, s, -1.0, 1.0) );
3288 }
3289 virtual const double eeffAFBbottom(const double pol_e, const double pol_p, const double s) const
3290 {
3291 return ( ( eeffsigma(quarks[BOTTOM], pol_e, pol_p, s, 0.0, 1.0) - eeffsigma(quarks[BOTTOM], pol_e, pol_p, s, -1.0, 0.0) ) / eeffsigma(quarks[BOTTOM], pol_e, pol_p, s, -1.0, 1.0) );
3292 }
3293
3294/* BEGIN: REMOVE FROM THE PACKAGE */
3296 //LEP2 Inclusive Observables
3297
3298 virtual const double LEP2sigmaE(const double s) const;
3299 virtual const double LEP2sigmaMu(const double s) const;
3300 virtual const double LEP2sigmaTau(const double s) const;
3301 virtual const double LEP2sigmaHadron(const double s) const;
3302 virtual const double LEP2sigmaCharm(const double s) const;
3303 virtual const double LEP2sigmaBottom(const double s) const;
3304
3305 virtual const double LEP2AFBe(const double s) const;
3306 virtual const double LEP2AFBmu(const double s) const;
3307 virtual const double LEP2AFBtau(const double s) const;
3308 virtual const double LEP2AFBcharm(const double s) const;
3309 virtual const double LEP2AFBbottom(const double s) const;
3310 virtual const double LEP2Rcharm(const double s) const;
3311 virtual const double LEP2Rbottom(const double s) const;
3312
3313 //LEP2 Differential Observables
3314 virtual const double LEP2dsigmadcosE(const double s, const double cos) const;
3315 virtual const double LEP2dsigmadcosMu(const double s, const double cos) const;
3316 virtual const double LEP2dsigmadcosTau(const double s, const double cos) const;
3317
3318 virtual const double LEP2dsigmadcosBinE(const double s, const double cos, const double cosmin, const double cosmax) const;
3319 virtual const double LEP2dsigmadcosBinMu(const double s, const double cos, const double cosmin, const double cosmax) const;
3320 virtual const double LEP2dsigmadcosBinTau(const double s, const double cos, const double cosmin, const double cosmax) const;
3321
3322/* END: REMOVE FROM THE PACKAGE */
3323
3324bool setFlagSigmaForAFB(const bool flagSigmaForAFB_i)
3325{
3326 bSigmaForAFB = flagSigmaForAFB_i;
3327 return true;
3328}
3329
3330bool setFlagSigmaForR(const bool flagSigmaForR_i)
3331{
3332 bSigmaForR = flagSigmaForR_i;
3333 return true;
3334}
3335
3342virtual const double getmq(const QCD::quark q, const double mu) const
3343{
3344 return m_q(q, mu, FULLNLO);
3345}
3346
3352 {
3353 this->requireCKM = requireCKM;
3354 }
3355
3360 void setCKM(const CKM& CKMMatrix)
3361 {
3362 myCKM = CKMMatrix;
3363 }
3364
3370 const gslpp::matrix<gslpp::complex>& getYu() const
3371 {
3372 return Yu;
3373 }
3374
3380 void setYu(const gslpp::matrix<gslpp::complex>& Yu)
3381 {
3382 this->Yu = Yu;
3383 }
3384
3390 const gslpp::matrix<gslpp::complex>& getYd() const
3391 {
3392 return Yd;
3393 }
3394
3400 void setYd(const gslpp::matrix<gslpp::complex>& Yd)
3401 {
3402 this->Yd = Yd;
3403 }
3404
3410 const gslpp::matrix<gslpp::complex>& getYe() const
3411 {
3412 return Ye;
3413 }
3414
3420 void setYe(const gslpp::matrix<gslpp::complex>& Ye)
3421 {
3422 this->Ye = Ye;
3423 }
3424
3429 const bool isSMSuccess() const
3430 {
3431 return SMSuccess;
3432 }
3433
3438 void setSMSuccess(bool success) const
3439 {
3440 SMSuccess = success;
3441 }
3442
3443
3445protected:
3446
3452 virtual void setParameter(const std::string name, const double& value);
3453
3457 virtual void computeCKM();
3458
3462 virtual void computeYukawas();
3463
3467// gslpp::matrix<gslpp::complex> VCKM; ///< The %CKM matrix.
3468// gslpp::matrix<gslpp::complex> UPMNS; ///< The %PMNS matrix.
3470 gslpp::matrix<gslpp::complex> Yu;
3471 gslpp::matrix<gslpp::complex> Yd;
3472 gslpp::matrix<gslpp::complex> Yn;
3473 gslpp::matrix<gslpp::complex> Ye;
3474
3475
3476 // model parameters
3477 double AlsMz;
3478 double Mz;
3479 double Mw_inp;
3480 double GF;
3481 double ale;
3482 double dAle5Mz;
3483 double mHl;
3484 double delMw;
3488 double delGammaZ;
3489 double delsigma0H;
3490 double delR0l;
3491 double delR0c;
3492 double delR0b;
3495 double lambda;
3496 double A;
3497 double rhob;
3498 double etab;
3499 double Vus;
3500 double Vud;
3501 double Vcb;
3502 double Vub;
3503 double gamma;
3504 double muw;
3506
3507 // non-model parameters
3508 double dAl5hMz;
3509
3511 // For EWPO
3512
3520
3531 double SchemeToDouble(const std::string scheme) const
3532 {
3533 if (scheme.compare("NORESUM") == 0)
3534 return 0.0;
3535 else if (scheme.compare("OMSI") == 0)
3536 return 1.0;
3537 else if (scheme.compare("INTERMEDIATE") == 0)
3538 return 2.0;
3539 else if (scheme.compare("OMSII") == 0)
3540 return 3.0;
3541 else if (scheme.compare("APPROXIMATEFORMULA") == 0)
3542 return 4.0;
3543 else
3544 throw std::runtime_error("EWSM::SchemeToDouble: bad scheme");
3545 }
3546
3552 bool checkEWPOscheme(const std::string scheme) const
3553 {
3554 if (scheme.compare("NORESUM") == 0
3555 || scheme.compare("OMSI") == 0
3556 || scheme.compare("INTERMEDIATE") == 0
3557 || scheme.compare("OMSII") == 0
3558 || scheme.compare("APPROXIMATEFORMULA") == 0)
3559 return true;
3560 else
3561 return false;
3562 }
3563
3564
3565 double m_q(const QCD::quark q, const double mu, const orders order=FULLNLO) const
3566 {
3567 switch(q) {
3568 case QCD::UP:
3569 case QCD::DOWN:
3570 case QCD::STRANGE:
3571 return Mrun(mu, getQuarks(q).getMass_scale(),
3572 getQuarks(q).getMass(), q, order);
3573 case QCD::CHARM:
3574 case QCD::BOTTOM:
3575 case QCD::TOP:
3576 return Mrun(mu, getQuarks(q).getMass(), q, order);
3577 default:
3578 throw std::runtime_error("Error in StandardModel::m_q()");
3579 }
3580 }
3581
3614 double resumMw(const double Mw_i, const double DeltaRho[orders_EW_size],
3615 const double DeltaR_rem[orders_EW_size]) const;
3616
3643 double resumRhoZ(const double DeltaRho[orders_EW_size],
3644 const double deltaRho_rem[orders_EW_size],
3645 const double DeltaRbar_rem, const bool bool_Zbb) const;
3646
3673 double resumKappaZ(const double DeltaRho[orders_EW_size],
3674 const double deltaKappa_rem[orders_EW_size],
3675 const double DeltaRbar_rem, const bool bool_Zbb) const;
3676
3699 double taub() const;
3700
3709 double Delta_EWQCD(const QCD::quark q) const;
3710
3720 double RVq(const QCD::quark q) const;
3721
3731 double RAq(const QCD::quark q) const;
3732
3746 double RVh() const;
3747
3751
3754
3755
3756 /* BEGIN: REMOVE FROM THE PACKAGE */
3758 //Migrated from LEP2ThObservables.h
3759
3760
3761
3763 mutable bool bSigmaForAFB;
3764 mutable bool bSigmaForR;
3765
3766
3767
3768 const double sigma_NoISR_l(const QCD::lepton l_flavor, const double s) const;
3769 const double sigma_NoISR_q(const QCD::quark q_flavor, const double s) const;
3770 const double AFB_NoISR_l(const QCD::lepton l_flavor, const double s) const;
3771 const double AFB_NoISR_q(const QCD::quark q_flavor, const double s) const;
3772
3773 const double Integrand_sigmaWithISR_l(double x, const QCD::lepton l_flavor, const double s) const;
3774
3775 const double getIntegrand_sigmaWithISR_mu130(double x) const;
3776 const double getIntegrand_sigmaWithISR_mu136(double x) const;
3777 const double getIntegrand_sigmaWithISR_mu161(double x) const;
3778 const double getIntegrand_sigmaWithISR_mu172(double x) const;
3779 const double getIntegrand_sigmaWithISR_mu183(double x) const;
3780 const double getIntegrand_sigmaWithISR_mu189(double x) const;
3781 const double getIntegrand_sigmaWithISR_mu192(double x) const;
3782 const double getIntegrand_sigmaWithISR_mu196(double x) const;
3783 const double getIntegrand_sigmaWithISR_mu200(double x) const;
3784 const double getIntegrand_sigmaWithISR_mu202(double x) const;
3785 const double getIntegrand_sigmaWithISR_mu205(double x) const;
3786 const double getIntegrand_sigmaWithISR_mu207(double x) const;
3787
3788
3789 const double getIntegrand_sigmaWithISR_tau130(double x) const;
3790 const double getIntegrand_sigmaWithISR_tau136(double x) const;
3791 const double getIntegrand_sigmaWithISR_tau161(double x) const;
3792 const double getIntegrand_sigmaWithISR_tau172(double x) const;
3793 const double getIntegrand_sigmaWithISR_tau183(double x) const;
3794 const double getIntegrand_sigmaWithISR_tau189(double x) const;
3795 const double getIntegrand_sigmaWithISR_tau192(double x) const;
3796 const double getIntegrand_sigmaWithISR_tau196(double x) const;
3797 const double getIntegrand_sigmaWithISR_tau200(double x) const;
3798 const double getIntegrand_sigmaWithISR_tau202(double x) const;
3799 const double getIntegrand_sigmaWithISR_tau205(double x) const;
3800 const double getIntegrand_sigmaWithISR_tau207(double x) const;
3801
3802
3803 const double Integrand_sigmaWithISR_q(double x, const QCD::quark q_flavor, const double s) const;
3804
3805 const double getIntegrand_sigmaWithISR_up130(double x) const;
3806 const double getIntegrand_sigmaWithISR_up133(double x) const;
3807 const double getIntegrand_sigmaWithISR_up136(double x) const;
3808 const double getIntegrand_sigmaWithISR_up161(double x) const;
3809 const double getIntegrand_sigmaWithISR_up167(double x) const;
3810 const double getIntegrand_sigmaWithISR_up172(double x) const;
3811 const double getIntegrand_sigmaWithISR_up183(double x) const;
3812 const double getIntegrand_sigmaWithISR_up189(double x) const;
3813 const double getIntegrand_sigmaWithISR_up192(double x) const;
3814 const double getIntegrand_sigmaWithISR_up196(double x) const;
3815 const double getIntegrand_sigmaWithISR_up200(double x) const;
3816 const double getIntegrand_sigmaWithISR_up202(double x) const;
3817 const double getIntegrand_sigmaWithISR_up205(double x) const;
3818 const double getIntegrand_sigmaWithISR_up207(double x) const;
3819
3820 const double getIntegrand_sigmaWithISR_down130(double x) const;
3821 const double getIntegrand_sigmaWithISR_down133(double x) const;
3822 const double getIntegrand_sigmaWithISR_down136(double x) const;
3823 const double getIntegrand_sigmaWithISR_down161(double x) const;
3824 const double getIntegrand_sigmaWithISR_down167(double x) const;
3825 const double getIntegrand_sigmaWithISR_down172(double x) const;
3826 const double getIntegrand_sigmaWithISR_down183(double x) const;
3827 const double getIntegrand_sigmaWithISR_down189(double x) const;
3828 const double getIntegrand_sigmaWithISR_down192(double x) const;
3829 const double getIntegrand_sigmaWithISR_down196(double x) const;
3830 const double getIntegrand_sigmaWithISR_down200(double x) const;
3831 const double getIntegrand_sigmaWithISR_down202(double x) const;
3832 const double getIntegrand_sigmaWithISR_down205(double x) const;
3833 const double getIntegrand_sigmaWithISR_down207(double x) const;
3834
3835 const double getIntegrand_sigmaWithISR_charm130(double x) const;
3836 const double getIntegrand_sigmaWithISR_charm133(double x) const;
3837 const double getIntegrand_sigmaWithISR_charm136(double x) const;
3838 const double getIntegrand_sigmaWithISR_charm161(double x) const;
3839 const double getIntegrand_sigmaWithISR_charm167(double x) const;
3840 const double getIntegrand_sigmaWithISR_charm172(double x) const;
3841 const double getIntegrand_sigmaWithISR_charm183(double x) const;
3842 const double getIntegrand_sigmaWithISR_charm189(double x) const;
3843 const double getIntegrand_sigmaWithISR_charm192(double x) const;
3844 const double getIntegrand_sigmaWithISR_charm196(double x) const;
3845 const double getIntegrand_sigmaWithISR_charm200(double x) const;
3846 const double getIntegrand_sigmaWithISR_charm202(double x) const;
3847 const double getIntegrand_sigmaWithISR_charm205(double x) const;
3848 const double getIntegrand_sigmaWithISR_charm207(double x) const;
3849
3850 const double getIntegrand_sigmaWithISR_strange130(double x) const;
3851 const double getIntegrand_sigmaWithISR_strange133(double x) const;
3852 const double getIntegrand_sigmaWithISR_strange136(double x) const;
3853 const double getIntegrand_sigmaWithISR_strange161(double x) const;
3854 const double getIntegrand_sigmaWithISR_strange167(double x) const;
3855 const double getIntegrand_sigmaWithISR_strange172(double x) const;
3856 const double getIntegrand_sigmaWithISR_strange183(double x) const;
3857 const double getIntegrand_sigmaWithISR_strange189(double x) const;
3858 const double getIntegrand_sigmaWithISR_strange192(double x) const;
3859 const double getIntegrand_sigmaWithISR_strange196(double x) const;
3860 const double getIntegrand_sigmaWithISR_strange200(double x) const;
3861 const double getIntegrand_sigmaWithISR_strange202(double x) const;
3862 const double getIntegrand_sigmaWithISR_strange205(double x) const;
3863 const double getIntegrand_sigmaWithISR_strange207(double x) const;
3864
3865 const double getIntegrand_sigmaWithISR_bottom130(double x) const;
3866 const double getIntegrand_sigmaWithISR_bottom133(double x) const;
3867 const double getIntegrand_sigmaWithISR_bottom136(double x) const;
3868 const double getIntegrand_sigmaWithISR_bottom161(double x) const;
3869 const double getIntegrand_sigmaWithISR_bottom167(double x) const;
3870 const double getIntegrand_sigmaWithISR_bottom172(double x) const;
3871 const double getIntegrand_sigmaWithISR_bottom183(double x) const;
3872 const double getIntegrand_sigmaWithISR_bottom189(double x) const;
3873 const double getIntegrand_sigmaWithISR_bottom192(double x) const;
3874 const double getIntegrand_sigmaWithISR_bottom196(double x) const;
3875 const double getIntegrand_sigmaWithISR_bottom200(double x) const;
3876 const double getIntegrand_sigmaWithISR_bottom202(double x) const;
3877 const double getIntegrand_sigmaWithISR_bottom205(double x) const;
3878 const double getIntegrand_sigmaWithISR_bottom207(double x) const;
3879
3880 const double Integrand_dsigmaBox_l(double cosTheta, const QCD::lepton l_flavor, const double s) const;
3881
3882 const double getIntegrand_dsigmaBox_mu130(double x) const;
3883 const double getIntegrand_dsigmaBox_mu133(double x) const;
3884 const double getIntegrand_dsigmaBox_mu136(double x) const;
3885 const double getIntegrand_dsigmaBox_mu161(double x) const;
3886 const double getIntegrand_dsigmaBox_mu167(double x) const;
3887 const double getIntegrand_dsigmaBox_mu172(double x) const;
3888 const double getIntegrand_dsigmaBox_mu183(double x) const;
3889 const double getIntegrand_dsigmaBox_mu189(double x) const;
3890 const double getIntegrand_dsigmaBox_mu192(double x) const;
3891 const double getIntegrand_dsigmaBox_mu196(double x) const;
3892 const double getIntegrand_dsigmaBox_mu200(double x) const;
3893 const double getIntegrand_dsigmaBox_mu202(double x) const;
3894 const double getIntegrand_dsigmaBox_mu205(double x) const;
3895 const double getIntegrand_dsigmaBox_mu207(double x) const;
3896
3897 const double getIntegrand_dsigmaBox_tau130(double x) const;
3898 const double getIntegrand_dsigmaBox_tau133(double x) const;
3899 const double getIntegrand_dsigmaBox_tau136(double x) const;
3900 const double getIntegrand_dsigmaBox_tau161(double x) const;
3901 const double getIntegrand_dsigmaBox_tau167(double x) const;
3902 const double getIntegrand_dsigmaBox_tau172(double x) const;
3903 const double getIntegrand_dsigmaBox_tau183(double x) const;
3904 const double getIntegrand_dsigmaBox_tau189(double x) const;
3905 const double getIntegrand_dsigmaBox_tau192(double x) const;
3906 const double getIntegrand_dsigmaBox_tau196(double x) const;
3907 const double getIntegrand_dsigmaBox_tau200(double x) const;
3908 const double getIntegrand_dsigmaBox_tau202(double x) const;
3909 const double getIntegrand_dsigmaBox_tau205(double x) const;
3910 const double getIntegrand_dsigmaBox_tau207(double x) const;
3911
3912
3913 const double Integrand_dsigmaBox_q(double cosTheta, const QCD::quark q_flavor, const double s) const;
3914
3915 const double getIntegrand_dsigmaBox_up130(double x) const;
3916 const double getIntegrand_dsigmaBox_up133(double x) const;
3917 const double getIntegrand_dsigmaBox_up136(double x) const;
3918 const double getIntegrand_dsigmaBox_up161(double x) const;
3919 const double getIntegrand_dsigmaBox_up167(double x) const;
3920 const double getIntegrand_dsigmaBox_up172(double x) const;
3921 const double getIntegrand_dsigmaBox_up183(double x) const;
3922 const double getIntegrand_dsigmaBox_up189(double x) const;
3923 const double getIntegrand_dsigmaBox_up192(double x) const;
3924 const double getIntegrand_dsigmaBox_up196(double x) const;
3925 const double getIntegrand_dsigmaBox_up200(double x) const;
3926 const double getIntegrand_dsigmaBox_up202(double x) const;
3927 const double getIntegrand_dsigmaBox_up205(double x) const;
3928 const double getIntegrand_dsigmaBox_up207(double x) const;
3929
3930 const double getIntegrand_dsigmaBox_down130(double x) const;
3931 const double getIntegrand_dsigmaBox_down133(double x) const;
3932 const double getIntegrand_dsigmaBox_down136(double x) const;
3933 const double getIntegrand_dsigmaBox_down161(double x) const;
3934 const double getIntegrand_dsigmaBox_down167(double x) const;
3935 const double getIntegrand_dsigmaBox_down172(double x) const;
3936 const double getIntegrand_dsigmaBox_down183(double x) const;
3937 const double getIntegrand_dsigmaBox_down189(double x) const;
3938 const double getIntegrand_dsigmaBox_down192(double x) const;
3939 const double getIntegrand_dsigmaBox_down196(double x) const;
3940 const double getIntegrand_dsigmaBox_down200(double x) const;
3941 const double getIntegrand_dsigmaBox_down202(double x) const;
3942 const double getIntegrand_dsigmaBox_down205(double x) const;
3943 const double getIntegrand_dsigmaBox_down207(double x) const;
3944
3945 const double getIntegrand_dsigmaBox_charm130(double x) const;
3946 const double getIntegrand_dsigmaBox_charm133(double x) const;
3947 const double getIntegrand_dsigmaBox_charm136(double x) const;
3948 const double getIntegrand_dsigmaBox_charm161(double x) const;
3949 const double getIntegrand_dsigmaBox_charm167(double x) const;
3950 const double getIntegrand_dsigmaBox_charm172(double x) const;
3951 const double getIntegrand_dsigmaBox_charm183(double x) const;
3952 const double getIntegrand_dsigmaBox_charm189(double x) const;
3953 const double getIntegrand_dsigmaBox_charm192(double x) const;
3954 const double getIntegrand_dsigmaBox_charm196(double x) const;
3955 const double getIntegrand_dsigmaBox_charm200(double x) const;
3956 const double getIntegrand_dsigmaBox_charm202(double x) const;
3957 const double getIntegrand_dsigmaBox_charm205(double x) const;
3958 const double getIntegrand_dsigmaBox_charm207(double x) const;
3959
3960 const double getIntegrand_dsigmaBox_strange130(double x) const;
3961 const double getIntegrand_dsigmaBox_strange133(double x) const;
3962 const double getIntegrand_dsigmaBox_strange136(double x) const;
3963 const double getIntegrand_dsigmaBox_strange161(double x) const;
3964 const double getIntegrand_dsigmaBox_strange167(double x) const;
3965 const double getIntegrand_dsigmaBox_strange172(double x) const;
3966 const double getIntegrand_dsigmaBox_strange183(double x) const;
3967 const double getIntegrand_dsigmaBox_strange189(double x) const;
3968 const double getIntegrand_dsigmaBox_strange192(double x) const;
3969 const double getIntegrand_dsigmaBox_strange196(double x) const;
3970 const double getIntegrand_dsigmaBox_strange200(double x) const;
3971 const double getIntegrand_dsigmaBox_strange202(double x) const;
3972 const double getIntegrand_dsigmaBox_strange205(double x) const;
3973 const double getIntegrand_dsigmaBox_strange207(double x) const;
3974
3975 const double getIntegrand_dsigmaBox_bottom130(double x) const;
3976 const double getIntegrand_dsigmaBox_bottom133(double x) const;
3977 const double getIntegrand_dsigmaBox_bottom136(double x) const;
3978 const double getIntegrand_dsigmaBox_bottom161(double x) const;
3979 const double getIntegrand_dsigmaBox_bottom167(double x) const;
3980 const double getIntegrand_dsigmaBox_bottom172(double x) const;
3981 const double getIntegrand_dsigmaBox_bottom183(double x) const;
3982 const double getIntegrand_dsigmaBox_bottom189(double x) const;
3983 const double getIntegrand_dsigmaBox_bottom192(double x) const;
3984 const double getIntegrand_dsigmaBox_bottom196(double x) const;
3985 const double getIntegrand_dsigmaBox_bottom200(double x) const;
3986 const double getIntegrand_dsigmaBox_bottom202(double x) const;
3987 const double getIntegrand_dsigmaBox_bottom205(double x) const;
3988 const double getIntegrand_dsigmaBox_bottom207(double x) const;
3989
3990
3991 const double Integrand_AFBnumeratorWithISR_l(double x, const QCD::lepton l_flavor, const double s) const;
3992
3993 const double getIntegrand_AFBnumeratorWithISR_mu130(double x) const;
3994 const double getIntegrand_AFBnumeratorWithISR_mu136(double x) const;
3995 const double getIntegrand_AFBnumeratorWithISR_mu161(double x) const;
3996 const double getIntegrand_AFBnumeratorWithISR_mu172(double x) const;
3997 const double getIntegrand_AFBnumeratorWithISR_mu183(double x) const;
3998 const double getIntegrand_AFBnumeratorWithISR_mu189(double x) const;
3999 const double getIntegrand_AFBnumeratorWithISR_mu192(double x) const;
4000 const double getIntegrand_AFBnumeratorWithISR_mu196(double x) const;
4001 const double getIntegrand_AFBnumeratorWithISR_mu200(double x) const;
4002 const double getIntegrand_AFBnumeratorWithISR_mu202(double x) const;
4003 const double getIntegrand_AFBnumeratorWithISR_mu205(double x) const;
4004 const double getIntegrand_AFBnumeratorWithISR_mu207(double x) const;
4005
4006 const double getIntegrand_AFBnumeratorWithISR_tau130(double x) const;
4007 const double getIntegrand_AFBnumeratorWithISR_tau136(double x) const;
4008 const double getIntegrand_AFBnumeratorWithISR_tau161(double x) const;
4009 const double getIntegrand_AFBnumeratorWithISR_tau172(double x) const;
4010 const double getIntegrand_AFBnumeratorWithISR_tau183(double x) const;
4011 const double getIntegrand_AFBnumeratorWithISR_tau189(double x) const;
4012 const double getIntegrand_AFBnumeratorWithISR_tau192(double x) const;
4013 const double getIntegrand_AFBnumeratorWithISR_tau196(double x) const;
4014 const double getIntegrand_AFBnumeratorWithISR_tau200(double x) const;
4015 const double getIntegrand_AFBnumeratorWithISR_tau202(double x) const;
4016 const double getIntegrand_AFBnumeratorWithISR_tau205(double x) const;
4017 const double getIntegrand_AFBnumeratorWithISR_tau207(double x) const;
4018
4019 const double Integrand_AFBnumeratorWithISR_q(double x, const QCD::quark q_flavor, const double s) const;
4020
4021 const double getIntegrand_AFBnumeratorWithISR_charm133(double x) const;
4022 const double getIntegrand_AFBnumeratorWithISR_charm167(double x) const;
4023 const double getIntegrand_AFBnumeratorWithISR_charm172(double x) const;
4024 const double getIntegrand_AFBnumeratorWithISR_charm183(double x) const;
4025 const double getIntegrand_AFBnumeratorWithISR_charm189(double x) const;
4026 const double getIntegrand_AFBnumeratorWithISR_charm192(double x) const;
4027 const double getIntegrand_AFBnumeratorWithISR_charm196(double x) const;
4028 const double getIntegrand_AFBnumeratorWithISR_charm200(double x) const;
4029 const double getIntegrand_AFBnumeratorWithISR_charm202(double x) const;
4030 const double getIntegrand_AFBnumeratorWithISR_charm205(double x) const;
4031 const double getIntegrand_AFBnumeratorWithISR_charm207(double x) const;
4032
4033 const double getIntegrand_AFBnumeratorWithISR_bottom133(double x) const;
4034 const double getIntegrand_AFBnumeratorWithISR_bottom167(double x) const;
4035 const double getIntegrand_AFBnumeratorWithISR_bottom172(double x) const;
4036 const double getIntegrand_AFBnumeratorWithISR_bottom183(double x) const;
4037 const double getIntegrand_AFBnumeratorWithISR_bottom189(double x) const;
4038 const double getIntegrand_AFBnumeratorWithISR_bottom192(double x) const;
4039 const double getIntegrand_AFBnumeratorWithISR_bottom196(double x) const;
4040 const double getIntegrand_AFBnumeratorWithISR_bottom200(double x) const;
4041 const double getIntegrand_AFBnumeratorWithISR_bottom202(double x) const;
4042 const double getIntegrand_AFBnumeratorWithISR_bottom205(double x) const;
4043 const double getIntegrand_AFBnumeratorWithISR_bottom207(double x) const;
4044
4045 /* END: REMOVE FROM THE PACKAGE */
4047private:
4057 /* BEGIN: REMOVE FROM THE PACKAGE */
4059 /* END: REMOVE FROM THE PACKAGE */
4060
4063 std::string FlagMw;
4064 std::string FlagRhoZ;
4065 std::string FlagKappaZ;
4068
4071
4072 mutable bool SMSuccess;
4073
4075 // Caches for EWPO
4076
4079 mutable double DeltaAlphaLepton_cache;
4080 mutable double DeltaAlpha_cache;
4081 mutable double Mw_cache;
4082 mutable double GammaW_cache;
4083 mutable gslpp::complex rhoZ_f_cache[12];
4084 mutable gslpp::complex kappaZ_f_cache[12];
4087 mutable bool useMw_cache;
4088 mutable bool useGammaW_cache;
4089 mutable bool useRhoZ_f_cache[12];
4090 mutable bool useKappaZ_f_cache[12];
4091
4092/* BEGIN: REMOVE FROM THE PACKAGE */
4093 // caches for the SM prediction of LEP2 Obs
4094// mutable double SMparams_cache[NumSMParamsForEWPO+3];
4095 mutable double SMresult_cache;
4096// mutable bool flag_cache[NUMofLEP2RCs];
4097// mutable double ml_cache, mq_cache, mqForHad_cache[6];
4098
4099 mutable double average;
4100 mutable double error;
4101 mutable gsl_function f_GSL;
4102 gsl_integration_workspace * w_GSL1;
4103/* END: REMOVE FROM THE PACKAGE */
4104
4106
4107 const double AlsE(double mu, orders order, bool Nf_thr) const;
4108 const double AlsEByOrder(double mu, orders order, bool Nf_thr) const;
4109
4110 const double AlsEWithInit(double mu, double alsi, double mu_i, const int nf_i, orders order) const;
4111 const double AleWithInit(double mu, double alsi, double mu_i, orders order) const;
4112 static const int CacheSize = 5;
4113 mutable double als_cache[11][CacheSize];
4114 mutable double ale_cache[10][CacheSize];
4116
4117
4118};
4119
4120#endif /* STANDARDMODEL_H */
std::map< std::string, double > DPars
Definition: Minimal.cpp:11
@ FULLNNNLO
Definition: OrderScheme.h:40
@ FULLNLO
Definition: OrderScheme.h:38
A class for the CKM matrix elements.
Definition: CKM.h:23
const gslpp::matrix< gslpp::complex > getCKM() const
A member for returning the CKM matrix.
Definition: CKM.h:59
A class for approximate formulae of the EW precision observables.
A class for one-loop corrections to the EW precision observables.
A class for three-loop corrections to the EW precision observables.
A class for three-loop corrections to the EW precision observables.
A class for three-loop corrections to the EW precision observables.
A class for the form factors , and in the processes at LEP-II.
A class for two-loop corrections to the EW precision observables.
Definition: EWSMTwoLoopEW.h:57
A class for two-loop corrections to the EW precision observables.
A class for cache variables used in computing radiative corrections to the EW precision observables.
Definition: EWSMcache.h:40
The parent class in Flavour for calculating all the Wilson coefficients for various Flavor Violating ...
Definition: Flavour.h:37
The parent class in LeptonFlavour for calculating all the Wilson coefficients for various Lepton Flav...
Definition: LeptonFlavour.h:26
std::string name
The name of the model.
Definition: Model.h:285
An observable class for the -boson mass.
Definition: Mw.h:22
A class for the PMNS matrix elements.
Definition: PMNS.h:23
gslpp::matrix< gslpp::complex > getPMNS() const
A member for returning the PMNS matrix.
Definition: PMNS.h:42
A class for particles.
Definition: Particle.h:26
A class for parameters related to QCD, hadrons and quarks.
Definition: QCD.h:304
const double Mrun(const double mu, const double m, const quark q, const orders order=FULLNNLO) const
Computes a running quark mass from .
Definition: QCD.cpp:1353
quark
An enum type for quarks.
Definition: QCD.h:323
@ UP
Definition: QCD.h:324
@ BOTTOM
Definition: QCD.h:329
@ TOP
Definition: QCD.h:328
@ DOWN
Definition: QCD.h:325
@ STRANGE
Definition: QCD.h:327
@ CHARM
Definition: QCD.h:326
lepton
An enum type for leptons.
Definition: QCD.h:310
@ MU
Definition: QCD.h:314
@ ELECTRON
Definition: QCD.h:312
@ TAU
Definition: QCD.h:316
const Particle & getQuarks(const QCD::quark q) const
A get method to access a quark as an object of the type Particle.
Definition: QCD.h:536
Particle quarks[6]
The vector of all SM quarks.
Definition: QCD.h:1027
const double Als(const double mu, const orders order=FULLNLO, const bool Nf_thr=true) const
Definition: QCD.cpp:762
A model class for the Standard Model.
const double getIntegrand_sigmaWithISR_down183(double x) const
virtual const double LEP2sigmaCharm(const double s) const
double dAle5Mz
The five-flavour hadronic contribution to the electromagnetic coupling, , used as input for FlagMWinp...
const double getMuw() const
A get method to retrieve the matching scale around the weak scale.
EWSMThreeLoopEW * myThreeLoopEW
A pointer to an object of type EWSMThreeLoopEW.
const double computeBrHtomumu() const
The Br in the Standard Model.
const double getIntegrand_dsigmaBox_strange167(double x) const
double dAl5hMz
The five-flavour hadronic contribution to the electromagnetic coupling, . (Non-input parameter)
double Delta_EWQCD(const QCD::quark q) const
The non-factorizable EW-QCD corrections to the partial widths for , denoted as .
const double getIntegrand_sigmaWithISR_strange167(double x) const
double delSin2th_b
The theoretical uncertainty in , denoted as .
const double computeGammaHTotal() const
The Higgs total width in the Standard Model.
virtual const double GammaHtoZZstar() const
The in the Standard Model.
const gslpp::matrix< gslpp::complex > & getYe() const
A get method to retrieve the Yukawa matrix of the charged leptons, .
virtual const double getPhiBs() const
Half the relative phase of the $B_s$ mixing amplitude w.r.t. the Standard Model one.
virtual const double SigmaeeHvv(const double sqrt_s, const double Pe, const double Pp) const
The in the Standard Model.
double Vub
used as an input for FlagWolfenstein = FALSE
virtual const double BrHtocc() const
The Br in the Standard Model.
double taub() const
Top-mass corrections to the vertex, denoted by .
static std::string SMvars[NSMvars]
A string array containing the labels of the model parameters in StandardModel.
const double MRL2eeff(const Particle f, const double s) const
const double getIntegrand_sigmaWithISR_strange189(double x) const
const double getIntegrand_dsigmaBox_mu200(double x) const
double A
The CKM parameter in the Wolfenstein parameterization.
const double getIntegrand_sigmaWithISR_strange161(double x) const
virtual const double sin2thetaEff(const Particle f) const
The effective weak mixing angle for at the the -mass scale.
EWSMcache * getMyEWSMcache() const
A get method to retrieve the member pointer of type EWSMcache.
bool requireCKM
An internal flag to control whether the CKM matrix has to be recomputed.
virtual const double GammaHtotautau() const
The in the Standard Model.
const double getIntegrand_AFBnumeratorWithISR_tau130(double x) const
virtual const double GammaZ(const Particle f) const
The partial decay width, .
const double uovers2(const double cosmin, const double cosmax) const
const gslpp::matrix< gslpp::complex > & getYu() const
A get method to retrieve the Yukawa matrix of the up-type quarks, .
virtual void computeCKM()
The method to compute the CKM matrix.
const double computeBrHtoZZ() const
The Br in the Standard Model.
const double getDelMw() const
A get method to retrieve the theoretical uncertainty in , denoted as .
const double getIntegrand_AFBnumeratorWithISR_charm133(double x) const
const double getIntegrand_sigmaWithISR_charm192(double x) const
const gslpp::matrix< gslpp::complex > getUPMNS() const
A get method to retrieve the object of the PMNS matrix.
const double getIntegrand_AFBnumeratorWithISR_tau200(double x) const
double gamma
used as an input for FlagWolfenstein = FALSE
double Vud
used as an input for FlagWolfenstein = FALSE and FlagUseVud = TRUE
const double getIntegrand_dsigmaBox_up130(double x) const
const double getIntegrand_dsigmaBox_down207(double x) const
virtual const double getCCC1() const
A virtual implementation for the RealWeakEFTCC class.
virtual const double BrHtoZZstar() const
The Br in the Standard Model.
const double getIntegrand_dsigmaBox_tau136(double x) const
const bool IsFlagNoApproximateGammaZ() const
A method to retrieve the model flag NoApproximateGammaZ.
virtual const double eeffRelectron(const double pol_e, const double pol_p, const double s) const
const gslpp::matrix< gslpp::complex > & getYd() const
A get method to retrieve the Yukawa matrix of the down-type quarks, .
virtual bool PreUpdate()
The pre-update method for StandardModel.
const double getIntegrand_AFBnumeratorWithISR_mu192(double x) const
const double getIntegrand_sigmaWithISR_tau207(double x) const
const Particle & getLeptons(const QCD::lepton p) const
A get method to retrieve the member object of a lepton.
const double getIntegrand_AFBnumeratorWithISR_bottom192(double x) const
const double getIntegrand_dsigmaBox_down200(double x) const
const double tovers2(const double cosmin, const double cosmax) const
const double getDelR0c() const
A get method to retrieve the theoretical uncertainty in , denoted as .
const double computeBrHto4q() const
The Br in the Standard Model.
gslpp::complex AHZga_W(const double tau, const double lambda) const
W loop function entering in the calculation of the effective coupling.
const double getIntegrand_dsigmaBox_strange172(double x) const
const double sW2_MSbar_Approx() const
The (approximated formula for the) square of the sine of the weak mixing angle in the MSbar scheme,...
const double getIntegrand_dsigmaBox_mu207(double x) const
const double computeSigmattH(const double sqrt_s) const
The ttH production cross section in the Standard Model.
virtual const StandardModel & getTrueSM() const
virtual ~StandardModel()
The default destructor.
virtual const double LEP2sigmaHadron(const double s) const
const double getIntegrand_dsigmaBox_up202(double x) const
const double getIntegrand_sigmaWithISR_down192(double x) const
const double computeSigmaggH(const double sqrt_s) const
The ggH cross section in the Standard Model.
std::string FlagRhoZ
A string for the model flag RhoZ.
const double getIntegrand_AFBnumeratorWithISR_mu207(double x) const
virtual const double BrHtomumu() const
The Br in the Standard Model.
const double getIntegrand_dsigmaBox_strange130(double x) const
const double getIntegrand_AFBnumeratorWithISR_charm183(double x) const
virtual const double LEP2sigmaTau(const double s) const
virtual const double eeffsigmaCharm(const double pol_e, const double pol_p, const double s) const
double Mz
The mass of the boson in GeV.
EWSMThreeLoopQCD * myThreeLoopQCD
A pointer to an object of type EWSMThreeLoopQCD.
virtual const double TauLFU_gtaugmuPi() const
The computation of the LFU ratio .
double delGammaWlv
The theoretical uncertainty in , denoted as .
virtual const double GammaHtogaga() const
The in the Standard Model.
const double getIntegrand_AFBnumeratorWithISR_tau172(double x) const
const double MLL2eeff(const Particle f, const double s, const double t) const
const double getIntegrand_sigmaWithISR_tau205(double x) const
double m_q(const QCD::quark q, const double mu, const orders order=FULLNLO) const
const double getIntegrand_AFBnumeratorWithISR_bottom196(double x) const
const double getIntegrand_sigmaWithISR_strange136(double x) const
const double getIntegrand_sigmaWithISR_down207(double x) const
const double computeBrHtocc() const
The Br in the Standard Model.
const double getIntegrand_sigmaWithISR_tau196(double x) const
const double getIntegrand_dsigmaBox_down136(double x) const
const double getIntegrand_dsigmaBox_strange189(double x) const
virtual const double getCCC3() const
A virtual implementation for the RealWeakEFTCC class.
const double getIntegrand_sigmaWithISR_bottom200(double x) const
const double MwFromMwbar(const double Mwbar) const
A method to convert the -boson mass in the complex-pole/fixed-width scheme to that in the experimenta...
const double getIntegrand_sigmaWithISR_up183(double x) const
virtual const double GammaHtobb() const
The in the Standard Model.
const double getIntegrand_dsigmaBox_mu136(double x) const
const double getIntegrand_sigmaWithISR_strange183(double x) const
EWSMThreeLoopEW2QCD * getMyThreeLoopEW2QCD() const
virtual const double AFB(const Particle f) const
const double getIntegrand_dsigmaBox_down183(double x) const
double RAq(const QCD::quark q) const
The radiator factor associated with the final-state QED and QCD corrections to the the axial-vector-c...
const double getIntegrand_dsigmaBox_charm130(double x) const
const double getMz() const
A get method to access the mass of the boson .
const double getIntegrand_AFBnumeratorWithISR_tau192(double x) const
const double getIntegrand_sigmaWithISR_strange200(double x) const
void setYu(const gslpp::matrix< gslpp::complex > &Yu)
A set method to set the Yukawa matrix of the up-type quarks, .
const double getIntegrand_AFBnumeratorWithISR_charm205(double x) const
void setFlagCacheInStandardModel(bool FlagCacheInStandardModel)
A set method to change the model flag CacheInStandardModel of StandardModel.
const double getIntegrand_sigmaWithISR_tau161(double x) const
const double getDelSin2th_q() const
A get method to retrieve the theoretical uncertainty in , denoted as .
const double getMw() const
A get method to access the input value of the mass of the boson .
const double getIntegrand_sigmaWithISR_up167(double x) const
double GammaW_cache
A cache of the value of .
double delMw
The theoretical uncertainty in , denoted as , in GeV.
const double getIntegrand_dsigmaBox_tau133(double x) const
virtual const double LEP2AFBtau(const double s) const
const double getIntegrand_AFBnumeratorWithISR_tau196(double x) const
bool flag_order[orders_EW_size]
An array of internal flags controlling the inclusions of higher-order corrections.
virtual const gslpp::complex rhoZ_f(const Particle f) const
The effective leptonic neutral-current coupling in the SM.
const double getIntegrand_AFBnumeratorWithISR_charm200(double x) const
virtual const double gLnuN2() const
The effective neutrino nucleon LH coupling: gLnuN2.
const double getIntegrand_sigmaWithISR_strange207(double x) const
const double getIntegrand_AFBnumeratorWithISR_mu196(double x) const
const double getIntegrand_dsigmaBox_mu189(double x) const
const double Ale(double mu, orders order, bool Nf_thr=true) const
The running electromagnetic coupling in the scheme.
const double getIntegrand_dsigmaBox_down161(double x) const
virtual const double Gamma_tau_l_nunu(const Particle l) const
The computation of the leptonic tau decays.
StandardModel()
The default constructor.
const double getIntegrand_dsigmaBox_strange207(double x) const
const double getIntegrand_sigmaWithISR_down202(double x) const
gsl_function f_GSL
const double Integrand_dsigmaBox_q(double cosTheta, const QCD::quark q_flavor, const double s) const
const double computeSigmaZWF(const double sqrt_s) const
The Z W interference fusion contribution to higgs-production cross section in the Standard Model.
EWSMTwoLoopQCD * getMyTwoLoopQCD() const
virtual const double LEP2dsigmadcosBinTau(const double s, const double cos, const double cosmin, const double cosmax) const
const double getIntegrand_dsigmaBox_bottom189(double x) const
EWSMTwoLoopEW * myTwoLoopEW
A pointer to an object of type EWSMTwoLoopEW.
const double computeSigmaVBF(const double sqrt_s) const
The VBF cross section in the Standard Model.
const double getIntegrand_AFBnumeratorWithISR_mu136(double x) const
const std::string getFlagMw() const
A method to retrieve the model flag Mw.
const double getIntegrand_AFBnumeratorWithISR_mu130(double x) const
virtual bool CheckParameters(const std::map< std::string, double > &DPars)
A method to check if all the mandatory parameters for StandardModel have been provided in model initi...
virtual const double getCBd() const
The ratio of the absolute value of the $B_d$ mixing amplitude over the Standard Model value.
virtual bool CheckFlags() const
A method to check the sanity of the set of model flags.
const double getIntegrand_dsigmaBox_charm161(double x) const
const double getIntegrand_dsigmaBox_up192(double x) const
void ComputeDeltaRho(const double Mw_i, double DeltaRho[orders_EW_size]) const
A method to collect computed via subclasses.
const double getIntegrand_sigmaWithISR_strange202(double x) const
double RVq(const QCD::quark q) const
The radiator factor associated with the final-state QED and QCD corrections to the the vector-current...
virtual const double LEP2Rbottom(const double s) const
const double getIntegrand_sigmaWithISR_tau183(double x) const
const double getIntegrand_dsigmaBox_mu133(double x) const
bool FlagFixMuwMut
A boolean for the model flag FixMuwMut.
const double AlsEByOrder(double mu, orders order, bool Nf_thr) const
const CKM & getCKM() const
A get method to retrieve the member object of type CKM.
const double getAlsMz() const
A get method to access the value of .
const double getIntegrand_sigmaWithISR_tau130(double x) const
const double getIntegrand_AFBnumeratorWithISR_bottom172(double x) const
virtual const double getPhiBd() const
Half the relative phase of the $B_d$ mixing amplitude w.r.t. the Standard Model one.
virtual const double RWlilj(const Particle li, const Particle lj) const
The lepton universality ratio .
gslpp::complex AH_f(const double tau) const
Fermionic loop function entering in the calculation of the effective and couplings.
const double computeSigmaWH(const double sqrt_s) const
The WH production cross section in the Standard Model.
virtual const double Gamma_inv() const
The invisible partial decay width of the boson, .
const double getIntegrand_dsigmaBox_bottom207(double x) const
const double getIntegrand_AFBnumeratorWithISR_mu172(double x) const
virtual const double Qwp() const
The computation of the proton weak charge: Qwp.
const double computeBrHtotautau() const
The Br in the Standard Model.
const double getIntegrand_dsigmaBox_tau183(double x) const
virtual const double getMHl() const
A get method to retrieve the Higgs mass .
const double getIntegrand_AFBnumeratorWithISR_mu202(double x) const
const double getIntegrand_AFBnumeratorWithISR_charm172(double x) const
const double computeSigmabbH(const double sqrt_s) const
The bbH production cross section in the Standard Model.
const double getIntegrand_dsigmaBox_down202(double x) const
virtual const double R_inv() const
The ratio of the invisible and leptonic (electron) decay widths of the boson, .
virtual bool Init(const std::map< std::string, double > &DPars)
A method to initialize the model parameters.
CKM myCKM
An object of type CKM.
virtual const double LEP2AFBe(const double s) const
double SMresult_cache
virtual const double eeffAFBmu(const double pol_e, const double pol_p, const double s) const
const double getIntegrand_sigmaWithISR_charm183(double x) const
const double getIntegrand_dsigmaBox_bottom202(double x) const
EWSMOneLoopEW * getMyOneLoopEW() const
A get method to retrieve the member pointer of type EWSMOneLoopEW,.
const double getIntegrand_dsigmaBox_bottom200(double x) const
virtual const double eeffAFBtau(const double pol_e, const double pol_p, const double s) const
bool checkEWPOscheme(const std::string scheme) const
A method to check if a given scheme name in string form is valid.
virtual const double eeffsigmaMu(const double pol_e, const double pol_p, const double s) const
const double getIntegrand_dsigmaBox_bottom205(double x) const
bool useDeltaAlpha_cache
const double getIntegrand_dsigmaBox_tau196(double x) const
const double computeBrHto4f() const
The Br in the Standard Model.
virtual const double getCCC4() const
A virtual implementation for the RealWeakEFTCC class.
const double getIntegrand_sigmaWithISR_tau200(double x) const
const double getIntegrand_sigmaWithISR_bottom207(double x) const
const double getIntegrand_dsigmaBox_bottom183(double x) const
virtual bool setFlagStr(const std::string name, const std::string value)
A method to set a flag of StandardModel.
bool requireYn
An internal flag to control whether the neutrino Yukawa matrix has to be recomputed.
const double getIntegrand_dsigmaBox_mu161(double x) const
virtual const double eeffsigma(const Particle f, const double pol_e, const double pol_p, const double s, const double cosmin, const double cosmax) const
virtual const double BrHtotautau() const
The Br in the Standard Model.
virtual const double Gamma_muon() const
The computation of the muon decay.
const double getIntegrand_dsigmaBox_up200(double x) const
const double intMLR2eeeets2(const double s, const double t0, const double t1) const
const double getIntegrand_dsigmaBox_tau192(double x) const
const double getIntegrand_AFBnumeratorWithISR_bottom207(double x) const
const double getIntegrand_AFBnumeratorWithISR_tau136(double x) const
const double DeltaAlphaTop(const double s) const
Top-quark contribution to the electromagnetic coupling , denoted as .
const double getIntegrand_sigmaWithISR_charm205(double x) const
const double getIntegrand_dsigmaBox_tau189(double x) const
virtual const double BrHtogaga() const
The Br in the Standard Model.
const double getIntegrand_AFBnumeratorWithISR_charm207(double x) const
gsl_integration_workspace * w_GSL1
virtual const double SigmaeeZH(const double sqrt_s, const double Pe, const double Pp) const
The in the Standard Model.
virtual const gslpp::complex kappaZ_f(const Particle f) const
The effective leptonic neutral-current coupling in the SM.
gslpp::matrix< gslpp::complex > Yn
The Yukawa matrix of the neutrinos.
const double getIntegrand_AFBnumeratorWithISR_charm202(double x) const
const double computeBrHtobb() const
The Br in the Standard Model.
EWSMTwoFermionsLEP2 * myTwoFermionsLEP2
A pointer to an object of type EWSMTwoFermionsLEP2.
const double computeGammaHgaga_tt() const
The top loop contribution to in the Standard Model.
const std::string getFlagRhoZ() const
A method to retrieve the model flag RhoZ.
Matching< StandardModelMatching, StandardModel > SMM
An object of type Matching.
virtual const double Gamma_had() const
The hadronic decay width of the boson, .
const double getIntegrand_sigmaWithISR_charm207(double x) const
const double ale_OS(const double mu, orders order=FULLNLO) const
The running electromagnetic coupling in the on-shell scheme.
bool setFlagSigmaForAFB(const bool flagSigmaForAFB_i)
virtual void computeYukawas()
The method to compute the Yukawas matrix.
virtual const double rho_GammaW(const Particle fi, const Particle fj) const
EW radiative corrections to the width of , denoted as .
const double AleWithInit(double mu, double alsi, double mu_i, orders order) const
const double getIntegrand_sigmaWithISR_bottom136(double x) const
virtual const double eeffRmuon(const double pol_e, const double pol_p, const double s) const
const double getIntegrand_sigmaWithISR_charm202(double x) const
gslpp::matrix< gslpp::complex > Yu
The Yukawa matrix of the up-type quarks.
virtual const double LEP2AFBcharm(const double s) const
std::string FlagMw
A string for the model flag Mw.
virtual const double eeffAFBcharm(const double pol_e, const double pol_p, const double s) const
const double Integrand_sigmaWithISR_l(double x, const QCD::lepton l_flavor, const double s) const
virtual const double eeffRstrange(const double pol_e, const double pol_p, const double s) const
const double getIntegrand_AFBnumeratorWithISR_tau183(double x) const
virtual const double GammaHtogg() const
The in the Standard Model.
double delsigma0H
The theoretical uncertainty in , denoted as in nb.
virtual const double BrHtobb() const
The Br in the Standard Model.
virtual const gslpp::complex gA_f(const Particle f) const
The effective leptonic neutral-current axial-vector coupling in the SM.
const double getIntegrand_sigmaWithISR_bottom130(double x) const
const int getIterationNo() const
const double getIntegrand_sigmaWithISR_charm161(double x) const
const double getIntegrand_sigmaWithISR_mu207(double x) const
virtual const double GammaHtoZga() const
The in the Standard Model.
const double getDelR0l() const
A get method to retrieve the theoretical uncertainty in , denoted as .
const double getIntegrand_sigmaWithISR_bottom192(double x) const
void setSMSuccess(bool success) const
A set method to change the success status of the Standard Model update and matching.
const double computeBrHto4l3() const
The Br in the Standard Model.
virtual const double alrmoller(const double q2, const double y) const
The computation of the parity violating asymmetry in Moller scattering.
double rhob
The CKM parameter in the Wolfenstein parameterization.
virtual const double eeffsigmaEtsub(const double pol_e, const double pol_p, const double s) const
Particle leptons[6]
An array of Particle objects for the leptons.
const double Integrand_AFBnumeratorWithISR_q(double x, const QCD::quark q_flavor, const double s) const
const double getIntegrand_sigmaWithISR_charm172(double x) const
const double getIntegrand_sigmaWithISR_charm200(double x) const
const std::string getFlagKappaZ() const
A method to retrieve the model flag KappaZ.
const double getIntegrand_dsigmaBox_up167(double x) const
const double getIntegrand_dsigmaBox_mu183(double x) const
const double computeSigmaggH_tt(const double sqrt_s) const
The square of the top-quark contribution to the ggH cross section in the Standard Model.
const double eeffsigmaEbin(const double pol_e, const double pol_p, const double s, const double cosmin, const double cosmax) const
double delSin2th_l
The theoretical uncertainty in , denoted as .
const double getIntegrand_dsigmaBox_charm167(double x) const
virtual const double RWc() const
The ratio .
const double getIntegrand_sigmaWithISR_mu172(double x) const
const double getIntegrand_sigmaWithISR_mu136(double x) const
const double sigma_NoISR_q(const QCD::quark q_flavor, const double s) const
bool FlagMWinput
A boolean for the model flag MWinput.
const double computeSigmaggH_bb(const double sqrt_s) const
The square of the bottom-quark contribution to the ggH cross section in the Standard Model.
EWSMThreeLoopQCD * getMyThreeLoopQCD() const
virtual const double Ruc() const
const double computeGammaHgg_tb() const
The top-bottom interference contribution to in the Standard Model.
const double getDelSin2th_b() const
A get method to retrieve the theoretical uncertainty in , denoted as .
const double getIntegrand_dsigmaBox_strange200(double x) const
const double getIntegrand_dsigmaBox_mu172(double x) const
const bool isSMSuccess() const
A get method to retrieve the success status of the Standard Model update and matching.
const double getIntegrand_sigmaWithISR_bottom196(double x) const
virtual const double sigma0_had() const
The hadronic cross section for at the -pole, .
virtual const double GammaW() const
The total width of the boson, .
const double intMLRtilde2eeeest2(const double s, const double t0, const double t1) const
virtual const double getCDMK() const
The ratio of the real part of the $K$ mixing amplitude over the Standard Model value.
const double getIntegrand_dsigmaBox_strange196(double x) const
const double s02() const
The square of the sine of the weak mixing angle defined without weak radiative corrections.
Flavour SMFlavour
An object of type Flavour.
gslpp::complex g_triangle(const double tau) const
Loop function entering in the calculation of the effective coupling.
const double getIntegrand_dsigmaBox_up207(double x) const
const double getIntegrand_sigmaWithISR_charm133(double x) const
bool FlagWithoutNonUniversalVC
A boolean for the model flag WithoutNonUniversalVC.
const double computeBrHto4l2() const
The Br in the Standard Model.
const double getIntegrand_sigmaWithISR_down200(double x) const
const double computeBrHtollvv3() const
The Br in the Standard Model.
bool setFlagSigmaForR(const bool flagSigmaForR_i)
bool FlagSMAux
A boolean for the model flag SMAux.
const double computeBrHtollvv2() const
The Br in the Standard Model.
const double getIntegrand_sigmaWithISR_charm130(double x) const
const double getIntegrand_sigmaWithISR_down196(double x) const
const double getIntegrand_sigmaWithISR_bottom172(double x) const
const double computeBrHtogg() const
The Br in the Standard Model.
gslpp::matrix< gslpp::complex > Yd
The Yukawa matrix of the down-type quarks.
const gslpp::matrix< gslpp::complex > getVCKM() const
A get method to retrieve the CKM matrix.
const double getIntegrand_sigmaWithISR_up133(double x) const
const double computeSigmatHq(const double sqrt_s) const
The tHq production cross section in the Standard Model.
virtual const double eeffAFBbottom(const double pol_e, const double pol_p, const double s) const
virtual const double LEP2AFBmu(const double s) const
virtual const double BrW(const Particle fi, const Particle fj) const
The branching ratio of the boson decaying into a SM fermion pair, .
const double getIntegrand_dsigmaBox_down189(double x) const
const double getDelR0b() const
A get method to retrieve the theoretical uncertainty in , denoted as .
const double getIntegrand_AFBnumeratorWithISR_bottom200(double x) const
const double getIntegrand_sigmaWithISR_charm167(double x) const
const double Integrand_dsigmaBox_l(double cosTheta, const QCD::lepton l_flavor, const double s) const
EWSMOneLoopEW * myOneLoopEW
A pointer to an object of type EWSMOneLoopEW.
double delR0c
The theoretical uncertainty in , denoted as .
const double computeGammaHZga_WW() const
The loop contribution to in the Standard Model. Currently it returns the value of tab 41 in ref....
virtual const double epsilon2() const
The SM contribution to the epsilon parameter .
const double getIntegrand_dsigmaBox_charm205(double x) const
const double getIntegrand_sigmaWithISR_up172(double x) const
virtual const double LEP2dsigmadcosMu(const double s, const double cos) const
const double getDelSin2th_l() const
A get method to retrieve the theoretical uncertainty in , denoted as .
const double getIntegrand_dsigmaBox_up161(double x) const
virtual const double Qwemoller(const double q2, const double y) const
The computation of the electron's weak charge.
const double getIntegrand_dsigmaBox_tau202(double x) const
virtual const double DeltaR() const
The SM prediction for derived from that for the boson mass.
std::string FlagKappaZ
A string for the model flag KappaZ.
const double getIntegrand_dsigmaBox_charm183(double x) const
const double getIntegrand_dsigmaBox_bottom196(double x) const
virtual const double LEP2dsigmadcosBinE(const double s, const double cos, const double cosmin, const double cosmax) const
const double getIntegrand_sigmaWithISR_strange205(double x) const
const double getIntegrand_dsigmaBox_bottom172(double x) const
virtual const double Gamma_Z() const
The total decay width of the boson, .
LeptonFlavour * getMyLeptonFlavour() const
const double getIntegrand_dsigmaBox_down133(double x) const
const double getIntegrand_dsigmaBox_mu192(double x) const
const double getIntegrand_dsigmaBox_down167(double x) const
const double getIntegrand_dsigmaBox_mu196(double x) const
virtual const double gRnuN2() const
The effective neutrino nucleon RH coupling: gRnuN2.
const double getIntegrand_sigmaWithISR_strange133(double x) const
const bool IsFlagWithoutNonUniversalVC() const
A method to retrieve the model flag WithoutNonUniversalVC.
const double getIntegrand_sigmaWithISR_tau189(double x) const
const double getIntegrand_AFBnumeratorWithISR_bottom205(double x) const
double GF
The Fermi constant in .
const double getIntegrand_sigmaWithISR_down130(double x) const
const double getIntegrand_sigmaWithISR_strange192(double x) const
virtual const double epsilonb() const
The SM contribution to the epsilon parameter .
EWSMApproximateFormulae * myApproximateFormulae
A pointer to an object of type EWSMApproximateFormulae.
virtual const gslpp::complex deltaRhoZ_f(const Particle f) const
Flavour non-universal vertex corrections to , denoted by .
const double getIntegrand_dsigmaBox_bottom130(double x) const
virtual const double eeffAFBetsub(const double pol_e, const double pol_p, const double s) const
virtual const double GammaHtocc() const
The in the Standard Model.
const double getIntegrand_sigmaWithISR_strange130(double x) const
const double getIntegrand_AFBnumeratorWithISR_bottom189(double x) const
const double getIntegrand_dsigmaBox_tau172(double x) const
const double intMRR2eeeeus2(const double s, const double t0, const double t1) const
double resumKappaZ(const double DeltaRho[orders_EW_size], const double deltaKappa_rem[orders_EW_size], const double DeltaRbar_rem, const bool bool_Zbb) const
A method to compute the real part of the effetvive coupling from , and .
virtual const double TauLFU_gtaugmuK() const
The computation of the LFU ratio .
const double getIntegrand_dsigmaBox_strange161(double x) const
virtual const double eeffRbottom(const double pol_e, const double pol_p, const double s) const
virtual const double LEP2sigmaE(const double s) const
const double getIntegrand_sigmaWithISR_down136(double x) const
const double Beta_s(int nm, unsigned int nf) const
QCD beta function coefficients including QED corrections - eq. (36) hep-ph/0512066.
const double AFB_NoISR_q(const QCD::quark q_flavor, const double s) const
virtual const double BrHtoWWstar() const
The Br in the Standard Model.
const double getIntegrand_sigmaWithISR_bottom205(double x) const
const double getIntegrand_dsigmaBox_tau200(double x) const
const double computeGammaHgaga_WW() const
The loop contribution to in the Standard Model.
virtual const double Mw() const
The SM prediction for the -boson mass in the on-shell scheme, .
virtual bool setFlag(const std::string name, const bool value)
A method to set a flag of StandardModel.
const double getIntegrand_sigmaWithISR_up192(double x) const
virtual const double A_f(const Particle f) const
The left-right asymmetry in at the -pole, .
const double MRR2eeff(const Particle f, const double s, const double t) const
double Mw_cache
A cache of the value of .
double DeltaAlphaLepton_cache
A cache of the value of .
const double getIntegrand_AFBnumeratorWithISR_mu205(double x) const
const double getIntegrand_dsigmaBox_bottom133(double x) const
const double computeGammaHZga_tt() const
The top loop contribution to in the Standard Model.
const double getIntegrand_dsigmaBox_strange205(double x) const
virtual const gslpp::complex deltaKappaZ_f(const Particle f) const
Flavour non-universal vertex corrections to , denoted by .
const double getIntegrand_AFBnumeratorWithISR_tau202(double x) const
const double getIntegrand_dsigmaBox_bottom167(double x) const
virtual const double amuon() const
The computation of the anomalous magnetic moment of the muon .
const Flavour & getFlavour() const
virtual const double BrHtoZga() const
The Br in the Standard Model.
double delGammaWqq
The theoretical uncertainty in , denoted as .
gslpp::complex I_triangle_1(const double tau, const double lambda) const
Loop function entering in the calculation of the effective coupling.
virtual const double epsilon1() const
The SM contribution to the epsilon parameter .
EWSMTwoLoopEW * getMyTwoLoopEW() const
const double computeBrHtoZga() const
The Br in the Standard Model.
double delSin2th_q
The theoretical uncertainty in , denoted as .
bool FlagUseVud
A boolean for the model flag UseVud.
const double getIntegrand_sigmaWithISR_bottom133(double x) const
const double getIntegrand_sigmaWithISR_mu161(double x) const
const double computeSigmaZH(const double sqrt_s) const
The ZH production cross section in the Standard Model.
const double computeSigmaWF(const double sqrt_s) const
The W fusion contribution to higgs-production cross section in the Standard Model.
double Mzbar() const
The -boson mass in the complex-pole/fixed-width scheme.
static const double Mw_error
The target accuracy of the iterative calculation of the -boson mass in units of GeV.
static const int NSMvars
The number of the model parameters in StandardModel.
const double getIntegrand_dsigmaBox_mu202(double x) const
const double getIntegrand_sigmaWithISR_up200(double x) const
const double sigma_NoISR_l(const QCD::lepton l_flavor, const double s) const
const double Integrand_sigmaWithISR_q(double x, const QCD::quark q_flavor, const double s) const
const double getIntegrand_dsigmaBox_up183(double x) const
virtual const double LEP2Rcharm(const double s) const
double SMparamsForEWPO_cache[NumSMParamsForEWPO]
const double computeBrHtogaga() const
The Br in the Standard Model.
bool FlagNoApproximateGammaZ
A boolean for the model flag NoApproximateGammaZ.
const double getIntegrand_sigmaWithISR_bottom202(double x) const
const double getIntegrand_sigmaWithISR_mu130(double x) const
const double getIntegrand_dsigmaBox_strange183(double x) const
gslpp::complex kappaZ_f_cache[12]
A cache of the value of .
virtual const double getCCC2() const
A virtual implementation for the RealWeakEFTCC class.
const double getIntegrand_dsigmaBox_up133(double x) const
virtual const double GammaHtoss() const
The in the Standard Model.
double lambda
The CKM parameter in the Wolfenstein parameterization.
const double getIntegrand_sigmaWithISR_strange172(double x) const
virtual StandardModelMatching & getMatching() const
A get method to access the member reference of type StandardModelMatching.
virtual const double Qwn() const
The computation of the neutron weak charge: Qwn.
const double getGF() const
A get method to retrieve the Fermi constant .
double als_cache[11][CacheSize]
Cache for .
const double getIntegrand_AFBnumeratorWithISR_mu189(double x) const
const double computeGammaHgg_tt() const
The top loop contribution to in the Standard Model.
void ComputeDeltaR_rem(const double Mw_i, double DeltaR_rem[orders_EW_size]) const
A method to collect computed via subclasses.
const double getIntegrand_AFBnumeratorWithISR_bottom202(double x) const
virtual const double eeffAFBe(const double pol_e, const double pol_p, const double s) const
const double getIntegrand_dsigmaBox_strange202(double x) const
void setRequireCKM(bool requireCKM)
A set method to change the value of requireCKM.
const double getIntegrand_sigmaWithISR_mu202(double x) const
const double Integrand_AFBnumeratorWithISR_l(double x, const QCD::lepton l_flavor, const double s) const
const double getIntegrand_sigmaWithISR_up189(double x) const
const double getIntegrand_dsigmaBox_down172(double x) const
const double getIntegrand_dsigmaBox_tau161(double x) const
const double getIntegrand_sigmaWithISR_tau192(double x) const
virtual const double RZlilj(const Particle li, const Particle lj) const
The lepton universality ratio .
const double getIntegrand_dsigmaBox_strange136(double x) const
const double computeSigmaZF(const double sqrt_s) const
The Z fusion contribution to higgs-production cross section in the Standard Model.
const double Als(const double mu, const orders order, const bool Nf_thr, const bool qed_flag) const
The running QCD coupling in the scheme including QED corrections.
const double getIntegrand_AFBnumeratorWithISR_charm189(double x) const
const double getIntegrand_sigmaWithISR_up136(double x) const
const double Beta_e(int nm, unsigned int nf) const
QED beta function coefficients - eq. (36) hep-ph/0512066.
const double getIntegrand_dsigmaBox_bottom192(double x) const
double Vcb
used as an input for FlagWolfenstein = FALSE
const double c02() const
The square of the cosine of the weak mixing angle defined without weak radiative corrections.
virtual const double LEP2sigmaMu(const double s) const
const double getIntegrand_dsigmaBox_mu205(double x) const
const double getIntegrand_sigmaWithISR_strange196(double x) const
const double getIntegrand_sigmaWithISR_tau136(double x) const
virtual const double R0_f(const Particle f) const
The ratio .
const double getIntegrand_dsigmaBox_up205(double x) const
const double Mw_tree() const
The tree-level mass of the boson, .
const double getDelSigma0H() const
A get method to retrieve the theoretical uncertainty in , denoted as .
const double AlsE(double mu, orders order, bool Nf_thr) const
const double getIntegrand_AFBnumeratorWithISR_mu200(double x) const
virtual const double eeffRcharm(const double pol_e, const double pol_p, const double s) const
const double getIntegrand_sigmaWithISR_up202(double x) const
bool useDeltaAlphaLepton_cache
virtual const double TauLFU_gtauge() const
The computation of the LFU ratio .
const double getIntegrand_dsigmaBox_mu167(double x) const
double etab
The CKM parameter in the Wolfenstein parameterization.
const double getIntegrand_sigmaWithISR_down167(double x) const
virtual const double eeffsigmaBottom(const double pol_e, const double pol_p, const double s) const
const double getIntegrand_sigmaWithISR_mu205(double x) const
virtual const double getCCC5() const
A virtual implementation for the RealWeakEFTCC class.
virtual const double eeffsigmaE(const double pol_e, const double pol_p, const double s) const
virtual const double SigmaeeHee(const double sqrt_s, const double Pe, const double Pp) const
The in the Standard Model.
const double getIntegrand_sigmaWithISR_down133(double x) const
virtual const double BrHtogg() const
The Br in the Standard Model.
EWSMcache * myEWSMcache
A pointer to an object of type EWSMcache.
gslpp::complex rhoZ_f_cache[12]
A cache of the value of .
const double getIntegrand_AFBnumeratorWithISR_tau207(double x) const
const double AFB_NoISR_l(const QCD::lepton l_flavor, const double s) const
virtual const double eeffsigmaStrange(const double pol_e, const double pol_p, const double s) const
virtual const double epsilon3() const
The SM contribution to the epsilon parameter .
LeptonFlavour * myLeptonFlavour
A pointer to an object of the type LeptonFlavour.
EWSMApproximateFormulae * getMyApproximateFormulae() const
A get method to retrieve the member pointer of type EWSMApproximateFormulae.
double ale_cache[10][CacheSize]
Cache for .
const double getIntegrand_sigmaWithISR_charm136(double x) const
const double getIntegrand_dsigmaBox_down130(double x) const
virtual const double Dalpha5hMz() const
The 5-quark contribution to the running of the em constant to the pole. .
double resumRhoZ(const double DeltaRho[orders_EW_size], const double deltaRho_rem[orders_EW_size], const double DeltaRbar_rem, const bool bool_Zbb) const
A method to compute the real part of the effective coupling from , and .
gslpp::complex I_triangle_2(const double tau, const double lambda) const
Loop function entering in the calculation of the effective coupling.
virtual const double TauLFU_gtaugmu() const
The computation of the LFU ratio .
const double getIntegrand_dsigmaBox_strange192(double x) const
void setCKM(const CKM &CKMMatrix)
A set method to change the CKM matrix.
virtual const double N_nu() const
The number of neutrinos obtained indirectly from the measurements at the Z pole, .
const double getIntegrand_sigmaWithISR_down189(double x) const
virtual const double TauLFU_gmuge() const
The computation of the LFU ratio .
virtual const double ThetaLnuN() const
The effective neutrino nucleon LH parameter: ThetaLnuN.
const double getIntegrand_dsigmaBox_charm196(double x) const
bool checkSMparamsForEWPO()
A method to check whether the parameters relevant to the EWPO are updated.
const double getIntegrand_dsigmaBox_bottom161(double x) const
virtual const double BrHtoss() const
The Br in the Standard Model.
const double getIntegrand_sigmaWithISR_up207(double x) const
const double Alstilde5(const double mu) const
The value of at any scale with the number of flavours and full EW corrections.
const double computeSigmaggH_tb(const double sqrt_s) const
The top-bottom interference contribution to the ggH cross section in the Standard Model.
const double getIntegrand_sigmaWithISR_charm196(double x) const
const double getIntegrand_AFBnumeratorWithISR_tau205(double x) const
const double computeBrHtoss() const
The Br in the Standard Model.
const double getIntegrand_sigmaWithISR_mu189(double x) const
const double getIntegrand_sigmaWithISR_mu196(double x) const
const double getIntegrand_AFBnumeratorWithISR_bottom167(double x) const
const double getIntegrand_sigmaWithISR_tau202(double x) const
double Mw_inp
The mass of the boson in GeV used as input for FlagMWinput = TRUE.
const double getIntegrand_dsigmaBox_charm172(double x) const
orders_EW
An enumerated type representing perturbative orders of radiative corrections to EW precision observab...
@ EW1
One-loop of .
@ EW2QCD1
Three-loop of .
@ EW2
Two-loop of .
@ orders_EW_size
The size of this enum.
@ EW3
Three-loop of .
@ EW1QCD1
Two-loop of .
@ EW1QCD2
Three-loop of .
const double getIntegrand_dsigmaBox_charm136(double x) const
double mHl
The Higgs mass in GeV.
const double getIntegrand_sigmaWithISR_up205(double x) const
const double getIntegrand_AFBnumeratorWithISR_charm196(double x) const
gslpp::complex AH_W(const double tau) const
W loop function entering in the calculation of the effective coupling.
double SchemeToDouble(const std::string scheme) const
A method to convert a given scheme name in string form into a floating-point number with double preci...
const double getIntegrand_sigmaWithISR_bottom167(double x) const
const double getIntegrand_dsigmaBox_charm133(double x) const
const double getIntegrand_dsigmaBox_charm192(double x) const
const double getIntegrand_dsigmaBox_tau207(double x) const
double ale
The fine-structure constant .
bool flagLEP2[NUMofLEP2RCs]
const double DeltaAlphaL5q() const
The sum of the leptonic and the five-flavour hadronic corrections to the electromagnetic coupling at...
const double getIntegrand_sigmaWithISR_bottom183(double x) const
double AlsMz
The strong coupling constant at the Z-boson mass, .
double delGammaZ
The theoretical uncertainty in , denoted as , in GeV.
const double getIntegrand_AFBnumeratorWithISR_mu161(double x) const
gslpp::matrix< gslpp::complex > Ye
The Yukawa matrix of the charged leptons.
const double getIntegrand_dsigmaBox_down196(double x) const
const double computeBrHtoevmuv() const
The Br in the Standard Model.
virtual const double eeffsigmaTau(const double pol_e, const double pol_p, const double s) const
virtual bool PostUpdate()
The post-update method for StandardModel.
const double getIntegrand_dsigmaBox_charm189(double x) const
const double sW2_ND() const
The square of the sine of the weak mixing angle in the MSbar-ND scheme (w/o decoupling $\alpha\ln(m_t...
double muw
A matching scale around the weak scale in GeV.
double RVh() const
The singlet vector corrections to the hadronic -boson width, denoted as .
bool useRhoZ_f_cache[12]
virtual const double alphaMz() const
The electromagnetic coupling at the -mass scale, .
virtual const double eeffsigmaHadron(const double pol_e, const double pol_p, const double s) const
double delR0b
The theoretical uncertainty in , denoted as .
const double getIntegrand_AFBnumeratorWithISR_bottom133(double x) const
const double DeltaAlpha() const
The total corrections to the electromagnetic coupling at the -mass scale, denoted as .
const double getIntegrand_sigmaWithISR_up161(double x) const
virtual void setParameter(const std::string name, const double &value)
A method to set the value of a parameter of StandardModel.
gslpp::complex f_triangle(const double tau) const
Loop function entering in the calculation of the effective and couplings.
const double getIntegrand_dsigmaBox_mu130(double x) const
const double getIntegrand_sigmaWithISR_mu200(double x) const
void setYd(const gslpp::matrix< gslpp::complex > &Yd)
A set method to set the Yukawa matrix of the down-type quarks, .
const double getIntegrand_dsigmaBox_charm207(double x) const
static const int CacheSize
Defines the depth of the cache.
const double getDAle5Mz() const
A get method to retrieve the five-flavour hadronic contribution to the electromagnetic coupling,...
const double getIntegrand_AFBnumeratorWithISR_tau189(double x) const
static const int NumSMParamsForEWPO
The number of the SM parameters that are relevant to the EW precision observables.
virtual const double eeffAFBstrange(const double pol_e, const double pol_p, const double s) const
const double getIntegrand_sigmaWithISR_tau172(double x) const
double Vus
used as an input for FlagWolfenstein = FALSE
virtual const double GammaHtoWWstar() const
The in the Standard Model.
const double getIntegrand_AFBnumeratorWithISR_charm167(double x) const
const double computeGammaHgaga_tW() const
The mixed loop contribution to in the Standard Model.
const double computeGammaHgg_bb() const
The bottom loop contribution to in the Standard Model.
const double getIntegrand_sigmaWithISR_mu183(double x) const
EWSMTwoLoopQCD * myTwoLoopQCD
A pointer to an object of type EWSMTwoLoopQCD.
virtual const double LEP2dsigmadcosTau(const double s, const double cos) const
virtual const double GammaHtomumu() const
The in the Standard Model.
const double getIntegrand_AFBnumeratorWithISR_tau161(double x) const
static const double GeVminus2_to_nb
void setFlagNoApproximateGammaZ(bool FlagNoApproximateGammaZ)
bool FlagWolfenstein
A boolean for the model flag Wolfenstein.
const double getIntegrand_sigmaWithISR_down172(double x) const
virtual const double LEP2sigmaBottom(const double s) const
virtual const double getmq(const QCD::quark q, const double mu) const
The MSbar running quark mass computed at NLO.
const double getIntegrand_dsigmaBox_tau130(double x) const
virtual const double DeltaRbar() const
The SM prediction for derived from that for the -boson mass.
const double computeBrHto4v() const
The Br in the Standard Model.
const double getIntegrand_dsigmaBox_strange133(double x) const
const double getIntegrand_dsigmaBox_tau205(double x) const
virtual const double getCBs() const
The ratio of the absolute value of the $B_s$ mixing amplitude over the Standard Model value.
bool requireYe
An internal flag to control whether the charged-lepton Yukawa matrix has to be recomputed.
EWSMThreeLoopEW2QCD * myThreeLoopEW2QCD
A pointer to an object of type EWSMThreeLoopEW2QCD.
const double getDelGammaWqq() const
A get method to retrieve the theoretical uncertainty in , denoted as .
const double getIntegrand_dsigmaBox_charm200(double x) const
const double v() const
The Higgs vacuum expectation value.
const double getIntegrand_sigmaWithISR_up196(double x) const
const double getIntegrand_dsigmaBox_up196(double x) const
const double getIntegrand_sigmaWithISR_bottom189(double x) const
const double getIntegrand_AFBnumeratorWithISR_bottom183(double x) const
virtual const double LEP2dsigmadcosE(const double s, const double cos) const
virtual bool Update(const std::map< std::string, double > &DPars)
The update method for StandardModel.
const gslpp::matrix< gslpp::complex > & getYn() const
A get method to retrieve the Yukawa matrix of the neutrinos, .
const double getIntegrand_sigmaWithISR_bottom161(double x) const
const double getIntegrand_dsigmaBox_tau167(double x) const
gslpp::complex AHZga_f(const double tau, const double lambda) const
Fermionic loop function entering in the calculation of the effective coupling.
const double getIntegrand_sigmaWithISR_mu192(double x) const
const double getIntegrand_dsigmaBox_up189(double x) const
const double MwbarFromMw(const double Mw) const
A method to convert the -boson mass in the experimental/running-width scheme to that in the complex-p...
const double getIntegrand_AFBnumeratorWithISR_mu183(double x) const
const double getIntegrand_dsigmaBox_up136(double x) const
const double getDelGammaWlv() const
A get method to retrieve the theoretical uncertainty in , denoted as .
bool useKappaZ_f_cache[12]
const double getDelGammaZ() const
A get method to retrieve the theoretical uncertainty in , denoted as .
virtual const double eeffRelectrontsub(const double pol_e, const double pol_p, const double s) const
const double getIntegrand_dsigmaBox_charm202(double x) const
const double computeGammaHZga_tW() const
The mixed loop contribution to in the Standard Model.
virtual const double ThetaRnuN() const
The effective neutrino nucleon RH parameter: ThetaRnuN.
const double intMLL2eeeeus2(const double s, const double t0, const double t1) const
EWSMTwoFermionsLEP2 * getMyTwoFermionsLEP2() const
A get method to retrieve the member pointer of type EWSMTwoFermionsLEP2.
const double DeltaAlphaLepton(const double s) const
Leptonic contribution to the electromagnetic coupling , denoted as .
const double getIntegrand_sigmaWithISR_charm189(double x) const
const double getIntegrand_sigmaWithISR_down205(double x) const
virtual const gslpp::complex gV_f(const Particle f) const
The effective leptonic neutral-current vector coupling in the SM.
const double getIntegrand_dsigmaBox_up172(double x) const
const double getIntegrand_sigmaWithISR_up130(double x) const
const double MLR2eeff(const Particle f, const double s) const
virtual const double gAnue() const
The effective (muon) neutrino-electron axial-vector coupling: gAnue.
void setYe(const gslpp::matrix< gslpp::complex > &Ye)
A set method to set the Yukawa matrix of the charged leptons, .
bool FlagCacheInStandardModel
A flag for caching (true by default).
bool SMSuccess
A boolean for the success of the Standard Model update and matching.
virtual bool InitializeModel()
A method to initialize the model.
double resumMw(const double Mw_i, const double DeltaRho[orders_EW_size], const double DeltaR_rem[orders_EW_size]) const
A method to compute the -boson mass from and .
const double getAle() const
A get method to retrieve the fine-structure constant .
const double getIntegrand_AFBnumeratorWithISR_charm192(double x) const
const double getIntegrand_dsigmaBox_down205(double x) const
virtual const double cW2() const
const double getIntegrand_dsigmaBox_down192(double x) const
virtual const double GammaHTot() const
The total Higgs width in the Standard Model.
virtual const double getCepsK() const
The ratio of the imaginary part of the $K$ mixing amplitude over the Standard Model value.
EWSMThreeLoopEW * getMyThreeLoopEW() const
const double computeBrHtoWW() const
The Br in the Standard Model.
virtual const double LEP2AFBbottom(const double s) const
double DeltaAlpha_cache
A cache of the value of .
virtual const double gVnue() const
The effective (muon) neutrino-electron vector coupling: gVnue.
const double AlsEWithInit(double mu, double alsi, double mu_i, const int nf_i, orders order) const
double delR0l
The theoretical uncertainty in , denoted as .
virtual const double LEP2dsigmadcosBinMu(const double s, const double cos, const double cosmin, const double cosmax) const
const double getIntegrand_sigmaWithISR_down161(double x) const
virtual const double eeffRtau(const double pol_e, const double pol_p, const double s) const
const double getIntegrand_dsigmaBox_bottom136(double x) const
A class for the matching in the Standard Model.
Test Observable.
Test Observable.
orders
An enum type for orders in QCD.
Definition: OrderScheme.h:33