Associated t \bar{t} H production at the LHC: theoretical predictions at NLO+NNLL accuracy
Anna Kulesza, Leszek Motyka, Tomasz Stebel, Vincent Theeuwes
MMS-TP-17-06
Associated t ¯ tH production at the LHC: theoreticalpredictions at NLO+NNLL accuracy Anna Kulesza a, , Leszek Motyka b, , Tomasz Stebel b,c, and Vincent Theeuwes d, a Institute for Theoretical Physics, WWU Münster, D-48149 Münster, Germany b Institute of Physics, Jagellonian University, S.Łojasiewicza 11, 30-348 Kraków, Poland c Institute of Nuclear Physics PAN, Radzikowskiego 152, 31-342 Kraków, Poland d Department of Physics, SUNY Buffalo, 261 Fronczak Hall, Buffalo, NY 14260-1500,USA
Abstract
We perform threshold resummation of soft gluon corrections to the total crosssection and the invariant mass distribution for the process pp → t ¯ tH . The resum-mation is carried out at next-to-next-to-leading-logarithmic (NNLL) accuracy usingthe direct QCD Mellin space technique in the three-particle invariant mass kinemat-ics. After presenting analytical expressions we discuss the impact of resummation onthe numerical predictions for the associated Higgs boson production with top quarksat the LHC. We find that NLO+NNLL resummation leads to predictions for whichthe central values are remarkably stable with respect to scale variation and for whichtheoretical uncertainties are reduced in comparison to NLO predictions. [email protected] [email protected] [email protected] [email protected] a r X i v : . [ h e p - ph ] O c t ontents → processes with two massive coloured particles in the final state 33 Numerical results for the pp → t ¯ tH process at NLO+NNLL accuracy 84 Summary 15 The measurement of Higgs boson production rates in the pp → t ¯ tH process is of the centralimportance to the LHC research program. The process has been intensively searched for inRun 1 [1–5] and its measurement is among the highest priorities of the LHC Run 2 physicsprogramme [6–8]. The associated production process offers a direct way to probe thestrength of the top–Higgs Yukawa coupling without making any assumptions regarding itsnature. As the top–Higgs Yukawa coupling is especially sensitive to the underlying physics, t ¯ tH production provides a vital test of the Standard Model (SM) and possibly a means toprobe the beyond the SM physics indirectly. It is thus highly important that precise andreliable theoretical predictions are available for this process.For these reasons, a large amount of effort has been invested to improve theoretical de-scription of the t ¯ tH production. While the next-to-leading-order (NLO) QCD, i.e. O ( α α ) predictions were obtained some time ago [9, 10], they have been newly recalculated andmatched to parton showers in [11–14]. In the last years, the mixed QCD-weak correc-tions [15] and QCD-EW corrections [16, 17] of O ( α α ) are also available. Furthermore,the NLO EW and QCD corrections to the hadronic t ¯ tH production with off-shell top andantitop quarks have been recently obtained [18, 19]. For the most part, the NLO QCDcorrections are ∼ % at the Run 2 LHC energies, whereas the size of the (electro)weak cor-rection is more than ten times smaller. The scale uncertainty of the NLO QCD correctionsis estimated to be ∼ % [9, 10, 20].In general, if for a given process one expects that a significant part of higher ordercorrections originates from emission of soft and/or collinear gluons, it is possible to improvethe accuracy of theoretical predictions by employing methods of resummation. Relyingon principles of factorization between various dynamical modes, they allow an all-ordercalculation of dominant logarithmic corrections originating from a certain kinematical limit.Supplementing fixed-order results with resummation leads not only to change in the valueof the cross section but also offers a better control over the theoretical error, in particulardue to cancellations of factorization scale dependence between parton distribution functions(pdfs) and the partonic cross sections. The universality of resummation concepts warrantstheir applications to scattering processes with arbitrary many partons in the final state [21,22], thus also to a class → processes and in particular the associated t ¯ tH productionat the LHC.The first step in this direction was performed by us in [23], where we presented thefirst calculation of the resummed total cross section for the t ¯ tH production at the next-to-leading-logarithmic (NLL) accuracy. The calculation relied on application of the traditionalMellin-space resummation formalism in the absolute threshold limit, i.e. in the limit of the2artonic energy √ ˆ s approaching the production threshold M = 2 m t + m H , ˆ s → M , where m t is the top quark mass and m H is the Higgs boson mass. In [23], we have achieved anall-order improvement of the theoretical predictions by taking into account a well-definedsubclass of higher order corrections. However, due to the suppression of the available 3-particle phase-space in the absolute production threshold limit, it is to be expected that thenumerical impact of formally large logarithmic corrections resummed in these kinematicswill be somewhat diminished and that contributions prevailing numerically might comefrom regions further away from the absolute threshold scale M .Subsequently we have performed [24] resummation of NLL corrections arising in thelimit of √ ˆ s approaching the invariant mass threshold Q with Q = ( p t + p ¯ t + p H ) . Wehave considered cross sections differential in the invariant mass Q , as well as the total crosssections obtained after integration over Q . For a → process, this type of resummationis often referred to as threshold resummation in the pair-invariant mass (PIM) kinemat-ics. Threshold resummation can be also performed in the framework of the soft-collineareffective theory (SCET). The first application of this technique to a → process, morespecifically to the process pp → t ¯ tW ± , was presented in [25]. The SCET framework wasalso used to obtain an approximation of the next-to-next-to-leading order (NNLO) t ¯ tH cross section and distributions [26], following from the expansion of the NNLL resum-mation formula. Recently, NNLL results for t ¯ tW [27], t ¯ tH [28] and t ¯ tZ [29] associatedproduction processes appeared, based on expressing the SCET formulas in Mellin space.In this paper, we continue the work presented in [24] and perform threshold resumma-tion in the invariant mass limit at the NNLL accuracy using the direct QCD [30] Mellin-space approach. Compared to NLL calculations, the anomalous dimensions governingresummation need to be implemented with accuracy higher by one order. In contrast tothe absolute threshold limit considered in [23], the soft anomalous dimension is a matrixin the colour space containing non-zero off-diagonal elements, thus requiring an implemen-tation of the diagonalization procedure. We then match our NNLL cross section with thefixed-order cross section at NLO. The invariant mass kinematics also offers an opportunityto perform resummation for the differential distributions in Q .The paper is structured as follows. In Section 2 we review threshold resummation inMellin space, stressing the difference between the resummation in the invariant mass andthe absolute threshold limits. The numerical results and their discussion is presented inSection 3, where we also compare our results to those in [28]. We summarize our mostimportant findings in Section 4. → processes with two massive coloured particles inthe final state The resummation of soft gluon corrections to the differential cross section dσ pp → t ¯ tH /dQ is performed in Mellin space, where the Mellin moments are taken w.r.t. the variable ρ = Q /S . At the partonic level, the Mellin moments for the process ij → klB , where i, j denote massless coloured partons, k, l two massive quarks and B a massive colour-singletparticle, are given by d ˜ˆ σ ij → klB dQ ( N, Q , { m } , µ , µ ) = (cid:90) d ˆ ρ ˆ ρ N − d ˆ σ ij → klB dQ ( ˆ ρ, Q , { m } , µ , µ ) , (1)3ith ˆ ρ = Q / ˆ s and { m } denoting all masses entering the calculations.Taking the Mellin transform allows one to systematically treat the logarithmic termsof the form α n s [log m (1 − z ) / (1 − z )] + , with m ≤ n − and z = Q / ˆ s , appearing in theperturbative expansion of the partonic cross section to all orders in α s . In Mellin spacethese logarithms turn into logarithms of the variable N , and the threshold limit z → corresponds to the limit N → ∞ .The resummed cross section in the N -space has the form [31, 32] d ˜ˆ σ (res) ij → klB dQ ( N, Q , { m } , µ , µ ) = (2) = Tr (cid:2) H ij → klB ( Q , { m } , µ , µ ) S ij → klB ( N + 1 , Q , { m } , µ ) (cid:3) × ∆ i ( N + 1 , Q , µ , µ )∆ j ( N + 1 , Q , µ , µ ) , where the trace is taken over colour space. The appearance of colour dependence in Eq. (2)is inherently related to the fact that soft radiation is coherently sensitive to the colourstructure of the hard process from which it is emitted. The matrix H ij → klB indicatesthe hard-scattering contributions, absorbing the off-shell effects, projected onto the chosencolour basis. The colour matrix S ij → klB represents the soft wide-angle emission. Thefunctions ∆ i and ∆ j sum the logarithmic contributions due to (soft-)collinear radiationfrom the incoming partons. The radiative factors are thus universal for a specific initialstate parton, i.e. they depend neither on the underlying colour structure nor the process.At LO the t ¯ tH production receives contributions from the q ¯ q and gg channels. Weanalyze the colour structure of the underlying processes in the s -channel colour bases, { c qI } and { c gI } , with c q = δ α i α j δ α k α l , c q = T aα i α j T aα k α l ,c g = δ a i a j δ α k α l , c g = T bα l α k d ba i a j , c g = iT bα l α k f ba i a j . The hard function H ij → klB carries no dependence on N and is given by the perturbativeexpansion H ij → klB = H (0) ij → klB + α s π H (1) ij → klB + . . . (3)In order to perform resummation at NLL accuracy one needs to know H (0) ij → klB , whereasNNLL accuracy requires the knowledge of the H (1) ij → klB coefficient.The soft function, on the other hand, resums logarithms of N at the rate of one powerof the logarithm per power of the strong coupling. These single logarithms due to thesoft emission can be confronted with double logarithms due to soft and collinear emissionsresummed by the jet factors ∆ i and ∆ j . As a dimensionless function, S ij → klB dependsonly on the ratio of the scales. At the same time, the dependence on N enters only via Q/N [33], making S ij → klB dependent on the ratio Q/ ( N µ R ) . The soft function is givenby a solution of the renormalization group equation [31, 34] and has the form S ij → klB ( N, Q , { m } , µ ) = ¯U ij → klB ( N, Q , { m } , µ ) ˜S ij → klB ( α s ( Q / ¯ N )) × U ij → klB ( N, Q , { m } , µ ) , (4) In fact, the soft function S ij → klB as well as the radiative factors ∆ i , ∆ j are dimensionless functionsof the ratios of the scales and the coupling constant at the renormalization scale. The current notationindicating dependence on the scales is introduced for brevity. ˜S ij → klB plays a role of a boundary condition and is obtained by taking S ij → klB at Q / ( ¯ N µ ) = 1 with ¯ N = N e γ E and γ E denoting the Euler constant. It is a purely eikonalfunction [31, 34, 35] and can be calculated perturbatively ˜S ij → klB = ˜S (0) ij → klB + α s π ˜S (1) ij → klB + . . . (5)At the lowest order the colour matrix is given by: (cid:16) ˜S (0) ij → kl (cid:17) IJ = Tr (cid:104) c † I c J (cid:105) . (6)Similarly to the hard function, knowledge of S (0) ij → klB is required in order to perform re-summation at NLL accuracy and a result for S (1) ij → klB at the NNLL accuracy. The scaleof α s in ˜S ij → kl , equal to Q / ¯ N , results in an order α log ¯ N term if we expand ˜S ij → kl in α s ( µ ) .The soft function evolution matrices ¯U ij → klB , U ij → klB contain logarithmic enhance-ments due to soft wide-angle emissions [31, 36] ¯U ij → klB (cid:0) N, Q , { m } , µ (cid:1) = ¯P exp (cid:34)(cid:90) Q/ ¯ Nµ R dqq Γ † ij → klB (cid:0) α s (cid:0) q (cid:1)(cid:1)(cid:35) , (7) U ij → klB (cid:0) N, Q , { m } , µ (cid:1) = P exp (cid:34)(cid:90) Q/ ¯ Nµ R dqq Γ ij → klB (cid:0) α s (cid:0) q (cid:1)(cid:1)(cid:35) , where P and ¯P denote the path- and reverse path-ordering in the variable q , respectively.The soft anomalous dimension Γ ij → klB is a perturbative function in α s : Γ ij → klB ( α s ) = (cid:20)(cid:16) α s π (cid:17) Γ (1) ij → klB + (cid:16) α s π (cid:17) Γ (2) ij → klB + . . . (cid:21) . (8)In order to perform resummation at NLL accuracy we need to know Γ (1) ij → klB , whereasNNLL accuracy requires including Γ (2) ij → klB . The one-loop soft anomalous dimension forthe process ij → klB with k, l being heavy quarks can be found e.g. in [23]. The two-loopcontributions to the soft anomalous dimension were calculated in [37, 38] . In the triple-invariant mass (TIM) kinematics, the soft anomalous dimension matrix in general containsoff-diagonal terms, thus complicating the evaluation of the resummed cross section. AtNNLL additional difficulty arises because of non-commutativity of Γ (1) ij → klB and Γ (2) ij → klB matrices.We make use of the method of [31] in order to diagonalize the one-loop soft anomalousdimension matrix. Denoting the diagonalization matrix for Γ (1) ij → klB by R we have Γ (1) R = R − Γ (1) ij → klB R , (9)where the diagonalized matrix is given by eigenvalues λ (1) I of Γ (1) ij → klB Γ (1) R,IJ = λ (1) I δ IJ , For simplicity, the argument dependence of the soft anomalous dimension on the mass scales is sup-pressed in Eq. (7). Note that while using the radiative factors as given in [39, 40], we need to subtract the collinear softradiation already included in ∆ i , ∆ j from the eikonal cross section used to calculate the soft function. Γ (1) R = (cid:104) −→ λ (1) (cid:105) D with −→ λ (1) = (cid:110) λ (1)1 , . . . , λ (1) D (cid:111) . The othermatrices are transformed as: Γ (2) R = R − Γ (2) ij → klB R , H R = R − H ij → klB (cid:0) R − (cid:1) † , (10) ˜S R = R † ˜S ij → klB R . At NLL accuracy, by changing the colour basis to the one in which Γ (1) ij → klB is diagonal,the path ordered exponentials in Eq. (4) reduce to sum over simple exponentials. At NNLLaccuracy, to recast the path ordered exponential of the soft anomalous dimension matrixin a form containing simple exponential functions, we make use of a technique detailed ine.g. [41, 42] resulting in U R ( N, Q , { m } , µ ) = (cid:32) + α s (cid:0) Q / ¯ N (cid:1) π K (cid:33) (cid:32) α s (cid:0) µ (cid:1) α s (cid:0) Q / ¯ N (cid:1) (cid:33) −→ λ (1)2 πb D (cid:32) − α s (cid:0) µ (cid:1) π K (cid:33) , (11)with the subscript D indicating a diagonal matrix. The matrix K is given by K IJ = δ IJ λ (1) I b b − (cid:16) Γ (2) R (cid:17) IJ πb + λ (1) I − λ (1) J , (12)where b and b are the first two coefficients of expansion β QCD in α s : b = 11 C A − n f T R π , (13) b = 17 C − n f T R (10 C A + 6 C F )24 π . (14)In our calculation we set n f = 5 .In the diagonal basis of the one-loop soft anomalous dimension, up to NNLL accuracyEq. (2) can be written as d ˜ˆ σ (NNLL) ij → klB dQ ( N, Q , { m } , µ ) = Tr (cid:2) H R ( Q , { m } , µ , µ ) ¯U R ( N + 1 , Q , { m } , Q ) × ˜S R ( N + 1 , Q , { m } ) U R ( N + 1 , Q , { m } , Q ) (cid:105) × ∆ i ( N + 1 , Q , µ , µ )∆ j ( N + 1 , Q , µ , µ ) . (15)In the above equation, the H R and ˜S R are hard and soft function matrices projectedonto R colour basis. They are calculated at the NLO accuracy, i.e. including the O ( α s ) terms in Eqs. (3) and (5). The LO hard matrix is derived from the Born cross section.The NLO hard matrix contains non-logarithmic contributions which are independent of N . They consist of virtual loop contributions, real terms of collinear origin and the contri-butions from the evolution matrices U R and ¯U R , corresponding to evolution between µ R and Q , expanded up to O ( α s ) . The colour-decomposed virtual corrections are extractedfrom the calculations of the NLO cross section in the PowHel framework [13]. Aside fromevolution terms, the remaining terms in H (1) R are obtained from the infrared-limit of the6eal corrections [43] using the method initially proposed in [44, 45]. Additionally, we recal-culate the one-loop soft function ˜S (1) [25, 42]. The dependence on N in the soft function ˜S R enters only through the argument of α s in Eq. (5).Substituting the expression for the running coupling, we obtain up to NNLL accuracyfor the soft matrix evolution factors in Eq. (15) U R ( N, Q , { m } , Q ) = (cid:18) + α s ( µ ) π (1 − λ ) K (cid:19) (cid:104) e g s ( N ) −→ λ (1) (cid:105) D (cid:18) − α s ( µ ) π K (cid:19) , (16) ¯U R ( N, Q , { m } , Q ) = (cid:18) − α s ( µ ) π K † (cid:19) (cid:20) e g s ( N ) (cid:16) −→ λ (1) (cid:17) ∗ (cid:21) D (cid:18) + α s ( µ ) π (1 − λ ) K † (cid:19) , (17)where: g s ( N ) = 12 πb (cid:26) log(1 − λ ) + α s ( µ ) (cid:20) b b log(1 − λ )1 − λ − γ E b λ − λ + b log (cid:18) Q µ (cid:19) λ − λ (cid:21)(cid:27) , (18)and λ = α s ( µ ) b log N. (19)The U R and ¯U R factors in Eqs. (16), (17) correspond to evolution from Q/ ¯ N up to Q anddepend on µ R only through the argument α s . The N -independent evolution from µ R to Q is incorporated into the hard function, as noted earlier.The other factors contributing to the resummation of logarithms, i.e. the radiativefactors for the incoming partons, ∆ i and ∆ j are widely known. The results at NLLaccuracy can be found for example in [39, 46] and at the NNLL level in [40].As already noted, at NLL accuracy, by changing the colour basis to R -basis, the pathordered exponentials in Eq. (4) reduce to simple exponentials. Equivalently, the NLL accu-racy can be obtained by neglecting terms suppressed by a factor of α s in Eqs. (16), (17) and (18).This results in the soft matrix evolution factors turning into exponential functions for theeigenvalues of the soft anomalous dimension matrix. At NLL, it is also enough to onlyknow the LO contributions to the hard and soft function, which results in the followingexpression for the resummed cross section in the Mellin space d ˜ˆ σ (NLL) ij → klB dQ ( N, Q , { m } , µ , µ ) = H (0) R,IJ ( Q , { m } ) ˜S (0) R,JI × ∆ i ( N + 1 , Q , µ , µ )∆ j ( N + 1 , Q , µ , µ ) × exp (cid:20) log(1 − λ )2 πb (cid:16)(cid:16) λ (1) J (cid:17) ∗ + λ (1) I (cid:17)(cid:21) , (20)where the color indices I and J are implicitly summed over. The trace of the product oftwo matrices H R (0) and ˜S (0) R returns the LO cross section.The incoming parton radiativefactors ∆ i are now considered only at NLL accuracy.In order to improve the accuracy of the numerical approximation provided by resum-mation, it is customary to include terms up to O ( α s ) in the expansion of the hard and soft7unction leading to d ˜ˆ σ (NLL w C ) ij → klB dQ ( N, Q , { m } , µ , µ ) = H R,IJ ( Q , { m } , µ , µ ) ˜S R,JI ( Q , { m } ) × ∆ i ( N + 1 , Q , µ , µ )∆ j ( N + 1 , Q , µ , µ ) × exp (cid:20) log(1 − λ )2 πb (cid:16)(cid:16) λ (1) J (cid:17) ∗ + λ (1) I (cid:17)(cid:21) . (21)where H R ˜S R = H (0) R ˜S (0) R + α s π (cid:104) H (1) R ˜S (0) R + H (0) R ˜S (1) R (cid:105) . We will refer to this result as "NLL w C ", since the N -independent O ( α s ) terms in the hardand soft function are often collected together in one function, known as the hard matchingcoefficient, C . Although we choose to treat these terms as in Eq. (21), we keep the name"w C " ("w" standing for "with") as a useful shorthand.The resummation-improved cross sections for the pp → t ¯ tH process are obtainedthrough matching the resummed expressions with the full NLO cross sections dσ (matched) h h → klB dQ ( Q , { m } , µ , µ ) = dσ (NLO) h h → klB dQ ( Q , { m } , µ , µ ) (22) + dσ (res − exp) h h → klB dQ ( Q , { m } , µ , µ ) with dσ (res − exp) h h → klB dQ ( Q , { m } , µ , µ ) = (cid:88) i,j (cid:90) C dN πi ρ − N f ( N +1) i/h ( µ ) f ( N +1) j/h ( µ ) × d ˜ˆ σ (res) ij → klB dQ ( N, Q , { m } , µ , µ ) − d ˜ˆ σ (res) ij → klB dQ ( N, Q , { m } , µ , µ ) | (NLO) , (23)where "matched" can stand for "NLO+NNLL", "NLO+NLL" or "NLO+NLL w C " and"res" for "NNLL", "NLL" or "NLL w C ", correspondingly. The inclusive cross sectionis obtained by integrating the invariant mass distribution given in Eq. (15) over Q and ˆ σ (res) ij → klB ( N, µ , µ ) | (NLO) represents its perturbative expansion truncated at NLO. The mo-ments of the parton distribution functions (pdf) f i/h ( x, µ ) are defined in the standard way f ( N ) i/h ( µ ) ≡ (cid:90) dx x N − f i/h ( x, µ ) . The inverse Mellin transform (23) is evaluated numerically using a contour C in thecomplex- N space according to the “Minimal Prescription” method developed in Ref. [39]. pp → t ¯ tH process at NLO+NNLLaccuracy In this section we present and discuss our state-of-the-art NLO+NNLL predictions forthe t ¯ tH production process at the LHC for two collision energies √ S = 13 TeV and8 S = 14 TeV. The results for the total cross section which we present below are obtained byintegrating out the invariant mass distribution over invariant mass Q . The distribution in Q undergoes resummation of soft gluon corrections in the threshold limit ˆ s → Q , i.e. in theinvariant mass kinematics. This approach is different from directly resumming correctionsto the total cross section in the absolute threshold limit ˆ s → M , which we performedin [23]. Numerical results involving resummation are obtained using two independentlydeveloped in-house computer codes. Apart from NLL and NNLL predictions matched toNLO according to Eq. (22), we also show the NLL predictions supplemented with the O ( α s ) non-logarithmic contributions ("NLL w C "), also matched to NLO.In the phenomenological analysis we use m t = 173 GeV and m H = 125 GeV. TheNLO cross section is calculated using the aMC@NLO code [47]. We perform the currentanalysis employing PDF4LHC15_100 sets [48–53] and use the corresponding values of α s .In particular, for the NLO+NLL predictions we use the NLO sets, whereas the NLO+NNLLpredictions are calculated with NNLO sets. For the sake of comparison with Broggio etal. [28], we adopt the same choice of pdfs, i.e. MMHT2014 [50].We present most of our analysis for two choices of the central values of the renormal-ization and factorization scales: µ = µ F , = µ R , = Q and µ = µ F , = µ R , = M/ . Theformer choice is motivated by invariant mass Q being the natural scale for the invariantmass kinematics used in resummation. The latter choice of the scale is often made in theNLO calculations of the total cross section reported in the literature, see e.g. [20]. Bystudying results for these two relatively distant scales, we aim to cover a span of scalechoices relevant in the problem. The theoretical error due to scale variation is calculatedusing the so called 7-point method, where the minimum and maximum values obtainedwith ( µ F /µ , µ R /µ ) = (0 . , . , (0 . , , (1 , . , (1 , , (1 , , (2 , , (2 , are considered.For reasons of technical simplicity, the pdf error is calculated for the NLO predictions,however we expect that adding the soft gluon corrections only minimally influences thevalue of the pdf error.As discussed in the previous section for the evaluation of the first-order hard functionmatrix H (1) IJ we need to know one-loop virtual corrections to the process, decomposed intovarious colour transitions IJ . We extract them numerically using the publicly available PowHel implementation of the t ¯ tH process [13]. In particular, we use analytical rela-tions to translate between virtual corrections split into various colour configurations inthe colour flow basis used in [13] and our default singlet-octet(s) bases. We cross checkthe consistency of results obtained in this way by comparing the colour-summed one-loopvirtual contributions to Tr (cid:104) H ( ) ˜S ( ) (cid:105) with the full one-loop virtual correction given by the PowHel package [13], as well as the
POWHEG implementation of the t ¯ tH process [14] and thestandalone MadLoop implementation in aMC@NLO [11].We begin the discussion of numerical results by analyzing how well the full NLO resultfor the total cross section is approximated by the expansion of the resummed cross sectionup to the same accuracy in α s as in NLO. It was first pointed out in [25] in the contextof the t ¯ tW production and then later in [26] for the t ¯ tH process that the qg productionchannel carries a relatively large numerical significance, especially in relation to the scaleuncertainty. This is due to the fact that a non-zero contribution from the qg channelappears first at NLO, i.e. it is subleading w.r.t. contributions from the q ¯ q and gg channels.Correspondingly, no resummation is performed for this channel and it enters the matchedresummation-improved formula Eq. (22) only through a fixed order contribution at NLO.It is then clear that in order to estimate how much of the NLO result is constituted by9he terms accounted for in the resummed expression, Eq. (15), its expansion should not bedirectly compared with full NLO but with NLO cross section without a contribution fromthe qg channel. We obtain the latter result from the PowHel package [13]. Its comparisonwith the expansion of the resummed expression in Eq. (15) up to NLO accuracy in α s as a function of the scale µ = µ F = µ R is shown in Figure 1 for √ S = 14 TeV and twochoices of the central scale µ = Q and µ = M/ . While in both cases the expansion ofthe resummed cross section differs significantly from the full NLO, the NLO result withthe qg channel contribution subtracted is much better approximated by the expansion,especially for the dynamical scale choice µ = Q and for the fixed scale choice µ ≥ M/ ,for the physically motivated scale choices. Such good agreement lets us conclude that theNNLL resummed formula will indeed take into account a prevailing part of the higher-ordercontributions from the q ¯ q and gg channels to all orders in α s . µ/µ σ ( pp → t ¯ tH + X ) [fb] √ S = 14 TeV µ = Q NLONLL | NLO (w C )NLO (no qg) µ/µ σ ( pp → t ¯ tH + X ) [fb] √ S = 14 TeV µ = M/ NLONLL | NLO (w C )NLO (no qg) Figure 1: Comparison between the expansion of the resummed expression Eq. (15) upto NLO accuracy in α s , the full NLO result and the NLO result without the qg channelcontribution as a function of the scale µ = µ F = µ R .Our numerical predictions for the total cross sections at √ S = 13 TeV and √ S = 14 TeV are shown in Table 1. We report results obtained with our default scale choice µ = Q as well as the fixed scale µ = M/ . Additionally, we also provide results for the ‘in-between’ choice of µ = Q/ . While for these choices of central scale the NLO results varyby 20 %, the variation reduces as the accuracy of resumation increases. In particular, theNLO+NNLL results show a remarkable stability w.r.t. the scale choice. We also observethat the 7-point method scale uncertainty of the results gets reduced with the increasingaccuracy. In particular, for all scale choices, the scale uncertainty of NLO+NNLL crosssection is reduced compared to the NLO scale uncertainty calculated in the same way. Thedegree up to which the scale uncertainty is reduced depends on the specific choice of thecentral scale. For example, for µ = Q/ the theoretical precision of the NLO+NNLLprediction is improved by about 40% with respect to the NLO result, bringing the scale Although the expansion and the NLO results w/o the qg channel contribution agree very well at thislevel of accuracy in α s , since we do not know the second-order hard-matching coefficients we cannot expectan equally good approximation of the NNLO result by the expansion of the NNLL formula. The value of 10% scale error often quoted in the literature relates to a variation by factors of 0.5 or 2around µ = M/ , while here we consider a much wider range between M/ and Q . S [TeV] µ NLO [fb] NLO+NLL[fb] NLO+NLL with C [fb] NLO+NNLL[fb]13 Q +11 . − . +9 . − . +8 . − . +7 . − . Q/ +9 . − . +8 . − . +6 . − . +6 . − . M/ +5 . − . +8 . − . +5 . − . +5 . − . Q +11 . − . +9 . − . +8 . − . +7 . − . Q/ +9 . − . +8 . − . +6 . − . +6 . − . M/ +6 . − . +8 . − . +6 . − . +5 . − . Table 1: Total cross section predictions for pp → t ¯ tH at various LHC collision energiesand central scale choices. The listed error is the theoretical error due to scale variationcalculated using the 7-point method.error calculated with the 7-point method down to less than 6.5% of the central cross sectionvalue. The results shown in Table 1 are further graphically presented in Fig. 2. N L O N L O+ N LL N L O+ N LL ( w C ) N L O+ NN LL σ ( pp → t ¯ tH + X ) [fb] √ S = 13 TeV µ = M/ µ = Q/ µ = Q N L O N L O+ N LL N L O+ N LL ( w C ) N L O+ NN LL σ ( pp → t ¯ tH + X ) [fb] √ S = 14 TeV µ = M/ µ = Q/ µ = Q Figure 2: Graphical illustration of results presented in Table 1.The size of the K NNLL factor measuring the impact of the higher-order logarithmiccorrections, defined as the ratio of the NLO+NNLL to NLO cross sections, is shown inTable 2. It varies depending on the value of the central scale. The variation is almostentirely driven by the scale dependence of the NLO cross section. For the choice µ = Q the K NNLL -factor can be as high as 1.19.Given the conspicuous stability of the NLO+NNLL results, see Fig. 2, we are encour-aged to combine our results obtained for various scale choices. For this purpose we adoptthe method proposed by the Higgs Cross Section Working Group [20]. In this way, weobtain for the t ¯ tH cross section at 13 TeV σ NLO+NNLL = 500 +7 . . − . − . fb , and at 14 TeV σ NLO+NNLL = 604 +7 . . − . − . fb , where the first error is the theoretical uncertainty due to scale variation and the seconderror is the pdf uncertainty. 11 S [TeV] µ NLO+NNLL [fb] K NNLL factor13 Q +7 . . − . − . Q/ +6 . . − . − . M/ +5 . . − . − . Q +7 . . − . − . Q/ +6 . . − . − . M/ +5 . . − . − . pp → t ¯ tH at various LHCcollision energies and central scale choices. The first error is the theoretical error due toscale variation calculated using the 7-point method and the second is the pdf error.Our findings are further illustrated in the plots in Fig. 3 and Fig. 4. We show therethe scale dependence of t ¯ tH total cross sections calculated with the factorization andrenormalization scale kept equal, µ = µ F = µ R for two LHC collision energies √ S = 13 TeV and √ S = 14 TeV. As readily expected, apart from quantitative differences thereis no visible disparity between the qualitative behaviour of results for the two energies.For the central scale choice of µ = Q , we observe a steady increase in the stability ofthe cross section value w.r.t. scale variation as the accuracy of resummation improvesfrom NLO+NLL to NLO+NNLL. Our final NLO+NNLL prediction is characterised bya very low scale dependence if µ F = µ R choice is made. Correspondingly, if calculatedonly along the µ F = µ R direction, the theoretical error on the NLO+NNLL predictiondue to scale variation would be at the level of 1%, which is a significant reduction fromthe 10% variation of the NLO, c.f. Table 1. Results obtained with the scale choice of µ = M/ behave mostly in a similar way. Only in the very low scale regime, µ (cid:46) . M ,the NLO+NNLL cross section shows a stronger scale dependence. For this scale choice, therise of the matched resummed predictions with the diminishing scale is driven by the fallof the expanded resummed NLL | NLO results, cf. Fig.1, and a therefore is a consequence ofthe relatively large scale dependence of NLO contributions stemming from the qg channel.We further investigate the dependence on the scale but showing separately the renor-malization and factorization scale dependence while keeping the other scale fixed. Fig. 5shows the dependence on µ R and Fig. 6 on µ F for the √ S = 14 TeV. We conclude thatthe weak scale dependence present when the scales are varied simultaneously is a result ofthe opposite behaviour of the total cross section under µ F and µ R variations. The effect issimilar to the cancellations between renormalization and factorization scale dependenciesfor threshold resummation in the absolute threshold limit which we observed in [23]. Thetypical decrease of the cross section with increasing µ R originates from running of α s . Thebehaviour under variation of the factorization scale, on the other hand, is related to theeffect of scaling violation of pdfs at probed values of x . In this context, it is interesting toobserve that the NLO+NLL predictions in Fig. 6 show very little µ F dependence aroundthe central scale, in agreement with expectation of the factorization scale dependence inthe resummed exponential and in the pdfs cancelling each other, here up to NLL. Therelatively strong dependence on µ F of the NLO+NLL (w C ) predictions can be then easilyunderstood: the resummed expression will take into account higher-order scale-dependentterms which involve higher-order terms of both logarithmic (in N ) and non-logarithmicorigin. The latter terms do not have their equivalent in the pdf evolution since the pdfs12 .00100.00200.00300.00400.00500.00600.00700.00800.000.2 0.5 1 2 5 µ/µ σ ( pp → t ¯ tH + X ) [fb] √ S = 13 TeV µ = Q NLONLO+NLLNLO+NLL (w C )NLO+NNLL µ/µ σ ( pp → t ¯ tH + X ) [fb] √ S = 13 TeV µ = M/ NLONLO+NLLNLO+NLL (w C )NLO+NNLL Figure 3: Scale dependence of the total cross section for the process pp → t ¯ tH at the LHCwith √ S = 13 TeV. Results shown for the choice µ = µ F = µ R and two central scale values µ = Q (left plot) and µ = M/ (right plot).do not carry any process-specific information. Correspondingly, the µ F dependence doesnot cancel and can lead to strong effects if the non-logarithmic terms are numericallysignificant.Given apparent cancellations between µ R and µ F scale dependence, we believe thatthe 7-point method of estimating the scale error, allowing for an independent variationof µ R only (for µ F = µ ), is better suited here as an estimate of the theory error thanthe often used variation of µ = µ F = µ R . Another reason for our preference of thisconservative estimate is presence of the hard and soft functions in the resummation formula,Eq. (15), which involve virtual corrections and are known only up to the order α s . Dueto suppression of the LO phase space, they provide a relatively significant part of theNLO+NNLL corrections to the total cross sections, cf. Table 1. It is then justified tosuppose that a similar situation might take place also at higher logarithmic orders andthat the value of the yet unknown two-loop virtual corrections which feed into the second-order coefficients in Eq. (15) can have a non-negligible impact on the predictions. Withthe 7-point method error estimate, we expect that this effect is included within the size ofthe error.Our observation of stability of the predictions w.r.t. scale variation is also confirmed atthe differential level. In Fig. 7 we show the differential distribution in the invariant mass Q of the t ¯ tH system produced at √ S = 14 TeV. While the NLO distributions calculatedwith µ = Q and µ = M/ differ visibly, the NLO+NNLL distributions for these scalechoices are very close in shape and value. The stability of the NLO+NNLL distributionw.r.t. the scale choice is demonstrated explicitly in Fig. 8. Correspondingly, the ratios ofthe NLO+NNLL to NLO distributions differ. In particular for the choice of µ = Q theNNLL differential K-factor grows with the invariant mass and can be higher than 1.2 atlarge Q . The scale error for the invariant mass distribution is also calculated using the7-point method. The error bands are slightly narrower for the NLO+NNLL distributionsthan at NLO. If the scale errors were calculated by variation of µ = µ F = µ R by factors of0.5 and 2, the NLO+NNLL error bands would be considerably narrower.We complete this part of the discussion by comparing resummed results obtained using13 .00200.00400.00600.00800.001000.000.2 0.5 1 2 5 µ/µ σ ( pp → t ¯ tH + X ) [fb] √ S = 14 TeV µ = Q NLONLO+NLLNLO+NLL (w C )NLO+NNLL µ/µ σ ( pp → t ¯ tH + X ) [fb] √ S = 14 TeV µ = M/ NLONLO+NLLNLO+NLL (w C )NLO+NNLL Figure 4: Scale dependence of the total cross section for the process pp → t ¯ tH at the LHCwith √ S = 14 TeV. Results shown for the choice µ = µ F = µ R and two central scale values µ = Q (left plot) and µ = M/ (right plot).the invariant mass kinematics with those obtained earlier by us in the absolute massthreshold limit [23]. At 13 (14) TeV, our most accurate prediction in these kinematics, i.e.the NLO+NLL cross section including the first-order hard-matching coefficient, evaluatedwith PDF4LHC15_100 pdf sets, amounted to σ NLO+NLL w C = 530 +7 . − . (641 +7 . − . ) . Theabsolute mass threshold approach allows only for a fixed scale choice, which is taken tobe µ = µ F = µ R = M/ . Comparing this result with our NLO+NLL predictions for thesame scale choice in the invariant mass kinematics, cf. Table 1, we see that the resultscalculated using the two resummation methods agree within errors.In the remaining part of this section we comment on the relation of our results to theresults of Broggio et al. [28]. That work relies on a resummation formula derived in theSCET framework in [26], though for the purpose of numerical calculations the Mellinspace is adopted. In order to facilitate a comparison with results of [28] we recalculateour results as a function of the scale µ = µ F = µ R using MMHT2014 pdfs as in [28]. Theoutcome is presented in Fig. 9, where we show the NLO cross section and the matchedresummed cross sections at various accuracy as a function of µ = µ F = µ R for the rangeof scales same as in Fig. 1 of [28]. Comparing the two figures, we find a qualitativelysimilar behaviour of the NLO+NNLL cross sections as a function of the scale. Likewise,we obtain bigger NNLL corrections for the µ = µ F = µ R = Q scale choice than for µ = µ F = µ R = Q/ . However, our NLO+NNLL results appear to be more stable wrt.the scale variation, leading to very little difference between the predictions for µ = Q/ and µ = Q (cf. also Table 1). Fig. 9 additionally illustrates another feature of our results,namely that for physically relevant values of µ (cid:38) . Q the scale dependence diminishes asthe accuracy of the predictions increases, independently on the choice of the central scale µ . However, it has to be noted that the scale choices made to obtain results reported inthis paper and [28] are not equivalent. While our resummed expressions depend on µ F and µ R , the formulas used in [28] contain dependence on the hard and soft scales µ h and µ s , as well as µ F . The µ s scale in [28] is chosen in such a way as to mimic the scale ofsoft radiation in the Mellin-space framework, i.e. µ s = Q/ ¯ N . Furthermore, for a given14 .00200.00400.00600.00800.001000.000.2 0.5 1 2 5 µ R /µ σ ( pp → t ¯ tH + X ) [fb] √ S = 14 TeV µ F = µ = Q NLONLO+NLLNLO+NLL (w C )NLO+NNLL µ R /µ σ ( pp → t ¯ tH + X ) [fb] √ S = 14 TeV µ F = µ = M/ NLONLO+NLLNLO+NLL (w C )NLO+NNLL Figure 5: Renormalization scale dependence of the total cross section for the process pp → t ¯ tH at the LHC with √ S = 14 TeV and µ F = µ F , kept fixed. Results shown twocentral scale choces µ = µ F , = µ R , = Q (left plot) and µ = µ F , = µ R , = M/ (rightplot). µ F the resummed central results of [28] are obtained with a fixed hard scale µ h = Q ,while the exact and approximate NLO results are evaluated keeping all other scales equalto the factorization scale. There is one choice of factorization scale for which the scalesetting procedure of [28] corresponds to simultaneous variation of µ = µ F = µ R , that is µ F = Q . For this choice we obtain σ NLO+NNLL = 501 . +38 . − . fb, to be compared with . +42 . − . fb reported in [28], i.e. the central results of the two calculations agree within2.5%. (The scale errors given together with the central values are expected to vary due tothe different methods used for calculating them.) At NLO+NLL accuracy we do not findan agreement with [28]. We conclude that the differences in properties of the NLO+NNLLcross sections reported here and in [28] are likely related to handling of scale setting in thetwo resummation approaches. In this work, we have investigated the impact of the soft gluon emission effects on the totalcross section for the process pp → t ¯ tH at the LHC. The resummation of soft gluon emissionhas been performed using the Mellin-moment resummation technique at the NLO+NNLLaccuracy in the three particle invariant mass kinematics. We have considered the differ-ential distribution in the invariant mass as well as the total cross section, obtained byintegrating the distribution. Our NLO+NNLL predictions are very stable with respect toa choice of the central scale µ for the invariant mass distribution, and consequently alsofor the total cross section. As this is not the case for the NLO predictions, the NNLLcorrections vary in size, depending on the choice of the scale. In general, for the energiesand scale choices considered they provide a non-negative modification of the cross section,which for the scale choice of µ = Q can be even higher than 20% at larger values of Q .We estimate the theoretical error due to scale variation by using the 7-point method,allowing for independent variation of renormalization and factorization scales. The overall15 .00200.00400.00600.00800.001000.000.2 0.5 1 2 5 µ F /µ σ ( pp → t ¯ tH + X ) [fb] √ S = 14 TeV µ R = µ = Q NLONLO+NLLNLO+NLL (w C )NLO+NNLL µ F /µ σ ( pp → t ¯ tH + X ) [fb] √ S = 14 TeV µ R = µ = M/ NLONLO+NLLNLO+NLL (w C )NLO+NNLL Figure 6: Factorization scale dependence of the total cross section for the process pp → t ¯ tH at the LHC with √ S = 14 TeV and µ R = µ R, kept fixed. Results shown two central scalechoces µ = µ F , = µ R , = Q (left plot) and µ = µ F , = µ R , = M/ (right plot).size of the theoretical scale error becomes gradually smaller as the accuracy of resumma-tion increases, albeit the reduction is relatively modest. The reduction would have beenmuch more significant if the scale error had been estimated by simultaneous variation ofrenormalization and factorization scales, i.e. of µ = µ F = µ R . However, as it seems thatthe reduction in this case is a result of cancellations between factorization and renormal-ization scale dependencies, we choose a more conservative 7-point approach for estimatingthe error.The stability of NLO+NNLL results w.r.t. the scale choice allows us to derive our bestprediction for the pp → t ¯ tH total cross section at 13 TeV σ NLO+NNLL = 500 +7 . . − . − . fb , and at 14 TeV σ NLO+NNLL = 604 +7 . . − . − . fb , where the first error is the theoretical uncertainty due to scale variation and the seconderror is the pdf uncertainty. We note that the predictions are very close in their centralvalue to the corresponding NLO predictions obtained for the scale choice µ = M/ andare compatible with them within errors, vindicating the appropriatness of this commonlymade choice. However, in comparison with the NLO predictions obtained in this way, ourNLO+NNLL predictions are characterized by the overall smaller size of the theory errorrelated to scale variation. For an equivalent scale choice setup, our NLO+NNLL resultsfor the t ¯ tH production process at the LHC agree with the results previously obtained byBroggio et al. [28]. Acknowledgments
We are grateful to M. Krämer for providing us with a numerical code for NLO t ¯ tH crosssection calculations [9]. We would like to thank Daniel Schwartländer for his input in the16 .010.020.030.040.050.060.070.080.0 σ ( pp → t ¯ tH + X ) [fb per bin] √ S = 14 TeV µ = Q Q [GeV] σ X /σ NLO
NLONLO+NNLL σ ( pp → t ¯ tH + X ) [fb per bin] √ S = 14 TeV µ = M/ Q [GeV] σ X /σ NLO
NLONLO+NNLL
Figure 7: Comparison of the NLO+NNLL and NLO invariant mass distributions for theprocess pp → t ¯ tH at the LHC with √ S = 14 TeV. Results shown for two central scalechoices µ = Q (left plot) and µ = M/ (right plot). Lower panels show the ratio of thedistributions w.r.t. the NLO predictions.later stages of this work. AK gratefully acknowledges valuable exchanges with A. Kardosand Z. Trocsanyi. This work has been supported in part by the DFG grant KU 3103/1.Support of the Polish National Science Centre grant no. DEC-2014/13/B/ST2/02486 iskindly acknowledged. This work was also partially supported by the U.S. National ScienceFoundation, under grants PHY-0969510, the LHC Theory Initiative, PHY-1417317 andPHY-1619867. AK would like to thank the Theory Group at CERN, where part of thiswork was carried out, for its kind hospitality. TS acknowledges support in the form of aWestfälische Wilhelms-Universität Internationalisation scholarship. References [1] G. Aad et al. [ATLAS Collaboration], Phys. Lett. B (2015) 222[arXiv:1409.3122 [hep-ex]].[2] G. Aad et al. [ATLAS Collaboration], Eur. Phys. J. C (2015) no.7, 349[arXiv:1503.05066 [hep-ex]].[3] V. Khachatryan et al. [CMS Collaboration], JHEP (2014) 087 [JHEP (2014) 106] [arXiv:1408.1682 [hep-ex]].[4] V. Khachatryan et al. [CMS Collaboration], Eur. Phys. J. C (2015) no.6, 251[arXiv:1502.02485 [hep-ex]].[5] G. Aad et al. [ATLAS Collaboration], Phys. Lett. B (2015) 519[arXiv:1506.05988 [hep-ex]].[6] D. de Florian et al. [LHC Higgs Cross Section Working Group], arXiv:1610.07922[hep-ph]. 17 .010.020.030.040.050.060.070.080.0 σ ( pp → t ¯ tH + X ) [fb per bin] √ S = 14 TeV0.801.001.20 500 600 700 800 900 1000 1100 1200 Q [GeV] σ X /σ NLO+NNLL | µ = Q NLO+NNLL ( µ = M/ )NLO+NNLL ( µ = Q ) Figure 8: Comparison of the NLO+NNLL invariant mass distributions for the process pp → t ¯ tH at the LHC with √ S = 14 TeV calculated with µ = Q and µ = M/ . Lowerpanel shows the ratio to the NLO+NNLL distribution with µ = Q .[7] The ATLAS collaboration [ATLAS Collaboration], ATLAS-CONF-2016-058,ATLAS-CONF-2016-067, ATLAS-CONF-2016-068, ATLAS-CONF-2016-080.[8] CMS Collaboration [CMS Collaboration], CMS-PAS-HIG-15-008,CMS-PAS-HIG-16-004, CMS-PAS-HIG-16-020, CMS-PAS-HIG-16-022,CMS-PAS-HIG-16-038, CMS-PAS-HIG-17-004.[9] W. Beenakker, S. Dittmaier, M. Krämer, B. Plumper, M. Spira and P. M. Zerwas,Phys. Rev. Lett. (2001) 201805 [hep-ph/0107081]; W. Beenakker, S. Dittmaier,M. Krämer, B. Plumper, M. Spira and P. M. Zerwas, Nucl. Phys. B (2003) 151[hep-ph/0211352].[10] L. Reina and S. Dawson, Phys. Rev. Lett. (2001) 201804 [hep-ph/0107101];L. Reina, S. Dawson and D. Wackeroth, Phys. Rev. D (2002) 053017[hep-ph/0109066]; S. Dawson, L. H. Orr, L. Reina and D. Wackeroth, Phys. Rev. D (2003) 071503 [hep-ph/0211438]; S. Dawson, C. Jackson, L. H. Orr, L. Reina andD. Wackeroth, Phys. Rev. D (2003) 034022 [hep-ph/0305087].[11] V. Hirschi, R. Frederix, S. Frixione, M. V. Garzelli, F. Maltoni and R. Pittau, JHEP (2011) 044 [arXiv:1103.0621 [hep-ph]].[12] R. Frederix, S. Frixione, V. Hirschi, F. Maltoni, R. Pittau and P. Torrielli, Phys.Lett. B (2011) 427 [arXiv:1104.5613 [hep-ph]].[13] M. V. Garzelli, A. Kardos, C. G. Papadopoulos and Z. Trocsanyi, Europhys. Lett. (2011) 11001 [arXiv:1108.0387 [hep-ph]].[14] H. B. Hartanto, B. Jager, L. Reina and D. Wackeroth, Phys. Rev. D (2015) 9,094003 [arXiv:1501.04498 [hep-ph]].[15] S. Frixione, V. Hirschi, D. Pagani, H. S. Shao and M. Zaro, JHEP (2014) 065[arXiv:1407.0823 [hep-ph]]. 18 µ/µ σ ( pp → t ¯ tH + X ) [fb] √ S = 13 TeV µ = Q NLONLO+NLLNLO+NLL (w C )NLO+NNLL Figure 9: Scale dependence of the total cross section for the process pp → t ¯ tH at theLHC with √ S = 13 TeV, calculated using MMHT 2014 pdfs. Results shown for the choice µ = µ F = µ R and the central scale value µ = Q .[16] S. Frixione, V. Hirschi, D. Pagani, H.-S. Shao and M. Zaro, JHEP (2015) 184[arXiv:1504.03446 [hep-ph]].[17] Y. Zhang, W. G. Ma, R. Y. Zhang, C. Chen and L. Guo, Phys. Lett. B (2014) 1[arXiv:1407.1110 [hep-ph]].[18] A. Denner and R. Feger, JHEP (2015) 209 [arXiv:1506.07448 [hep-ph]].[19] A. Denner, J. N. Lang, M. Pellen and S. Uccirati, JHEP (2017) 053[arXiv:1612.07138 [hep-ph]].[20] S. Dittmaier et al. [LHC Higgs Cross Section Working Group Collaboration],arXiv:1101.0593 [hep-ph].[21] R. Bonciani, S. Catani, M. L. Mangano and P. Nason, Phys. Lett. B (2003) 268[hep-ph/0307035].[22] S. M. Aybat, L. J. Dixon and G. F. Sterman, Phys. Rev. Lett. (2006) 072001[hep-ph/0606254]; S. M. Aybat, L. J. Dixon and G. F. Sterman, Phys. Rev. D (2006) 074004 [hep-ph/0607309].[23] A. Kulesza, L. Motyka, T. Stebel and V. Theeuwes, JHEP (2016) 065[arXiv:1509.02780 [hep-ph]].[24] A. Kulesza, L. Motyka, T. Stebel and V. Theeuwes, arXiv:1609.01619 [hep-ph].[25] H. T. Li, C. S. Li and S. A. Li, Phys. Rev. D (2014) 9, 094009 [arXiv:1409.1460[hep-ph]].[26] A. Broggio, A. Ferroglia, B. D. Pecjak, A. Signer and L. L. Yang, JHEP (2016)124 [arXiv:1510.01914 [hep-ph]]. 1927] A. Broggio, A. Ferroglia, G. Ossola and B. D. Pecjak, JHEP (2016) 089[arXiv:1607.05303 [hep-ph]].[28] A. Broggio, A. Ferroglia, B. D. Pecjak and L. L. Yang, JHEP (2017) 126[arXiv:1611.00049 [hep-ph]].[29] A. Broggio, A. Ferroglia, G. Ossola, B. D. Pecjak and R. D. Sameshima,arXiv:1702.00800 [hep-ph].[30] G. Sterman and M. Zeng, JHEP (2014) 132 [arXiv:1312.5397 [hep-ph]].[31] N. Kidonakis, G. Oderda and G. Sterman, Nucl. Phys. B , 365 (1998).[arXiv:hep-ph/9803241].[32] H. Contopanagos, E. Laenen and G. F. Sterman, Nucl. Phys. B (1997) 303[hep-ph/9604313].[33] M. Czakon, A. Mitov and G. F. Sterman, Phys. Rev. D (2009) 074017[arXiv:0907.1790 [hep-ph]].[34] N. Kidonakis and G. Sterman, Nucl. Phys. B (1997) 321[arXiv:hep-ph/9705234].[35] L. J. Dixon, L. Magnea and G. F. Sterman, JHEP (2008) 022 [arXiv:0805.3515[hep-ph]].[36] N. Kidonakis, G. Oderda and G. Sterman, Nucl. Phys. B , 299 (1998)[arXiv:hep-ph/9801268].[37] A. Ferroglia, M. Neubert, B. D. Pecjak and L. L. Yang, Phys. Rev. Lett. (2009)201601 [arXiv:0907.4791 [hep-ph]].[38] A. Ferroglia, M. Neubert, B. D. Pecjak and L. L. Yang, JHEP (2009) 062[arXiv:0908.3676 [hep-ph]].[39] S. Catani, M. L. Mangano, P. Nason and L. Trentadue, Nucl. Phys. B , 273(1996) [arXiv:hep-ph/9604351].[40] S. Catani, D. de Florian, M. Grazzini and P. Nason, JHEP (2003) 028[hep-ph/0306211].[41] A. J. Buras, Rev. Mod. Phys. (1980) 199.[42] V. Ahrens, A. Ferroglia, M. Neubert, B. D. Pecjak and L. L. Yang, JHEP (2010) 097 [arXiv:1003.5827 [hep-ph]].[43] S. Catani, S. Dittmaier, M. H. Seymour and Z. Trocsanyi, Nucl. Phys. B (2002)189 [hep-ph/0201036].[44] W. Beenakker, S. Brensing, M. Kramer, A. Kulesza, E. Laenen and I. Niessen, JHEP (2012) 076 [arXiv:1110.2446 [hep-ph]].[45] W. Beenakker et al. , JHEP (2013) 120 [arXiv:1304.6354 [hep-ph]].2046] R. Bonciani, S. Catani, M. L. Mangano and P. Nason, Nucl. Phys. B , 424(1998) [arXiv:hep-ph/9801375].[47] J. Alwall et al. , JHEP (2014) 079 [arXiv:1405.0301 [hep-ph]].[48] J. Butterworth et al. , J. Phys. G (2016) 023001 [arXiv:1510.03865 [hep-ph]].[49] S. Dulat et al. , Phys. Rev. D (2016) no.3, 033006 [arXiv:1506.07443 [hep-ph]].[50] L. A. Harland-Lang, A. D. Martin, P. Motylinski and R. S. Thorne, Eur. Phys. J. C (2015) no.5, 204 [arXiv:1412.3989 [hep-ph]].[51] R. D. Ball et al. [NNPDF Collaboration], JHEP (2015) 040 [arXiv:1410.8849[hep-ph]].[52] J. Gao and P. Nadolsky, JHEP (2014) 035 [arXiv:1401.0013 [hep-ph]].[53] S. Carrazza, S. Forte, Z. Kassabov, J. I. Latorre and J. Rojo, Eur. Phys. J. C75