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
509class StandardModel : public QCD {
510public:
511
512
513 // Radiative Corrections for the LEP-II observables
514 enum LEP2RCs {
515 Weak = 0,
521 };
522
528 EW1 = 0,
535 };
536
537
538
539// static const int NSMvars = 39; ///< The number of the model parameters in %StandardModel. !!! PMNS INCLUDED
540 static const int NSMvars = 26;
544 static std::string SMvars[NSMvars];
545
546 static const double GeVminus2_to_nb;
547
552 static const double Mw_error;
553
558
562 virtual ~StandardModel();
563
564
566 // Initialization
567
574 virtual bool InitializeModel();
575
576
578 // Model parameters
579
586 virtual bool Init(const std::map<std::string, double>& DPars);
587
595 virtual bool PreUpdate();
596
604 virtual bool Update(const std::map<std::string, double>& DPars);
605
615 virtual bool PostUpdate();
616
617 const int getIterationNo() const
618 {
619 return iterationNo;
620 }
621
629 virtual bool CheckParameters(const std::map<std::string, double>& DPars);
630
631
633 // Flags
634
641 virtual bool setFlag(const std::string name, const bool value);
642
649 virtual bool setFlagStr(const std::string name, const std::string value);
650
655 virtual bool CheckFlags() const;
656
667 {
669 }
670
679 const bool IsFlagNoApproximateGammaZ() const
680 {
682 }
683
685 {
686 this->FlagNoApproximateGammaZ = FlagNoApproximateGammaZ;
687 }
688
694 const std::string getFlagMw() const
695 {
696 return FlagMw;
697 }
698
704 const std::string getFlagRhoZ() const
705 {
706 return FlagRhoZ;
707 }
708
714 const std::string getFlagKappaZ() const
715 {
716 return FlagKappaZ;
717 }
718
731 {
732 this->FlagCacheInStandardModel = FlagCacheInStandardModel;
733 }
734
735
737 // get and set methods for class members
738
744 const Particle& getLeptons(const QCD::lepton p) const
745 {
746 return leptons[p];
747 }
748
753 const double getMz() const
754 {
755 return Mz;
756 }
757
762 const double getMw() const
763 {
764 return Mw_inp;
765 }
766
771 const double getAlsMz() const
772 {
773 return AlsMz;
774 }
775
780 const double getGF() const
781 {
782 return GF;
783 }
784
789 const double getAle() const
790 {
791 return ale;
792 }
793
800 const double getDAle5Mz() const
801 {
802 return dAle5Mz;
803 }
804
809 virtual const double getMHl() const
810 {
811 return mHl;
812 }
813
819 const double getDelMw() const
820 {
821 return delMw;
822 }
823
830 const double getDelSin2th_l() const
831 {
832 return delSin2th_l;
833 }
834
841 const double getDelSin2th_q() const
842 {
843 return delSin2th_q;
844 }
845
852 const double getDelSin2th_b() const
853 {
854 return delSin2th_b;
855 }
856
862 const double getDelGammaZ() const
863 {
864 return delGammaZ;
865 }
866
872 const double getDelSigma0H() const
873 {
874 return delsigma0H;
875 }
876
882 const double getDelR0l() const
883 {
884 return delR0l;
885 }
886
892 const double getDelR0c() const
893 {
894 return delR0c;
895 }
896
902 const double getDelR0b() const
903 {
904 return delR0b;
905 }
906
911 const gslpp::matrix<gslpp::complex> getVCKM() const
912 {
913 return myCKM.getCKM();
914 }
915
920 const CKM& getCKM() const
921 {
922 return myCKM;
923 }
924
925
926
931 const gslpp::matrix<gslpp::complex> getUPMNS() const
932 {
933 return myPMNS.getPMNS();
934 }
935
941 const gslpp::matrix<gslpp::complex>& getYn() const
942 {
943 return Yn;
944 }
945
951 const double getMuw() const
952 {
953 return muw;
954 }
955
956 virtual const StandardModel& getTrueSM() const
957 {
958 throw std::runtime_error("StandardModel::getTrueSM() must be overridden by the NP extension.");
959 }
960
966 {
967 return SMM.getObj();
968 }
969
975 {
976 return myEWSMcache;
977 }
978
984 {
985 return myOneLoopEW;
986 }
987
993 {
995 }
996
997 /* BEGIN: REMOVE FROM THE PACKAGE */
1003 {
1004 return myTwoFermionsLEP2;
1005 }
1006 /* END: REMOVE FROM THE PACKAGE */
1007
1009 {
1010 return myThreeLoopEW;
1011 }
1012
1014 {
1015 return myThreeLoopEW2QCD;
1016 }
1017
1019 {
1020 return myThreeLoopQCD;
1021 }
1022
1024 {
1025 return myTwoLoopEW;
1026 }
1027
1029 {
1030 return myTwoLoopQCD;
1031 }
1032
1033 const Flavour& getFlavour() const
1034 {
1035 return SMFlavour;
1036 }
1037
1039 {
1040 return myLeptonFlavour;
1041 }
1043 // QED coupling
1044
1056 const double ale_OS(const double mu, orders order = FULLNLO) const;
1057
1064 const double Beta_s(int nm, unsigned int nf) const;
1065
1072 const double Beta_e(int nm, unsigned int nf) const;
1073
1074 //make the QCD implementation of Als available here
1075 using QCD::Als;
1076
1088 const double Als(const double mu, const orders order, const bool Nf_thr, const bool qed_flag) const
1089 {
1090 if (qed_flag && order == FULLNNNLO)
1091 return AlsE(mu, order, Nf_thr);
1092
1093 return Als(mu, order, Nf_thr);
1094 }
1095
1105 const double Ale(double mu, orders order, bool Nf_thr = true) const;
1106
1113 const double DeltaAlphaLepton(const double s) const;
1114
1127 const double DeltaAlphaL5q() const;
1128
1135 const double DeltaAlphaTop(const double s) const;
1136
1149 const double DeltaAlpha() const;
1150
1160 virtual const double alphaMz() const;
1161
1168 const double Alstilde5(const double mu) const;
1169
1173 virtual const double getCCC1() const { return 0.; };
1174
1178 virtual const double getCCC2() const { return 0.; };
1179
1183 virtual const double getCCC3() const { return 0.; };
1184
1188 virtual const double getCCC4() const { return 0.; };
1189
1193 virtual const double getCCC5() const { return 0.; };
1194
1195
1196
1198 // Higgs VEV
1199
1208 const double v() const;
1209
1210
1212 // The W-boson mass
1213
1218 const double Mw_tree() const;
1219
1235 const double s02() const;
1236
1249 const double c02() const;
1250
1277 virtual const double Mw() const;
1278
1279
1288 virtual const double Dalpha5hMz() const;
1289
1299 virtual const double cW2(const double Mw_i) const;
1300 virtual const double cW2() const;
1301
1311 virtual const double sW2(const double Mw_i) const;
1312 const double sW2() const;
1313
1320 const double sW2_MSbar_Approx() const;
1321
1328 const double sW2_ND() const;
1329
1330
1331
1353 virtual const double DeltaR() const;
1354
1363 void ComputeDeltaRho(const double Mw_i, double DeltaRho[orders_EW_size]) const;
1364
1373 void ComputeDeltaR_rem(const double Mw_i, double DeltaR_rem[orders_EW_size]) const;
1374
1375
1377 // The W and Z masses in the complex-pole/fixed-width scheme
1378
1409 double Mzbar() const;
1410
1431 const double MwbarFromMw(const double Mw) const;
1432
1456 const double MwFromMwbar(const double Mwbar) const;
1457
1475 virtual const double DeltaRbar() const;
1476
1477
1479 // The W-boson decay width
1480
1489 virtual const double rho_GammaW(const Particle fi, const Particle fj) const;
1490
1517 virtual const double GammaW(const Particle fi, const Particle fj) const;
1518
1523 virtual const double GammaW() const;
1524
1529 virtual const double BrW(const Particle fi, const Particle fj) const;
1530
1535 virtual const double RWlilj(const Particle li, const Particle lj) const;
1536
1537 virtual const double Ruc() const; //AG:added
1538
1543 virtual const double RWc() const;
1544
1546 // EWPO at Z-pole
1547
1566 virtual const double A_f(const Particle f) const;
1567
1573 virtual const double AFB(const Particle f) const;
1574
1594 virtual const double sin2thetaEff(const Particle f) const;
1595
1621 virtual const double GammaZ(const Particle f) const;
1622
1633 virtual const double Gamma_inv() const;
1634
1649 virtual const double Gamma_had() const;
1650
1664 virtual const double Gamma_Z() const;
1665
1670 virtual const double RZlilj(const Particle li, const Particle lj) const;
1671
1685 virtual const double sigma0_had() const;
1686
1701 virtual const double R0_f(const Particle f) const;
1702
1711 virtual const double R_inv() const;
1712
1722 virtual const double N_nu() const;
1723
1724
1726 // Zff effective couplings
1727
1737 virtual const gslpp::complex gV_f(const Particle f) const;
1738
1748 virtual const gslpp::complex gA_f(const Particle f) const;
1749
1765 virtual const gslpp::complex rhoZ_f(const Particle f) const;
1766
1794 virtual const gslpp::complex kappaZ_f(const Particle f) const;
1795
1816 virtual const gslpp::complex deltaRhoZ_f(const Particle f) const;
1817
1842 virtual const gslpp::complex deltaKappaZ_f(const Particle f) const;
1843
1844
1846 // Epsilon parameters for EWPO
1847
1859 virtual const double epsilon1() const;
1860
1881 virtual const double epsilon2() const;
1882
1900 virtual const double epsilon3() const;
1901
1919 virtual const double epsilonb() const;
1920
1921
1923 // EW low-energy observables: Parity violation
1924
1925
1931 virtual const double amuon() const;
1932
1940 virtual const double Qwemoller(const double q2, const double y) const;
1941
1942
1950 virtual const double alrmoller(const double q2, const double y) const;
1951
1952
1959 virtual const double Qwp() const;
1960
1961
1968 virtual const double Qwn() const;
1969
1971 // EW low-energy observables: neutrino-scattering
1972
1979 virtual const double gLnuN2() const;
1980
1981
1988 virtual const double gRnuN2() const;
1989
1990
1997 virtual const double ThetaLnuN() const;
1998
1999
2006 virtual const double ThetaRnuN() const;
2007
2008
2015 virtual const double gVnue() const;
2016
2017
2024 virtual const double gAnue() const;
2025
2026
2028 // Lepton decays
2029
2030 // Muon decays
2031
2037 virtual const double Gamma_muon() const;
2038
2039
2040 // Tau decays
2041
2047 virtual const double Gamma_tau_l_nunu(const Particle l) const;
2048
2049
2050 // Lepton Flavor universality tests in Tau decays
2051
2057 virtual const double TauLFU_gmuge() const;
2058
2059
2065 virtual const double TauLFU_gtaugmu() const;
2066
2067
2073 virtual const double TauLFU_gtauge() const;
2074
2080 virtual const double TauLFU_gtaugmuPi() const;
2081
2087 virtual const double TauLFU_gtaugmuK() const;
2088
2089
2091 // For EWPO caches
2092
2100 static const int NumSMParamsForEWPO = 33;
2101
2116 bool checkSMparamsForEWPO();
2117
2119 // Several Higgs-related quantities used in Higgs coupling analysis
2121
2127 gslpp::complex f_triangle(const double tau) const;
2133 gslpp::complex g_triangle(const double tau) const;
2139 gslpp::complex I_triangle_1(const double tau, const double lambda) const;
2145 gslpp::complex I_triangle_2(const double tau, const double lambda) const;
2153 gslpp::complex AH_f(const double tau) const;
2154
2162 gslpp::complex AH_W(const double tau) const;
2163
2170 gslpp::complex AHZga_f(const double tau, const double lambda) const;
2171
2178 gslpp::complex AHZga_W(const double tau, const double lambda) const;
2179
2181
2187 virtual const double SigmaeeZH(const double sqrt_s, const double Pe, const double Pp) const;
2188
2194 virtual const double SigmaeeHvv(const double sqrt_s, const double Pe, const double Pp) const;
2195
2201 virtual const double SigmaeeHee(const double sqrt_s, const double Pe, const double Pp) const;
2202
2204
2210 virtual const double GammaHtogg() const;
2211
2217 virtual const double GammaHtoZZstar() const;
2218
2224 virtual const double GammaHtoWWstar() const;
2225
2231 virtual const double GammaHtoZga() const;
2232
2238 virtual const double GammaHtogaga() const;
2239
2245 virtual const double GammaHtomumu() const;
2246
2252 virtual const double GammaHtotautau() const;
2253
2259 virtual const double GammaHtocc() const;
2260
2266 virtual const double GammaHtoss() const;
2267
2273 virtual const double GammaHtobb() const;
2274
2280 virtual const double GammaHTot() const;
2281
2283
2289 virtual const double BrHtogg() const;
2290
2296 virtual const double BrHtoZZstar() const;
2297
2303 virtual const double BrHtoWWstar() const;
2304
2310 virtual const double BrHtoZga() const;
2311
2317 virtual const double BrHtogaga() const;
2318
2324 virtual const double BrHtomumu() const;
2325
2331 virtual const double BrHtotautau() const;
2332
2338 virtual const double BrHtocc() const;
2339
2345 virtual const double BrHtoss() const;
2346
2352 virtual const double BrHtobb() const;
2353
2355
2371 const double computeSigmaggH(const double sqrt_s) const
2372 {
2373 if (sqrt_s == 7.0) {
2374 return 16.83; // in pb for Mh=125.1 GeV
2375 } else if (sqrt_s == 8.0) {
2376 return 21.40; // in pb for Mh=125.1 GeV
2377 } else if (sqrt_s == 13.0) {
2378 return 48.61; // in pb for Mh=125.09 GeV
2379 } else if (sqrt_s == 13.6) {
2380 return 52.17; // in pb for Mh=125.09 GeV
2381 } else if (sqrt_s == 14.0) {
2382 return 54.72; // in pb for Mh=125.09 GeV
2383 } else if (sqrt_s == 27.0) {
2384 return 146.65; // in pb for Mh=125.09 GeV
2385 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2386 return 740.3; // in pb for Mh=125. GeV
2387 } else if (sqrt_s == 1.96) {
2388 return 0.9493; // in pb for Mh=125 GeV
2389 } else
2390 throw std::runtime_error("Bad argument in StandardModel::computeSigmaggH()");
2391 }
2392
2403 const double computeSigmaggH_tt(const double sqrt_s) const
2404 {
2405 if (sqrt_s == 7.0) {
2406 return 16.69; // in pb for Mh=125.09 GeV
2407 } else if (sqrt_s == 8.0) {
2408 return 21.20; // in pb for Mh=125.09 GeV
2409 } else if (sqrt_s == 13.0) {
2410 return 47.94; // in pb for Mh=125.09 GeV
2411 } else if (sqrt_s == 13.6) {
2412 return 51.534; // in pb for Mh=125.09 GeV (interpolation between 13 and 14 TeV)
2413 } else if (sqrt_s == 14.0) {
2414 return 53.93; // in pb for Mh=125.09 GeV
2415 } else if (sqrt_s == 27.0) {
2416 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_tt(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2417 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2418 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_tt(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2419 } else
2420 throw std::runtime_error("Bad argument in StandardModel::computeSigmaggH_tt()");
2421 }
2422
2433 const double computeSigmaggH_bb(const double sqrt_s) const
2434 {
2435 if (sqrt_s == 7.0) {
2436 return 0.04; // in pb for Mh=125.09 GeV
2437 } else if (sqrt_s == 8.0) {
2438 return 0.05; // in pb for Mh=125.09 GeV
2439 } else if (sqrt_s == 13.0) {
2440 return 0.10; // in pb for Mh=125.09 GeV
2441 } else if (sqrt_s == 13.6) {
2442 return 0.106; // in pb for Mh=125.09 GeV (interpolation between 13 and 14 TeV)
2443 } else if (sqrt_s == 14.0) {
2444 return 0.11; // in pb for Mh=125.09 GeV
2445 } else if (sqrt_s == 27.0) {
2446 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_bb(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2447 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2448 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_bb(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2449 } else
2450 throw std::runtime_error("Bad argument in StandardModel::computeSigmaggH_bb()");
2451 }
2452
2463 const double computeSigmaggH_tb(const double sqrt_s) const
2464 {
2465 if (sqrt_s == 7.0) {
2466 return -0.66; // in pb for Mh=125.09 GeV
2467 } else if (sqrt_s == 8.0) {
2468 return -0.82; // in pb for Mh=125.09 GeV
2469 } else if (sqrt_s == 13.0) {
2470 return -1.73; // in pb for Mh=125.09 GeV
2471 } else if (sqrt_s == 13.6) {
2472 return -1.844; // in pb for Mh=125.09 GeV (interpolation between 13 and 14 TeV)
2473 } else if (sqrt_s == 14.0) {
2474 return -1.92; // in pb for Mh=125.09 GeV
2475 } else if (sqrt_s == 27.0) {
2476 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_tb(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2477 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2478 return computeSigmaggH(sqrt_s) / computeSigmaggH(14.) * computeSigmaggH_tb(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2479 } else
2480 throw std::runtime_error("Bad argument in StandardModel::computeSigmaggH_tb()");
2481 }
2482
2496 const double computeSigmaVBF(const double sqrt_s) const
2497 {
2498 if (sqrt_s == 7.0) {
2499 return 1.241; // in pb for Mh=125.09 GeV
2500 } else if (sqrt_s == 8.0) {
2501 return 1.601; // in pb for Mh=125.09 GeV
2502 } else if (sqrt_s == 13.0) {
2503 return 3.766; // in pb for Mh=125.09 GeV
2504 } else if (sqrt_s == 13.6) {
2505 return 4.075; // in pb for Mh=125.09 GeV
2506 } else if (sqrt_s == 14.0) {
2507 return 4.260; // in pb for Mh=125.09 GeV
2508 } else if (sqrt_s == 27.0) {
2509 return 11.838; // in pb for Mh=125.09 GeV
2510 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2511 return 82.0; // in pb for Mh=125. GeV
2512 } else if (sqrt_s == 1.96) {
2513 return 0.0653; // in pb for Mh=125 GeV
2514 } else
2515 throw std::runtime_error("Bad argument in StandardModel::computeSigmaVBF()");
2516 }
2517
2529 const double computeSigmaWF(const double sqrt_s) const
2530 {
2531 if (sqrt_s == 7.0) {
2532 return 0.946; // in pb for Mh=125 GeV
2533 } else if (sqrt_s == 8.0) {
2534 return 1.220; // in pb for Mh=125 GeV
2535 } else if (sqrt_s == 13.0) {
2536 return 2.882; // in pb for Mh=125 GeV
2537 } else if (sqrt_s == 13.6) {
2538 return 3.1088; // in pb for Mh=125 GeV (interpolation between 13 and 14 TeV)
2539 } else if (sqrt_s == 14.0) {
2540 return 3.260; // in pb for Mh=125 GeV
2541 } else if (sqrt_s == 27.0) {
2542 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(14.) * computeSigmaWF(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2543 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2544 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(14.) * computeSigmaWF(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2545 } else if (sqrt_s == 1.96) {
2546 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(7.) * computeSigmaWF(7.); // in the absence of individual cross sections for TeVatron we rescale the LHC ones
2547 } else
2548 throw std::runtime_error("Bad argument in StandardModel::computeSigmaWF()");
2549 }
2550
2562 const double computeSigmaZF(const double sqrt_s) const
2563 {
2564 if (sqrt_s == 7.0) {
2565 return 0.333; // in pb for Mh=125 GeV
2566 } else if (sqrt_s == 8.0) {
2567 return 0.432; // in pb for Mh=125 GeV
2568 } else if (sqrt_s == 13.0) {
2569 return 1.049; // in pb for Mh=125 GeV
2570 } else if (sqrt_s == 13.6) {
2571 return 1.1342; // in pb for Mh=125 GeV (interpolation between 13 and 14 TeV)
2572 } else if (sqrt_s == 14.0) {
2573 return 1.191; // in pb for Mh=125 GeV
2574 } else if (sqrt_s == 27.0) {
2575 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(14.) * computeSigmaZF(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2576 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2577 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(14.) * computeSigmaZF(14.); // in the absence of this value we rescale the LHC result at 14 TeV
2578 } else if (sqrt_s == 1.96) {
2579 return computeSigmaVBF(sqrt_s) / computeSigmaVBF(7.) * computeSigmaZF(7.); // in the absence of individual cross sections for TeVatron we rescale the LHC ones
2580 } else
2581 throw std::runtime_error("Bad argument in StandardModel::computeSigmaZF()");
2582 }
2583
2591 const double computeSigmaZWF(const double sqrt_s) const
2592 {
2593 return 0.;
2594 }
2595
2609 const double computeSigmaWH(const double sqrt_s) const
2610 {
2611 if (sqrt_s == 7.0) {
2612 return 0.577; // in pb for Mh=125.1 GeV
2613 //return 0.5688; // in pb for Mh=125.6 GeV
2614 } else if (sqrt_s == 8.0) {
2615 return 0.7027; // in pb for Mh=125.1 GeV
2616 //return 0.6931; // in pb for Mh=125.6 GeV
2617 } else if (sqrt_s == 13.0) {
2618 return 1.358; // in pb for Mh=125.09 GeV
2619 } else if (sqrt_s == 13.6) {
2620 return 1.453; // in pb for Mh=125.09 GeV
2621 } else if (sqrt_s == 14.0) {
2622 return 1.498; // in pb for Mh=125.09 GeV
2623 } else if (sqrt_s == 27.0) {
2624 return 3.397; // in pb for Mh=125.09 GeV
2625 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2626 return 15.9; // in pb for Mh=125. GeV
2627 } else if (sqrt_s == 1.96) {
2628 return 0.1295; // in pb for Mh=125 GeV
2629 } else
2630 throw std::runtime_error("Bad argument in StandardModel::computeSigmaWH()");
2631 }
2632
2646 const double computeSigmaZH(const double sqrt_s) const
2647 {
2648 if (sqrt_s == 7.0) {
2649 return 0.3341; // in pb for Mh=125.1 GeV
2650 //return 0.3299; // in pb for Mh=125.6 GeV
2651 } else if (sqrt_s == 8.0) {
2652 return 0.4142; // in pb for Mh=125.1 GeV
2653 //return 0.4091; // in pb for Mh=125.6 GeV
2654 } else if (sqrt_s == 13.0) {
2655 return 0.880; // in pb for Mh=125.09 GeV
2656 } else if (sqrt_s == 13.6) {
2657 return 0.9422; // in pb for Mh=125.09 GeV
2658 } else if (sqrt_s == 14.0) {
2659 return 0.981; // in pb for Mh=125.09 GeV
2660 } else if (sqrt_s == 27.0) {
2661 return 2.463; // in pb for Mh=125.09 GeV
2662 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2663 return 11.26; // in pb for Mh=125. GeV
2664 } else if (sqrt_s == 1.96) {
2665 return 0.0785; // in pb for Mh=125 GeV
2666 } else
2667 throw std::runtime_error("Bad argument in StandardModel::computeSigmaZH()");
2668 }
2669
2686 const double computeSigmattH(const double sqrt_s) const
2687 {
2688 if (sqrt_s == 7.0) {
2689 return 0.0861; // in pb for Mh=125.1 GeV
2690 //return 0.0851; // in pb for Mh=125.6 GeV
2691 } else if (sqrt_s == 8.0) {
2692 return 0.129; // in pb for Mh=125.1 GeV
2693 //return 0.1274; // in pb for Mh=125.6 GeV
2694 } else if (sqrt_s == 13.0) {
2695 return 0.5060; // in pb for Mh=125.1 GeV
2696 } else if (sqrt_s == 13.6) {
2697 return 0.5688; // in pb for Mh=125.09 GeV
2698 } else if (sqrt_s == 14.0) {
2699 return 0.6128; // in pb for Mh=125.09 GeV
2700 } else if (sqrt_s == 27.0) {
2701 return 2.86; // in pb for Mh=125.09 GeV
2702 } else if ( (sqrt_s == 84.0) || (sqrt_s == 100.0) ) {
2703 return 37.9; // in pb for Mh=125. GeV
2704 } else if (sqrt_s == 1.96) {
2705 return 0.0043; // in pb for Mh=125 GeV
2706 } else
2707 throw std::runtime_error("Bad argument in StandardModel::computeSigmattH()");
2708 }
2709
2716 const double computeSigmatHq(const double sqrt_s) const
2717 {
2718 if (sqrt_s == 13.0) {
2719 return 0.07426; // in pb for Mh=125.09 GeV
2720 } else
2721 throw std::runtime_error("Bad argument in StandardModel::computeSigmatHq()");
2722 }
2723
2724
2725
2735 const double computeSigmabbH(const double sqrt_s) const
2736 {
2737 if (sqrt_s == 13.0){
2738 return 0.4863; // in pb for Mh=125.09 GeV
2739 }
2740 else if (sqrt_s == 13.6) {
2741 return 0.566; // in pb for Mh=125.09 GeV (NLO+NNLLpart+ybyt matching)
2742 } else
2743 throw std::runtime_error("Bad argument in StandardModel::computeSigmabbH()");
2744 }
2745
2752 const double computeBrHtogg() const
2753 {
2754 return 8.179e-2; // Mh=125.1 GeV
2755 }
2756
2763 const double computeBrHtoWW() const
2764 {
2765 //return 2.23e-1; // Mh=125.5 GeV
2766 return 2.154e-1; // Mh=125.1 GeV
2767 }
2768
2775 const double computeBrHtoZZ() const
2776 {
2777 return 2.643e-2; // Mh=125.1 GeV
2778 //return 2.79e-2; // Mh=125.6 GeV
2779 }
2780
2787 const double computeBrHtoZga() const
2788 {
2789 return 1.541e-3; // Mh=125.1 GeV
2790 //return 1.59e-3; // Mh=125.6 GeV
2791 }
2792
2799 const double computeBrHtogaga() const
2800 {
2801 return 2.27e-3; // Mh=125.1 GeV
2802 }
2803
2810 const double computeBrHtomumu() const
2811 {
2812 return 2.17e-4; // Mh=125.1 GeV
2813 }
2814
2821 const double computeBrHtotautau() const
2822 {
2823 return 6.256e-2; // Mh=125.1 GeV
2824 //return 6.22e-2; // Mh=125.6 GeV
2825 }
2826
2833 const double computeBrHtocc() const
2834 {
2835 return 2.883e-2; // Mh=125.1 GeV
2836 //return 2.86e-2; // Mh=125.6 GeV
2837 }
2838
2845 const double computeBrHtoss() const
2846 {
2847 return 4.0e-4;
2848 }
2849
2856 const double computeBrHtobb() const
2857 {
2858 return 5.807e-1; // Mh=125.1 GeV
2859 //return 5.67e-1; // Mh=125.6 GeV
2860 }
2861
2862 // H decays into 4-fermions
2863
2869 const double computeBrHto4v() const
2870 {
2871 return 1.06e-3;
2872 }
2873
2880 const double computeBrHtoevmuv() const
2881 {
2882 return 2.539e-03; // Mh=125.1 GeV
2883 }
2884
2891 const double computeBrHtollvv2() const
2892 {
2893 return 1.063e-02; // Mh=125.1 GeV
2894 }
2895
2902 const double computeBrHtollvv3() const
2903 {
2904 return 2.356e-02; // Mh=125.1 GeV
2905 }
2906
2913 const double computeBrHto4l2() const
2914 {
2915 return 1.252e-04; // Mh=125.1 GeV
2916 }
2917
2924 const double computeBrHto4l3() const
2925 {
2926 return 2.771e-04; // Mh=125.1 GeV
2927 }
2928
2935 const double computeBrHto4q() const
2936 {
2937 return 1.098e-01; // Mh=125.1 GeV
2938 }
2939
2946 const double computeBrHto4f() const
2947 {
2948 return 2.406e-01; // Mh=125.1 GeV
2949 }
2950
2957 const double computeGammaHTotal() const
2958 {
2959 return 4.101e-3; // Mh=125.1 GeV
2960 //return 4.15e-3; // Mh=125.6 GeV
2961 }
2962
2968 const double computeGammaHgg_tt() const
2969 {
2970 return 380.8; // in keV for Mh=125 GeV
2971 //return 389.6; // in keV for Mh=126 GeV
2972 }
2973
2979 const double computeGammaHgg_bb() const
2980 {
2981 return 3.96; // in keV for Mh=125 GeV
2982 //return 3.95; // in keV for Mh=126 GeV
2983 }
2984
2990 const double computeGammaHgg_tb() const
2991 {
2992 return -42.1; // in keV for Mh=125 GeV
2993 //return -42.7; // in keV for Mh=126 GeV
2994 }
2995
3001 const double computeGammaHZga_tt() const
3002 {
3003 return 21.74; // in eV for Mh=125 GeV
3004 //return 23.51; // in eV for Mh=126 GeV
3005 }
3006
3012 const double computeGammaHZga_WW() const
3013 {
3014 return 7005.6; // in eV for Mh=125 GeV
3015 //return 7648.4; // in eV for Mh=126 GeV
3016 }
3017
3023 const double computeGammaHZga_tW() const
3024 {
3025 return -780.4; // in eV for Mh=125 GeV
3026 //return -848.1; // in eV for Mh=126 GeV
3027 }
3028
3034 const double computeGammaHgaga_tt() const
3035 {
3036 return 662.84; // in eV for Mh=125 GeV
3037 //return 680.39; // in eV for Mh=126 GeV
3038 }
3039
3045 const double computeGammaHgaga_WW() const
3046 {
3047 return 14731.86; // in eV for Mh=125 GeV
3048 //return 15221.98; // in eV for Mh=126 GeV
3049 }
3050
3056 const double computeGammaHgaga_tW() const
3057 {
3058 return -6249.93; // in eV for Mh=125 GeV
3059 //return -6436.35; // in eV for Mh=126 GeV
3060 }
3061
3066 virtual const double getCBd() const
3067 {
3068 return 1.;
3069 }
3070
3075 virtual const double getCBs() const
3076 {
3077 return 1.;
3078 }
3079
3084 virtual const double getCDMK() const
3085 {
3086 return 1.;
3087 }
3088
3093 virtual const double getCepsK() const
3094 {
3095 return 1.;
3096 }
3097
3102 virtual const double getPhiBs() const
3103 {
3104 return 0.;
3105 }
3106
3111 virtual const double getPhiBd() const
3112 {
3113 return 0.;
3114 }
3115
3117 //Generic e+e- -> ff Inclusive Observables
3118
3119// For f!=e
3120
3121// Helicity amplitudes squared
3122 const double MLR2eeff(const Particle f, const double s) const;
3123 const double MRL2eeff(const Particle f, const double s) const;
3124
3125 const double MLL2eeff(const Particle f, const double s, const double t) const;
3126 const double MRR2eeff(const Particle f, const double s, const double t) const;
3127
3128// Some simple functions for cos \theta integrals
3129 const double tovers2(const double cosmin, const double cosmax) const;
3130 const double uovers2(const double cosmin, const double cosmax) const;
3131
3132// For f=e
3133
3134// Integrals of the SM squared amplitudes x (t/s)^2, (s/t)^2, (u/s)^2 in [t0, t1]
3135 const double intMLR2eeeets2(const double s, const double t0, const double t1) const;
3136 const double intMLRtilde2eeeest2(const double s, const double t0, const double t1) const;
3137 const double intMLL2eeeeus2(const double s, const double t0, const double t1) const;
3138 const double intMRR2eeeeus2(const double s, const double t0, const double t1) const;
3139
3140// Cross sections
3141
3142 const double eeffsigmaEbin(const double pol_e, const double pol_p, const double s, const double cosmin, const double cosmax) const;
3143
3144 virtual const double eeffsigmaE(const double pol_e, const double pol_p, const double s) const
3145 {
3146 double cosmin = -0.90; // As in LEP2
3147 double cosmax = 0.90; // As in LEP2
3148
3149 return eeffsigmaEbin(pol_e, pol_p, s, cosmin, cosmax);
3150 }
3151
3152 const double eeffsigma(const Particle f, const double pol_e, const double pol_p, const double s, const double cosmin, const double cosmax) const;
3153
3154 virtual const double eeffsigmaEtsub(const double pol_e, const double pol_p, const double s) const
3155 {
3156 double cosmin = -0.90; // As in LEP2
3157 double cosmax = 0.90; // As in LEP2
3158
3159 return eeffsigma(leptons[ELECTRON], pol_e, pol_p, s, cosmin, cosmax);
3160 }
3161
3162 virtual const double eeffsigmaMu(const double pol_e, const double pol_p, const double s) const
3163 {
3164 return eeffsigma(leptons[MU], pol_e, pol_p, s, -1.0, 1.0);
3165 }
3166
3167 virtual const double eeffsigmaTau(const double pol_e, const double pol_p, const double s) const
3168 {
3169 return eeffsigma(leptons[TAU], pol_e, pol_p, s, -1.0, 1.0);
3170 }
3171
3172 virtual const double eeffsigmaHadron(const double pol_e, const double pol_p, const double s) const
3173 {
3174 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)
3175 + eeffsigma(quarks[CHARM], pol_e, pol_p, s, -1.0, 1.0) + eeffsigma(quarks[STRANGE], pol_e, pol_p, s, -1.0, 1.0)
3176 + eeffsigma(quarks[BOTTOM], pol_e, pol_p, s, -1.0, 1.0) );
3177 }
3178
3179 virtual const double eeffsigmaStrange(const double pol_e, const double pol_p, const double s) const
3180 {
3181 return eeffsigma(quarks[STRANGE], pol_e, pol_p, s, -1.0, 1.0);
3182 }
3183
3184 virtual const double eeffsigmaCharm(const double pol_e, const double pol_p, const double s) const
3185 {
3186 return eeffsigma(quarks[CHARM], pol_e, pol_p, s, -1.0, 1.0);
3187 }
3188
3189 virtual const double eeffsigmaBottom(const double pol_e, const double pol_p, const double s) const
3190 {
3191 return eeffsigma(quarks[BOTTOM], pol_e, pol_p, s, -1.0, 1.0);
3192 }
3193
3194// Ratios
3195
3196 virtual const double eeffRelectron(const double pol_e, const double pol_p, const double s) const
3197 {
3198 return ( eeffsigmaHadron(pol_e, pol_p, s) / eeffsigmaE(pol_e, pol_p, s) );
3199 }
3200 virtual const double eeffRelectrontsub(const double pol_e, const double pol_p, const double s) const
3201 {
3202 return ( eeffsigmaHadron(pol_e, pol_p, s) / eeffsigmaEtsub(pol_e, pol_p, s) );
3203 }
3204 virtual const double eeffRmuon(const double pol_e, const double pol_p, const double s) const
3205 {
3206 return ( eeffsigmaHadron(pol_e, pol_p, s) / eeffsigmaMu(pol_e, pol_p, s) );
3207 }
3208 virtual const double eeffRtau(const double pol_e, const double pol_p, const double s) const
3209 {
3210 return ( eeffsigmaHadron(pol_e, pol_p, s) / eeffsigmaTau(pol_e, pol_p, s) );
3211 }
3212 virtual const double eeffRstrange(const double pol_e, const double pol_p, const double s) const
3213 {
3214 return ( eeffsigmaStrange(pol_e, pol_p, s) / eeffsigmaHadron(pol_e, pol_p, s) );
3215 }
3216 virtual const double eeffRcharm(const double pol_e, const double pol_p, const double s) const
3217 {
3218 return ( eeffsigmaCharm(pol_e, pol_p, s) / eeffsigmaHadron(pol_e, pol_p, s) );
3219 }
3220 virtual const double eeffRbottom(const double pol_e, const double pol_p, const double s) const
3221 {
3222 return ( eeffsigmaBottom(pol_e, pol_p, s) / eeffsigmaHadron(pol_e, pol_p, s) );
3223 }
3224
3225// FB asymmetries
3226
3227 virtual const double eeffAFBe(const double pol_e, const double pol_p, const double s) const
3228 {
3229 double cosmin = -0.90; // As in LEP2
3230 double cosmax = 0.90; // As in LEP2
3231
3232 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));
3233 }
3234 virtual const double eeffAFBetsub(const double pol_e, const double pol_p, const double s) const
3235 {
3236 double cosmin = -0.90; // As in LEP2
3237 double cosmax = 0.90; // As in LEP2
3238
3239 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) );
3240 }
3241 virtual const double eeffAFBmu(const double pol_e, const double pol_p, const double s) const
3242 {
3243 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) );
3244 }
3245 virtual const double eeffAFBtau(const double pol_e, const double pol_p, const double s) const
3246 {
3247 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) );
3248 }
3249 virtual const double eeffAFBstrange(const double pol_e, const double pol_p, const double s) const
3250 {
3251 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) );
3252 }
3253 virtual const double eeffAFBcharm(const double pol_e, const double pol_p, const double s) const
3254 {
3255 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) );
3256 }
3257 virtual const double eeffAFBbottom(const double pol_e, const double pol_p, const double s) const
3258 {
3259 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) );
3260 }
3261
3262/* BEGIN: REMOVE FROM THE PACKAGE */
3264 //LEP2 Inclusive Observables
3265
3266 virtual const double LEP2sigmaE(const double s) const;
3267 virtual const double LEP2sigmaMu(const double s) const;
3268 virtual const double LEP2sigmaTau(const double s) const;
3269 virtual const double LEP2sigmaHadron(const double s) const;
3270 virtual const double LEP2sigmaCharm(const double s) const;
3271 virtual const double LEP2sigmaBottom(const double s) const;
3272
3273 virtual const double LEP2AFBe(const double s) const;
3274 virtual const double LEP2AFBmu(const double s) const;
3275 virtual const double LEP2AFBtau(const double s) const;
3276 virtual const double LEP2AFBcharm(const double s) const;
3277 virtual const double LEP2AFBbottom(const double s) const;
3278 virtual const double LEP2Rcharm(const double s) const;
3279 virtual const double LEP2Rbottom(const double s) const;
3280
3281 //LEP2 Differential Observables
3282 virtual const double LEP2dsigmadcosE(const double s, const double cos) const;
3283 virtual const double LEP2dsigmadcosMu(const double s, const double cos) const;
3284 virtual const double LEP2dsigmadcosTau(const double s, const double cos) const;
3285
3286 virtual const double LEP2dsigmadcosBinE(const double s, const double cos, const double cosmin, const double cosmax) const;
3287 virtual const double LEP2dsigmadcosBinMu(const double s, const double cos, const double cosmin, const double cosmax) const;
3288 virtual const double LEP2dsigmadcosBinTau(const double s, const double cos, const double cosmin, const double cosmax) const;
3289
3290/* END: REMOVE FROM THE PACKAGE */
3291
3292bool setFlagSigmaForAFB(const bool flagSigmaForAFB_i)
3293{
3294 bSigmaForAFB = flagSigmaForAFB_i;
3295 return true;
3296}
3297
3298bool setFlagSigmaForR(const bool flagSigmaForR_i)
3299{
3300 bSigmaForR = flagSigmaForR_i;
3301 return true;
3302}
3303
3310virtual const double getmq(const QCD::quark q, const double mu) const
3311{
3312 return m_q(q, mu, FULLNLO);
3313}
3314
3320 {
3321 this->requireCKM = requireCKM;
3322 }
3323
3328 void setCKM(const CKM& CKMMatrix)
3329 {
3330 myCKM = CKMMatrix;
3331 }
3332
3338 const gslpp::matrix<gslpp::complex>& getYu() const
3339 {
3340 return Yu;
3341 }
3342
3348 void setYu(const gslpp::matrix<gslpp::complex>& Yu)
3349 {
3350 this->Yu = Yu;
3351 }
3352
3358 const gslpp::matrix<gslpp::complex>& getYd() const
3359 {
3360 return Yd;
3361 }
3362
3368 void setYd(const gslpp::matrix<gslpp::complex>& Yd)
3369 {
3370 this->Yd = Yd;
3371 }
3372
3378 const gslpp::matrix<gslpp::complex>& getYe() const
3379 {
3380 return Ye;
3381 }
3382
3388 void setYe(const gslpp::matrix<gslpp::complex>& Ye)
3389 {
3390 this->Ye = Ye;
3391 }
3392
3397 const bool isSMSuccess() const
3398 {
3399 return SMSuccess;
3400 }
3401
3406 void setSMSuccess(bool success) const
3407 {
3408 SMSuccess = success;
3409 }
3410
3411
3413protected:
3414
3420 virtual void setParameter(const std::string name, const double& value);
3421
3425 virtual void computeCKM();
3426
3430 virtual void computeYukawas();
3431
3435// gslpp::matrix<gslpp::complex> VCKM; ///< The %CKM matrix.
3436// gslpp::matrix<gslpp::complex> UPMNS; ///< The %PMNS matrix.
3438 gslpp::matrix<gslpp::complex> Yu;
3439 gslpp::matrix<gslpp::complex> Yd;
3440 gslpp::matrix<gslpp::complex> Yn;
3441 gslpp::matrix<gslpp::complex> Ye;
3442
3443
3444 // model parameters
3445 double AlsMz;
3446 double Mz;
3447 double Mw_inp;
3448 double GF;
3449 double ale;
3450 double dAle5Mz;
3451 double mHl;
3452 double delMw;
3456 double delGammaZ;
3457 double delsigma0H;
3458 double delR0l;
3459 double delR0c;
3460 double delR0b;
3461 double lambda;
3462 double A;
3463 double rhob;
3464 double etab;
3465 double Vus;
3466 double Vud;
3467 double Vcb;
3468 double Vub;
3469 double gamma;
3470 double muw;
3472
3473 // non-model parameters
3474 double dAl5hMz;
3475
3477 // For EWPO
3478
3486
3497 double SchemeToDouble(const std::string scheme) const
3498 {
3499 if (scheme.compare("NORESUM") == 0)
3500 return 0.0;
3501 else if (scheme.compare("OMSI") == 0)
3502 return 1.0;
3503 else if (scheme.compare("INTERMEDIATE") == 0)
3504 return 2.0;
3505 else if (scheme.compare("OMSII") == 0)
3506 return 3.0;
3507 else if (scheme.compare("APPROXIMATEFORMULA") == 0)
3508 return 4.0;
3509 else
3510 throw std::runtime_error("EWSM::SchemeToDouble: bad scheme");
3511 }
3512
3518 bool checkEWPOscheme(const std::string scheme) const
3519 {
3520 if (scheme.compare("NORESUM") == 0
3521 || scheme.compare("OMSI") == 0
3522 || scheme.compare("INTERMEDIATE") == 0
3523 || scheme.compare("OMSII") == 0
3524 || scheme.compare("APPROXIMATEFORMULA") == 0)
3525 return true;
3526 else
3527 return false;
3528 }
3529
3530
3531 double m_q(const QCD::quark q, const double mu, const orders order=FULLNLO) const
3532 {
3533 switch(q) {
3534 case QCD::UP:
3535 case QCD::DOWN:
3536 case QCD::STRANGE:
3537 return Mrun(mu, getQuarks(q).getMass_scale(),
3538 getQuarks(q).getMass(), q, order);
3539 case QCD::CHARM:
3540 case QCD::BOTTOM:
3541 case QCD::TOP:
3542 return Mrun(mu, getQuarks(q).getMass(), q, order);
3543 default:
3544 throw std::runtime_error("Error in StandardModel::m_q()");
3545 }
3546 }
3547
3580 double resumMw(const double Mw_i, const double DeltaRho[orders_EW_size],
3581 const double DeltaR_rem[orders_EW_size]) const;
3582
3609 double resumRhoZ(const double DeltaRho[orders_EW_size],
3610 const double deltaRho_rem[orders_EW_size],
3611 const double DeltaRbar_rem, const bool bool_Zbb) const;
3612
3639 double resumKappaZ(const double DeltaRho[orders_EW_size],
3640 const double deltaKappa_rem[orders_EW_size],
3641 const double DeltaRbar_rem, const bool bool_Zbb) const;
3642
3665 double taub() const;
3666
3675 double Delta_EWQCD(const QCD::quark q) const;
3676
3686 double RVq(const QCD::quark q) const;
3687
3697 double RAq(const QCD::quark q) const;
3698
3712 double RVh() const;
3713
3717
3720
3721
3722 /* BEGIN: REMOVE FROM THE PACKAGE */
3724 //Migrated from LEP2ThObservables.h
3725
3726
3727
3729 mutable bool bSigmaForAFB;
3730 mutable bool bSigmaForR;
3731
3732
3733
3734 const double sigma_NoISR_l(const QCD::lepton l_flavor, const double s) const;
3735 const double sigma_NoISR_q(const QCD::quark q_flavor, const double s) const;
3736 const double AFB_NoISR_l(const QCD::lepton l_flavor, const double s) const;
3737 const double AFB_NoISR_q(const QCD::quark q_flavor, const double s) const;
3738
3739 const double Integrand_sigmaWithISR_l(double x, const QCD::lepton l_flavor, const double s) const;
3740
3741 const double getIntegrand_sigmaWithISR_mu130(double x) const;
3742 const double getIntegrand_sigmaWithISR_mu136(double x) const;
3743 const double getIntegrand_sigmaWithISR_mu161(double x) const;
3744 const double getIntegrand_sigmaWithISR_mu172(double x) const;
3745 const double getIntegrand_sigmaWithISR_mu183(double x) const;
3746 const double getIntegrand_sigmaWithISR_mu189(double x) const;
3747 const double getIntegrand_sigmaWithISR_mu192(double x) const;
3748 const double getIntegrand_sigmaWithISR_mu196(double x) const;
3749 const double getIntegrand_sigmaWithISR_mu200(double x) const;
3750 const double getIntegrand_sigmaWithISR_mu202(double x) const;
3751 const double getIntegrand_sigmaWithISR_mu205(double x) const;
3752 const double getIntegrand_sigmaWithISR_mu207(double x) const;
3753
3754
3755 const double getIntegrand_sigmaWithISR_tau130(double x) const;
3756 const double getIntegrand_sigmaWithISR_tau136(double x) const;
3757 const double getIntegrand_sigmaWithISR_tau161(double x) const;
3758 const double getIntegrand_sigmaWithISR_tau172(double x) const;
3759 const double getIntegrand_sigmaWithISR_tau183(double x) const;
3760 const double getIntegrand_sigmaWithISR_tau189(double x) const;
3761 const double getIntegrand_sigmaWithISR_tau192(double x) const;
3762 const double getIntegrand_sigmaWithISR_tau196(double x) const;
3763 const double getIntegrand_sigmaWithISR_tau200(double x) const;
3764 const double getIntegrand_sigmaWithISR_tau202(double x) const;
3765 const double getIntegrand_sigmaWithISR_tau205(double x) const;
3766 const double getIntegrand_sigmaWithISR_tau207(double x) const;
3767
3768
3769 const double Integrand_sigmaWithISR_q(double x, const QCD::quark q_flavor, const double s) const;
3770
3771 const double getIntegrand_sigmaWithISR_up130(double x) const;
3772 const double getIntegrand_sigmaWithISR_up133(double x) const;
3773 const double getIntegrand_sigmaWithISR_up136(double x) const;
3774 const double getIntegrand_sigmaWithISR_up161(double x) const;
3775 const double getIntegrand_sigmaWithISR_up167(double x) const;
3776 const double getIntegrand_sigmaWithISR_up172(double x) const;
3777 const double getIntegrand_sigmaWithISR_up183(double x) const;
3778 const double getIntegrand_sigmaWithISR_up189(double x) const;
3779 const double getIntegrand_sigmaWithISR_up192(double x) const;
3780 const double getIntegrand_sigmaWithISR_up196(double x) const;
3781 const double getIntegrand_sigmaWithISR_up200(double x) const;
3782 const double getIntegrand_sigmaWithISR_up202(double x) const;
3783 const double getIntegrand_sigmaWithISR_up205(double x) const;
3784 const double getIntegrand_sigmaWithISR_up207(double x) const;
3785
3786 const double getIntegrand_sigmaWithISR_down130(double x) const;
3787 const double getIntegrand_sigmaWithISR_down133(double x) const;
3788 const double getIntegrand_sigmaWithISR_down136(double x) const;
3789 const double getIntegrand_sigmaWithISR_down161(double x) const;
3790 const double getIntegrand_sigmaWithISR_down167(double x) const;
3791 const double getIntegrand_sigmaWithISR_down172(double x) const;
3792 const double getIntegrand_sigmaWithISR_down183(double x) const;
3793 const double getIntegrand_sigmaWithISR_down189(double x) const;
3794 const double getIntegrand_sigmaWithISR_down192(double x) const;
3795 const double getIntegrand_sigmaWithISR_down196(double x) const;
3796 const double getIntegrand_sigmaWithISR_down200(double x) const;
3797 const double getIntegrand_sigmaWithISR_down202(double x) const;
3798 const double getIntegrand_sigmaWithISR_down205(double x) const;
3799 const double getIntegrand_sigmaWithISR_down207(double x) const;
3800
3801 const double getIntegrand_sigmaWithISR_charm130(double x) const;
3802 const double getIntegrand_sigmaWithISR_charm133(double x) const;
3803 const double getIntegrand_sigmaWithISR_charm136(double x) const;
3804 const double getIntegrand_sigmaWithISR_charm161(double x) const;
3805 const double getIntegrand_sigmaWithISR_charm167(double x) const;
3806 const double getIntegrand_sigmaWithISR_charm172(double x) const;
3807 const double getIntegrand_sigmaWithISR_charm183(double x) const;
3808 const double getIntegrand_sigmaWithISR_charm189(double x) const;
3809 const double getIntegrand_sigmaWithISR_charm192(double x) const;
3810 const double getIntegrand_sigmaWithISR_charm196(double x) const;
3811 const double getIntegrand_sigmaWithISR_charm200(double x) const;
3812 const double getIntegrand_sigmaWithISR_charm202(double x) const;
3813 const double getIntegrand_sigmaWithISR_charm205(double x) const;
3814 const double getIntegrand_sigmaWithISR_charm207(double x) const;
3815
3816 const double getIntegrand_sigmaWithISR_strange130(double x) const;
3817 const double getIntegrand_sigmaWithISR_strange133(double x) const;
3818 const double getIntegrand_sigmaWithISR_strange136(double x) const;
3819 const double getIntegrand_sigmaWithISR_strange161(double x) const;
3820 const double getIntegrand_sigmaWithISR_strange167(double x) const;
3821 const double getIntegrand_sigmaWithISR_strange172(double x) const;
3822 const double getIntegrand_sigmaWithISR_strange183(double x) const;
3823 const double getIntegrand_sigmaWithISR_strange189(double x) const;
3824 const double getIntegrand_sigmaWithISR_strange192(double x) const;
3825 const double getIntegrand_sigmaWithISR_strange196(double x) const;
3826 const double getIntegrand_sigmaWithISR_strange200(double x) const;
3827 const double getIntegrand_sigmaWithISR_strange202(double x) const;
3828 const double getIntegrand_sigmaWithISR_strange205(double x) const;
3829 const double getIntegrand_sigmaWithISR_strange207(double x) const;
3830
3831 const double getIntegrand_sigmaWithISR_bottom130(double x) const;
3832 const double getIntegrand_sigmaWithISR_bottom133(double x) const;
3833 const double getIntegrand_sigmaWithISR_bottom136(double x) const;
3834 const double getIntegrand_sigmaWithISR_bottom161(double x) const;
3835 const double getIntegrand_sigmaWithISR_bottom167(double x) const;
3836 const double getIntegrand_sigmaWithISR_bottom172(double x) const;
3837 const double getIntegrand_sigmaWithISR_bottom183(double x) const;
3838 const double getIntegrand_sigmaWithISR_bottom189(double x) const;
3839 const double getIntegrand_sigmaWithISR_bottom192(double x) const;
3840 const double getIntegrand_sigmaWithISR_bottom196(double x) const;
3841 const double getIntegrand_sigmaWithISR_bottom200(double x) const;
3842 const double getIntegrand_sigmaWithISR_bottom202(double x) const;
3843 const double getIntegrand_sigmaWithISR_bottom205(double x) const;
3844 const double getIntegrand_sigmaWithISR_bottom207(double x) const;
3845
3846 const double Integrand_dsigmaBox_l(double cosTheta, const QCD::lepton l_flavor, const double s) const;
3847
3848 const double getIntegrand_dsigmaBox_mu130(double x) const;
3849 const double getIntegrand_dsigmaBox_mu133(double x) const;
3850 const double getIntegrand_dsigmaBox_mu136(double x) const;
3851 const double getIntegrand_dsigmaBox_mu161(double x) const;
3852 const double getIntegrand_dsigmaBox_mu167(double x) const;
3853 const double getIntegrand_dsigmaBox_mu172(double x) const;
3854 const double getIntegrand_dsigmaBox_mu183(double x) const;
3855 const double getIntegrand_dsigmaBox_mu189(double x) const;
3856 const double getIntegrand_dsigmaBox_mu192(double x) const;
3857 const double getIntegrand_dsigmaBox_mu196(double x) const;
3858 const double getIntegrand_dsigmaBox_mu200(double x) const;
3859 const double getIntegrand_dsigmaBox_mu202(double x) const;
3860 const double getIntegrand_dsigmaBox_mu205(double x) const;
3861 const double getIntegrand_dsigmaBox_mu207(double x) const;
3862
3863 const double getIntegrand_dsigmaBox_tau130(double x) const;
3864 const double getIntegrand_dsigmaBox_tau133(double x) const;
3865 const double getIntegrand_dsigmaBox_tau136(double x) const;
3866 const double getIntegrand_dsigmaBox_tau161(double x) const;
3867 const double getIntegrand_dsigmaBox_tau167(double x) const;
3868 const double getIntegrand_dsigmaBox_tau172(double x) const;
3869 const double getIntegrand_dsigmaBox_tau183(double x) const;
3870 const double getIntegrand_dsigmaBox_tau189(double x) const;
3871 const double getIntegrand_dsigmaBox_tau192(double x) const;
3872 const double getIntegrand_dsigmaBox_tau196(double x) const;
3873 const double getIntegrand_dsigmaBox_tau200(double x) const;
3874 const double getIntegrand_dsigmaBox_tau202(double x) const;
3875 const double getIntegrand_dsigmaBox_tau205(double x) const;
3876 const double getIntegrand_dsigmaBox_tau207(double x) const;
3877
3878
3879 const double Integrand_dsigmaBox_q(double cosTheta, const QCD::quark q_flavor, const double s) const;
3880
3881 const double getIntegrand_dsigmaBox_up130(double x) const;
3882 const double getIntegrand_dsigmaBox_up133(double x) const;
3883 const double getIntegrand_dsigmaBox_up136(double x) const;
3884 const double getIntegrand_dsigmaBox_up161(double x) const;
3885 const double getIntegrand_dsigmaBox_up167(double x) const;
3886 const double getIntegrand_dsigmaBox_up172(double x) const;
3887 const double getIntegrand_dsigmaBox_up183(double x) const;
3888 const double getIntegrand_dsigmaBox_up189(double x) const;
3889 const double getIntegrand_dsigmaBox_up192(double x) const;
3890 const double getIntegrand_dsigmaBox_up196(double x) const;
3891 const double getIntegrand_dsigmaBox_up200(double x) const;
3892 const double getIntegrand_dsigmaBox_up202(double x) const;
3893 const double getIntegrand_dsigmaBox_up205(double x) const;
3894 const double getIntegrand_dsigmaBox_up207(double x) const;
3895
3896 const double getIntegrand_dsigmaBox_down130(double x) const;
3897 const double getIntegrand_dsigmaBox_down133(double x) const;
3898 const double getIntegrand_dsigmaBox_down136(double x) const;
3899 const double getIntegrand_dsigmaBox_down161(double x) const;
3900 const double getIntegrand_dsigmaBox_down167(double x) const;
3901 const double getIntegrand_dsigmaBox_down172(double x) const;
3902 const double getIntegrand_dsigmaBox_down183(double x) const;
3903 const double getIntegrand_dsigmaBox_down189(double x) const;
3904 const double getIntegrand_dsigmaBox_down192(double x) const;
3905 const double getIntegrand_dsigmaBox_down196(double x) const;
3906 const double getIntegrand_dsigmaBox_down200(double x) const;
3907 const double getIntegrand_dsigmaBox_down202(double x) const;
3908 const double getIntegrand_dsigmaBox_down205(double x) const;
3909 const double getIntegrand_dsigmaBox_down207(double x) const;
3910
3911 const double getIntegrand_dsigmaBox_charm130(double x) const;
3912 const double getIntegrand_dsigmaBox_charm133(double x) const;
3913 const double getIntegrand_dsigmaBox_charm136(double x) const;
3914 const double getIntegrand_dsigmaBox_charm161(double x) const;
3915 const double getIntegrand_dsigmaBox_charm167(double x) const;
3916 const double getIntegrand_dsigmaBox_charm172(double x) const;
3917 const double getIntegrand_dsigmaBox_charm183(double x) const;
3918 const double getIntegrand_dsigmaBox_charm189(double x) const;
3919 const double getIntegrand_dsigmaBox_charm192(double x) const;
3920 const double getIntegrand_dsigmaBox_charm196(double x) const;
3921 const double getIntegrand_dsigmaBox_charm200(double x) const;
3922 const double getIntegrand_dsigmaBox_charm202(double x) const;
3923 const double getIntegrand_dsigmaBox_charm205(double x) const;
3924 const double getIntegrand_dsigmaBox_charm207(double x) const;
3925
3926 const double getIntegrand_dsigmaBox_strange130(double x) const;
3927 const double getIntegrand_dsigmaBox_strange133(double x) const;
3928 const double getIntegrand_dsigmaBox_strange136(double x) const;
3929 const double getIntegrand_dsigmaBox_strange161(double x) const;
3930 const double getIntegrand_dsigmaBox_strange167(double x) const;
3931 const double getIntegrand_dsigmaBox_strange172(double x) const;
3932 const double getIntegrand_dsigmaBox_strange183(double x) const;
3933 const double getIntegrand_dsigmaBox_strange189(double x) const;
3934 const double getIntegrand_dsigmaBox_strange192(double x) const;
3935 const double getIntegrand_dsigmaBox_strange196(double x) const;
3936 const double getIntegrand_dsigmaBox_strange200(double x) const;
3937 const double getIntegrand_dsigmaBox_strange202(double x) const;
3938 const double getIntegrand_dsigmaBox_strange205(double x) const;
3939 const double getIntegrand_dsigmaBox_strange207(double x) const;
3940
3941 const double getIntegrand_dsigmaBox_bottom130(double x) const;
3942 const double getIntegrand_dsigmaBox_bottom133(double x) const;
3943 const double getIntegrand_dsigmaBox_bottom136(double x) const;
3944 const double getIntegrand_dsigmaBox_bottom161(double x) const;
3945 const double getIntegrand_dsigmaBox_bottom167(double x) const;
3946 const double getIntegrand_dsigmaBox_bottom172(double x) const;
3947 const double getIntegrand_dsigmaBox_bottom183(double x) const;
3948 const double getIntegrand_dsigmaBox_bottom189(double x) const;
3949 const double getIntegrand_dsigmaBox_bottom192(double x) const;
3950 const double getIntegrand_dsigmaBox_bottom196(double x) const;
3951 const double getIntegrand_dsigmaBox_bottom200(double x) const;
3952 const double getIntegrand_dsigmaBox_bottom202(double x) const;
3953 const double getIntegrand_dsigmaBox_bottom205(double x) const;
3954 const double getIntegrand_dsigmaBox_bottom207(double x) const;
3955
3956
3957 const double Integrand_AFBnumeratorWithISR_l(double x, const QCD::lepton l_flavor, const double s) const;
3958
3959 const double getIntegrand_AFBnumeratorWithISR_mu130(double x) const;
3960 const double getIntegrand_AFBnumeratorWithISR_mu136(double x) const;
3961 const double getIntegrand_AFBnumeratorWithISR_mu161(double x) const;
3962 const double getIntegrand_AFBnumeratorWithISR_mu172(double x) const;
3963 const double getIntegrand_AFBnumeratorWithISR_mu183(double x) const;
3964 const double getIntegrand_AFBnumeratorWithISR_mu189(double x) const;
3965 const double getIntegrand_AFBnumeratorWithISR_mu192(double x) const;
3966 const double getIntegrand_AFBnumeratorWithISR_mu196(double x) const;
3967 const double getIntegrand_AFBnumeratorWithISR_mu200(double x) const;
3968 const double getIntegrand_AFBnumeratorWithISR_mu202(double x) const;
3969 const double getIntegrand_AFBnumeratorWithISR_mu205(double x) const;
3970 const double getIntegrand_AFBnumeratorWithISR_mu207(double x) const;
3971
3972 const double getIntegrand_AFBnumeratorWithISR_tau130(double x) const;
3973 const double getIntegrand_AFBnumeratorWithISR_tau136(double x) const;
3974 const double getIntegrand_AFBnumeratorWithISR_tau161(double x) const;
3975 const double getIntegrand_AFBnumeratorWithISR_tau172(double x) const;
3976 const double getIntegrand_AFBnumeratorWithISR_tau183(double x) const;
3977 const double getIntegrand_AFBnumeratorWithISR_tau189(double x) const;
3978 const double getIntegrand_AFBnumeratorWithISR_tau192(double x) const;
3979 const double getIntegrand_AFBnumeratorWithISR_tau196(double x) const;
3980 const double getIntegrand_AFBnumeratorWithISR_tau200(double x) const;
3981 const double getIntegrand_AFBnumeratorWithISR_tau202(double x) const;
3982 const double getIntegrand_AFBnumeratorWithISR_tau205(double x) const;
3983 const double getIntegrand_AFBnumeratorWithISR_tau207(double x) const;
3984
3985 const double Integrand_AFBnumeratorWithISR_q(double x, const QCD::quark q_flavor, const double s) const;
3986
3987 const double getIntegrand_AFBnumeratorWithISR_charm133(double x) const;
3988 const double getIntegrand_AFBnumeratorWithISR_charm167(double x) const;
3989 const double getIntegrand_AFBnumeratorWithISR_charm172(double x) const;
3990 const double getIntegrand_AFBnumeratorWithISR_charm183(double x) const;
3991 const double getIntegrand_AFBnumeratorWithISR_charm189(double x) const;
3992 const double getIntegrand_AFBnumeratorWithISR_charm192(double x) const;
3993 const double getIntegrand_AFBnumeratorWithISR_charm196(double x) const;
3994 const double getIntegrand_AFBnumeratorWithISR_charm200(double x) const;
3995 const double getIntegrand_AFBnumeratorWithISR_charm202(double x) const;
3996 const double getIntegrand_AFBnumeratorWithISR_charm205(double x) const;
3997 const double getIntegrand_AFBnumeratorWithISR_charm207(double x) const;
3998
3999 const double getIntegrand_AFBnumeratorWithISR_bottom133(double x) const;
4000 const double getIntegrand_AFBnumeratorWithISR_bottom167(double x) const;
4001 const double getIntegrand_AFBnumeratorWithISR_bottom172(double x) const;
4002 const double getIntegrand_AFBnumeratorWithISR_bottom183(double x) const;
4003 const double getIntegrand_AFBnumeratorWithISR_bottom189(double x) const;
4004 const double getIntegrand_AFBnumeratorWithISR_bottom192(double x) const;
4005 const double getIntegrand_AFBnumeratorWithISR_bottom196(double x) const;
4006 const double getIntegrand_AFBnumeratorWithISR_bottom200(double x) const;
4007 const double getIntegrand_AFBnumeratorWithISR_bottom202(double x) const;
4008 const double getIntegrand_AFBnumeratorWithISR_bottom205(double x) const;
4009 const double getIntegrand_AFBnumeratorWithISR_bottom207(double x) const;
4010
4011 /* END: REMOVE FROM THE PACKAGE */
4013private:
4023 /* BEGIN: REMOVE FROM THE PACKAGE */
4025 /* END: REMOVE FROM THE PACKAGE */
4026
4029 std::string FlagMw;
4030 std::string FlagRhoZ;
4031 std::string FlagKappaZ;
4034
4037
4038 mutable bool SMSuccess;
4039
4041 // Caches for EWPO
4042
4045 mutable double DeltaAlphaLepton_cache;
4046 mutable double DeltaAlpha_cache;
4047 mutable double Mw_cache;
4048 mutable double GammaW_cache;
4049 mutable gslpp::complex rhoZ_f_cache[12];
4050 mutable gslpp::complex kappaZ_f_cache[12];
4053 mutable bool useMw_cache;
4054 mutable bool useGammaW_cache;
4055 mutable bool useRhoZ_f_cache[12];
4056 mutable bool useKappaZ_f_cache[12];
4057
4058/* BEGIN: REMOVE FROM THE PACKAGE */
4059 // caches for the SM prediction of LEP2 Obs
4060// mutable double SMparams_cache[NumSMParamsForEWPO+3];
4061 mutable double SMresult_cache;
4062// mutable bool flag_cache[NUMofLEP2RCs];
4063// mutable double ml_cache, mq_cache, mqForHad_cache[6];
4064
4065 mutable double average;
4066 mutable double error;
4067 mutable gsl_function f_GSL;
4068 gsl_integration_workspace * w_GSL1;
4069/* END: REMOVE FROM THE PACKAGE */
4070
4072
4073 const double AlsE(double mu, orders order, bool Nf_thr) const;
4074 const double AlsEByOrder(double mu, orders order, bool Nf_thr) const;
4075
4076 const double AlsEWithInit(double mu, double alsi, double mu_i, const int nf_i, orders order) const;
4077 const double AleWithInit(double mu, double alsi, double mu_i, orders order) const;
4078 static const int CacheSize = 5;
4079 mutable double als_cache[11][CacheSize];
4080 mutable double ale_cache[10][CacheSize];
4082
4083
4084};
4085
4086#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:36
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 .
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
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.
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 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
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