#include "METOOLS/Explicit/NLO_Counter_Terms.H" #include "ATOOLS/Org/Message.H" #include "ATOOLS/Phys/Flavour.H" #include "MODEL/Main/Running_AlphaS.H" using namespace METOOLS; using namespace ATOOLS; double METOOLS::Hab(const Flavour &a,const Flavour &b) { if (a.IsQuark()) { if (b.IsQuark()) return a==b?4.0/3.0*3.0/2.0:0.0; return 0.0; } else { if (b.IsQuark()) return 0.0; return 11.0/6.0*3.0-2.0/3.0*0.5*(Flavour(kf_jet).Size()/2); } } double METOOLS::FPab(const Flavour &a,const Flavour &b,const double &z) { if (a.IsQuark()) { if (b.IsQuark()) return a==b?-4.0/3.0*(1.0+z):0.0; return 4.0/3.0*(1.0+sqr(1.0-z))/z; } else { if (b.IsQuark()) return 1.0/2.0*(z*z+sqr(1.0-z)); return 3.0*2.0*((1.0-z)/z-1.0+z*(1.0-z)); } } double METOOLS::SPab(const Flavour &a,const Flavour &b,const double &z) { if (a.IsQuark()) { if (b.IsQuark()) return a==b?4.0/3.0*2.0/(1.0-z):0.0; return 0.0; } else { if (b.IsQuark()) return 0.0; return 3.0*2.0/(1.0-z); } } double METOOLS::IPab(const Flavour &a,const Flavour &b,const double &x) { if (a.IsQuark()) { if (b.IsQuark() && a==b) return 4.0/3.0*2.0*log(1.0/(1.0-x)); return 0.0; } else { if (b.IsQuark()) return 0.0; return 3.0*2.0*log(1.0/(1.0-x)); } } double METOOLS::AlphaSCounterTerm(const double& muR2,const double& muR2ref, const double& as, MODEL::One_Running_AlphaS * oras, const size_t oqcd) { DEBUG_FUNC(muR2<<" -> "< thrs(oras->Thresholds(muR2,muR2ref)); msg_Debugging()<<"Flavour thresholds in range ["<muR2ref)?muR2:muR2ref); thrs.insert(thrs.begin(),(muR2>muR2ref)?muR2ref:muR2); double betalog(0.); msg_Debugging()<<"\\sum_{\\mu_{th}} \\beta_0(nf(\\mu_i)) " <<"log(\\mu_{i+1}^2/\\mu_i^2) = "<Beta0(thrs[i+1]) <<" * "<Beta0((thrs[i+1]+thrs[i])/2.)*log(thrs[i+1]/thrs[i]); } if (muR2>muR2ref) betalog*=-1.; msg_Debugging()<<" = "<Contains(fl)) return 0.; pdf->Calculate(x,tref); double fb(pdf->GetXPDF(fl)/x); if (IsZero(fb,th)) { msg_Tracking()<Contains(jet[j])) continue; double Pf(METOOLS::FPab(jet[j],fl,z)); double Ps(METOOLS::SPab(jet[j],fl,z)); if (Pf+Ps==0.0) continue; double Pi(METOOLS::IPab(jet[j],fl,x)); double H(METOOLS::Hab(jet[j],fl)); pdf->Calculate(x/z,tref); double fa(pdf->GetXPDF(jet[j])/x); pdf->Calculate(x,tref); double fc(pdf->GetXPDF(jet[j])/x); msg_Debugging()<<" P_{"< "<