Phenomenology of combined resummation for Higgs and Drell-Yan

We study the phenomenological impact of a recently suggested formalism for the combination of a so-called consistent transverse momentum and threshold resummation, by using it to improve the fixed-order results. This formalism allows for a systematical improvement of the transverse momentum resummation that is valid in the entire range of $p_T$ by the inclusion of the threshold contribution. We use the Borel method as a suitable prescription for defining the inverse Mellin and Fourier transforms in the context of combined resummed expression. The study is applied to two QCD processes, namely the Higgs boson produced via gluon fusion and $Z$ boson production via the Drell--Yan mechanism. We compare our results to the standard transverse momentum resummation, as well as to the fixed-order results. We find that consistent transverse momentum resummation leads to faster perturbative convergence at small-$p_T$ while the inclusion of threshold resummation improves the agreement with fixed-order calculations. These effects are more pronounced in the case of Higgs which is known to have slower perturbative convergence.


Introduction
The current state-of-the-art accuracy in resummed calculations of transverse momentum spectra is N3LL matched with fixed-order predictions (NLO/NNLO), and has been achieved recently for a variety of processes including Higgs boson production [1][2][3], and Drell-Yan (DY) [4][5][6][7][8][9]. Since large logarithmic corrections present in transverse momentum and threshold resummation both originate from the emission of soft gluons, there have been attempts to construct a joint formalism that simultaneously resum logarithmic contributions that are enhanced at small-p T and at partonic threshold [10,11]. This joint resummation has been successful in producing phenomenological results at NLL accuracy for various processes including Higgs [12] and vector boson production via DY mechanism [13] which has been recently extended to NNLL [14]. While the aforementioned resummation were done in Fourier-Mellin space, joint resummation in direct space has been achieved up to NNLL accuracy using SCET [15]. In these references, however, the inclusion of soft contributions were not fully complete as will be discussed below.
Recently, a combination of threshold and transverse momentum resummation has been proposed in ref. [16] that has the following features: it reproduces the correct behaviour to any desired logarithmic order in the limit p T → 0 for fixed x; it reduces to threshold resummation in the soft limit x → 1 for fixed p T up to power correction in (1 − x); and it reproduces the total cross section at any given logarithmic order in the threshold limit upon integration over the transverse momentum. This combined resummation relies mainly on two main ingredients: (i) a modified transverse momentum resummation-that henceforth we call consistent resummation-which leads to threshold resummed expression upon integration over p T , (ii) the combination of the consistent expression with the pure threshold expression that takes into account all the logarithmic enhanced terms in the soft limit for finite-p T . Such a combined resummation is expected to allow for a systematic improvements of the transverse momentum that is valid for the entire range of p T but the phenomenological implications have never been studied.
Therefore, the main goal of this paper is to separately assess the effects of consistent and combined resummation on transverse momentum distributions. In order to do this, we have to deal with issues that arise from the particular construction of the consistent resummation. The first concerns the issue of performing the Fourier-Mellin inverse. In most resummation formulae [17][18][19][20][21][22][23][24][25][26][27][28][29][30], transverse momentum resummation is performed in Fourier-Mellin (N −b) space where closed expression can be found. However, in addition to the problem of Landau pole that prevents the existence of an inverse Mellin in standard resummation, consistent resummation presents extra-singularities due to the interplay between the Mellin moment N and the impact parameter b in the argument of the logarithms. We show that this issue can be addressed by slightly modifying the Borel prescription that was recently studied in ref. [31,32] in the context of threshold resummation for DY cross sections. The second has to do with logarithmic counting in which perturbative evolutions have to be treated in a different way to properly account for the threshold behaviours. Finally, in order to obtain valid predictions, we need to match the resummed expressions to the fixed-order results. This requires the computation of the inverse Fourier-Mellin transform order by order in the running of the coupling constant α s . This paper is organized as follows. We begin in Section 2 with a review of combined threshold and consistent resummation. We then describe the analytical procedures to address the issues introduced previously. Our phenomenological results are presented in Section 3 in which we separately study the impact of consistent and combined resummation to the Higgs boson production at LHC via gluon-gluon fusion and a DY process at Tevatron. Finally, conclusions are drawn in Section 4.

Analytic formulation of combined resummation
In this section, we first provide a brief review of the combined resummation. We then address specific issues that arise in the phenomenological studies of the consistent resummation: (i) the need for a Fourier-Mellin inversion procedure, (ii) the treatment of the PDF evolution, and (iii) the matching to the fixed-order calculations.
Regardless of the fact that we are here particularly interested in the phenomenological implication of the combined resummation to the case of Higgs and DY production, the expressions that follow are kept as general as possible such that they apply to general colour singlet hadronic production.

Combined consistent small-p T and threshold resummation
Consider the inclusive hard-scattering process h 1 + h 2 → F (Q, p T ) + X, where the collision of the hadrons h 1 and h 2 produces a final-state system F with an invariant mass Q and a transverse momentum p T accompanied by an arbitrary final state X. As mentioned earlier, in this paper we study the case where the final state F can be either a Higgs boson or a DY lepton pair via electroweak boson production. We denote by √ s and √ŝ the center of mass-energy of the colliding hadrons and partons respectively.
Transverse momentum distributions factorize into a convolution between a parton luminosity and a partonic cross section which can be written in the following way (2. 2) The differential p T distribution in Eq. (2.1) has been expressed in terms of the dimensionless variable ξ p = p 2 T /Q 2 and, for simplicity, the renormalization and factorization scale dependencies have been omitted. The sum is over the different partonic channels where a and b denote partons. As it turns out that in Mellin space, the convolution becomes a simple product, Eq. (2.1) can be re-written as: given that Notice here that the cross section and its Mellin transform are denoted with the same symbol and only distinguished trough their arguments. It is worth stressing that the Mellin transform of the hadronic cross section in Eq. (2.4a) is taken w.r.t. the scaling variable τ while the Mellin transform of the partonic one in Eq. (2.4b) is taken w.r.t.
x. Thanks to this choice of variables, one can take simultaneously the Mellin and Fourier transform of the cross section. Indeed, τ and p T are independent variables as τ ranges from 0 ≤ τ ≤ 1 and p T ranges from 0 ≤ p T ≤ ∞. The available phase induced by such choice of variables is depicted in the right-hand side of Fig. 1.
The partonic cross section in Eq. (2.3) is constructed by combining the consistent improved transverse momentum, to be discussed below, with threshold resummation. The two resummations are combined through a profile matching function which is chosen such that the combined result reproduces consistent and threshold resummation at small and large p T respectively. One possible expression of dσ ab /dξ p has been proposed in ref. [16] where the profile matching function is defined as Small-x Small-pT Threshold Figure 1. Available phase space for the production of a final-state system F with an invariant mass Q and a transverse momentum p T in the three kinematic limits (threshold, collinear, high-energy), parametrized in terms of the kinematic variable τ = Q 2 /s (left) and τ (right).
The values of k and m can be chosen arbitrarily provided that m < k. This is because combined resummation produces results that differ from consistent by O(ξ m p ) corrections when ξ p → 0, and from threshold by O(1/N k ) corrections when N → ∞.
The threshold part in Eq. (2.5) resums logarithms of the form α n s ln m (N ) that are enhanced in the large-N limit for fixed p T and is given in refs. [33,34] for Higgs and in refs. [23,24,26,28] for DY production. On the other hand, the consistent part represents the improved version of the standard transverse momentum resummation from which it differs by the inclusion of soft contributions that are power-suppressed terms in the small-p T limit for fixed N . It provides a modified transverse momentum expression that reproduces threshold resummation upon integration over p T .
While the phenomenological study of the threshold resummation poses no new problems, technical issues arise in the consistent resummation that require specific care. Hence, the following sections are mainly devoted to the consistent resummation. First, in order to fully take into account the kinematics constraint on the transverse momentum conservation, it is advantageous to work in Fourier space where the consistent part in Eq. (2.5) yields In the above expression, σ 0 represents the leading order Born cross section which is given in [29] for Higgs and [35] for DY; whileb = bQ represents the modified conjugate variable to ξ p in Fourier space. The resummed partonic cross section is represented by Σ ab where theb dependence is embodied in λ χ which is defined as with the modified Mellin variable defined asN = N exp(−γ E ) and b 0 = 2 exp(−γ E ). Notice that in comparison to the standard transverse momentum resummation, now the logarithm which is being resummed also depends on the Mellin variable N .
The logarithmically enhanced terms in Eq. (2.7) are contained in the resummed expression Σ ab which can be organized in the following way: where P denotes the different projectors of the LO anomalous dimension in the matrix flavour space. The function G contains all the large logarithms coming from the Sudakov form factor, the evolution of the PDFs, and the evolution of the coefficient functions. Its complete expression is derived in Section 2.3. On the other hand, the functionH contains all the terms that behave as a constant in the large-b limit.

Inverse Fourier-Mellin
In order to get physical cross sections, we have to compute the inverse Fourier-Mellin of the (N − b) space resummed expression. Mellin and/or Fourier resummed expressions generally give rise to technical complications when it comes to performing the inverse transform. Because the resummation corresponds to the asymptotic sum of a divergent series of p Tspace contributions [29], the Fourier-Mellin integral contains singularities. In a standard procedure, Fourier and Mellin inverse transform are computed using the Minimal prescription (MP) [36]. In the case of inverse Fourier transform, this amounts to finding ab contour that avoids the singularities; while for the inverse Mellin transform, this amounts to defining an integration path that passes to the left of the branch but to the right of all singularities as shown in Fig. 3. This prescription defines the resummed hadronic cross section as As a result, cross sections obtained through such a procedure is finite. However, in addition to the branch cut related to the Landau pole, the integrand in Eq. (2.7) has extrasingularities when This means that the location of the N -space singularities depends on b. Therefore, it is non-trivial to find a contour deformation that leaves all the singularities to the left and the branch cut to the right. An alternative approach consists on summing the divergent series using the Borel method which was first introduced in refs [ [36][37][38]. Recently, in ref. [31] it was shown in the context of resummation of transverse momentum distributions that the inverse Fourier-Mellin can be performed by first resumming the divergent series using the Borel method and then removing the divergence in the inverse Borel transform. This permits the computation of the asymptotic result the divergent series is tending to. Such a procedure can be slightly modified to work with transverse momentum resummation [39]. The steps are as follows. First, we expand the resummed component Σ ab (N, χ, α s ) as a series inᾱ s ln(χ) and tackle directly the divergence using the Borel summation. This then leads to a N -space expression that can be inverted using the MP. Hence, we have When k = 0 the zeroth order coefficient h 0 (N, α s ) is just a constant. Putting back the above equation into the partonic resummed expression in Eq. (2.7), exchanging the integral and the sum, and finally using the trick d k χ /d k | →0 = ln k χ, we get to the following result Remarks are in order concerning the above equations: (i) dσ /dξ p is just dσ cons /dξ p without the kinematic factor 1 + ξ p + ξ p 2N and the Born level cross section σ 0 , (ii) b 0 has been factorized out such that the integrand in Eq. (2.13a) is a function of N instead ofN . This allows us to perform the inverse Fourier transform analytically order by order. Hence, where the function M is defined as The function K denotes the modified Bessel function of the second kind. We can then get rid of the derivatives and the limit by performing a contour integral, where H represents a contour enclosing the singularity at ξ = 0.
Having expanded the function Σ allowed us to perform the inverse Fourier transform, but it leaves us with a series representation of the result. Thinking of the numerical limitations of series representations, we would rather trade the sum for an integral. This basically calls for an analytic continuation, which in this case is achieved using a Borel summation and a subsequent Borel back-transformation. First, writing down the series as The corresponding Borel sum is then given by where we resummed back the expansion of Σ using its definition in Eq. (2.12). We finally perform the Borel back-transformation, which represents our final resummed expression when multiplied by the factor σ 0 ( 1 + ξ p + ξ p ) 2N that we left out from the beginning. After cutting off the w-integration at some finite number C, we have the final expression dσ cons The Borel and Minimal Prescription differ in the way higher-twist behaviour of the resummed series is handled. However, in practice and as will be demonstrated in the next sections, this difference is negligible at collider energies we are interested in. The inclusion or exclusion of subleading terms in the Borel expression of Eq. (2.20) is regulated by the cutoff C which can be used to estimate the accuracy of the resummation procedure. The larger the value of C gets, the more we include power-correction terms in the series. However, as C gets bigger the integral in Eq. (2.20) becomes unstable and varies a lot for small changes of C eventually spoiling the accuracy of the result. This is due to the fact that the series in Eq. (2.12) is not Borel summable.
The following subsections deal with the numerical implementation of the procedure described above and, specifically, comment on the proper choice of the cutoff C regulating the truncation of the Borel integral making it convergent.

Contour and Borel integration
By examining the singularity structure of the integrand in Eq. (2.20) (see Fig. 2), we see that the contour integral can be written as an integral over the domain [0, 1], where f represents the integrand. The contour is parametrized by ξ(ϕ) = x 0 + R exp(2iπϕ) which is a circle around the real point ξ = x 0 ∈ R with a radius R > 0, in the positive direction around x 0 . The branch cut and poles of f are now mapped in terms of the variable ξ into a branch cut on the real axis from zero to w and poles in the negative integers (see Fig. 2). In practice, we choose x 0 = w/2 and R = w/2 + ζ (ζ ∈ R * + ). The value of ζ and the cutoff C must be chosen in such a way that the circle does not intersect the branch cut in ξ-space. Based on the twist-4 argument [40], the value of the cutoff should be chosen according to C ≥ 2. Throughout our implementation, we chose the minimal value for C: i.e C = 2.
The contour H is given by the blue circle.

Inverse Mellin
Thanks to the Borel procedure, the Minimal Prescription can now be applied to perform the Mellin back transformation. The inverse Mellin transform given by Eq. (2.10) can be rewritten in the following way 1 2πi where N 0 , as usual, must be located to the right of all the singularities of the integrand.
In the first line, we split the contour into two pieces, the first above the real axis and the second below. Since dσ/dξ p (τ * ) = [dσ/dξ p (τ )] * is real, the second piece turns out to be the negative of the complex conjugate of the first and hence can be combined to the former.
In the second line, we used the substitution N (u) = N 0 + (r + i) ln u with r > 0 being an arbitrary real and positive parameter which controls the slope of the path that enhances the numerical convergence of the integral. This parameter must be positive, so that for large imaginary values of N (i.e. small values of u), the prefactor converges to zero, requiring the real part of N (u) to be negative for u → 0 (i.e. N 0 − r ln 1/u < 0). Notice that the following is also an equivalent transformation:

PDF evolution for consistent resummation
In order to compute hadronic cross sections, we need to combine the partonic resummed cross section with the Parton Distribution Functions (PDFs) embodied in the parton luminosity L in Eq.
where hereH represents the hard factor which behaves as a constant in the large-b limit, C is the coefficient functions, U is the evolution operators that evolve the PDFs from the scale Q 2 /χ to µ F , and finally, S c is the Sudakov form factor. The subscripts define partonic indices where repeated indices are summed over. In the case of the Higgs production via gluon fusion, for instance, c = g. Finally, µ R and µ F represent the renormalization and factorization scale respectively. From Eq. (2.25a), we can see that the consistent resummed expression differs from the standard small-p T resummation in two ways. First, the argument of the logarithms in the Sudakov exponent includes the soft contribution. Second, as mentioned above, the argument of the running of coupling, both in the coefficient functions and in the evolution operators, is computed at the scale Q 2 /χ instead of b 2 0 /b 2 . In order to determine the appropriate logarithmic counting, we work out in this section the solution of the evolution equation explicitly. Let us first start with the coefficient functions C of Eq. (2.25a). One can set all the arguments of α s in C to be the same and take into account the difference in scales in an evolution factor R such that where (2.27) Let us now turn to the PDF evolution. Here, we can evolve the PDFs from Q 2 /χ to Q 2 and then from Q 2 to µ 2 F . In the case where Q = µ F , we have U ij N, α s (Q 2 ), α s (µ 2 F ) = 1. The expression of the PDF evolution can be written as [41]: where U (LO) and V respectively represent the lowest and higher perturbative order solutions to the evolution equation, andṼ ka denotes the k, a-element of the inverse matrix V −1 in the flavour space. Henceforth, we shall use the boldface notation to denote the representation in flavour matrix space. In the case of singlet, the lowest perturbative order solution is derived by diagonalizing the LO anomalous dimensions matrix γ (0) . Thus, following refs. [29,41], we have where λ ± and e ± represent the eigenvalues of the singlet matrix γ (0) and its projectors respectively whose expressions are given by Eq.(2.27) and Eq.(2.28) of [41]. On the other hand, the V functions which take into account the higher order solutions to the evolution equation can be perturbatively expanded as a series in α s , Each coefficient V (i) in Eq. (2.30) can be computed iteratively following [41]. For instance, the first order coefficient V (1) is given by where we have introduced the coefficients β i of the QCD β-function. Once again, one can take into account the difference in the argument of the running of coupling of V (N, α s (Q 2 /χ)) in Eq. (2.28a) by introducing an evolution factor where similarly to Eq. (2.26), (2.33) At this point, we can discuss the difference in counting between consistent and standard small-p T resummation. In the context of standard transverse momentum resummation, the leading-order solution to the evolution equation contributes with a single logarithm of b. However, from the point of view of the threshold resummation, the flavour diagonal anomalous dimensions carry an additional contribution of the form A p T g/q ln(N 2 ). To see this, let us compute the term in the exponent of the lowest-order solution to the evolution equation given by Eq. (2.29): Taking the large-N limit of the anomalous dimension leads to the following result where we have defined λ N = α s β 0 ln(N 2 ). Therefore, in order to compute N k LL consistent resummation, the large-N behaviour of the evolution has to be included up to N k LO, which is not the case in the standard transverse momentum resummation as N k−1 LO is enough.
That is, at NNLL, the solution of evolution is performed up to NNLO accuracy but with the NNLO anomalous dimension substituted by its large-N behaviour. Combining all the results, we can organize terms in such a way that we obtain the form presented in Eq. (2.9). The hard functionH is then given bỹ where here {P} = {r, p} = ± and we have definedC cl = C ci V il . Finally, the universal form factor in the exponent of Eq. (2.9) which contains all the logarithmic enhanced terms is expressed as where the evolution factorR cl is defined asR cl = R ci E il . The logarithmic expansion of S c starts at LL accuracy. While the standard small-p T resummation, the expansion of G r/p starts at NLL accuracy, their large-N behaviour already contribute at LL in the context of consistent in order to take into account for the soft behaviour. Similarly to the standard small-p T resummation, the flavour off-diagonal terms ofR cl (c = l) starts to contribute at NLL. The terms in the flavour diagonal, instead, start to contribute at NLL as opposed to the standard resummation procedure where they start to contribute at NNLL. To summarize, the truncation of the consistent resummation at a given logarithmic accuracy is defined in the following way: • At LL, we approximate the hard functionH by its lowest perturbative order and include g 1 in the Sudakov exponent. The G ± functions are included up to LO with γ (0) replaced by its large-N behaviour while theR functions are approximated to 1.
• At NLL, we includeH up to NLO with the function g 2 in the Sudakov exponent. In addition, we include the complete LO term in G ± together with the NLO term where γ (0) and γ (1) are replaced the their large-N behaviour. The full expression ofR gq is also included while we only include the large-N behaviour of V (1) (N ) inR gg .
• At NNLL, we includeH up to NNLO and the function g 3 in the Sudakov exponent. Similarly, we also include the complete NLO term in G ± together with the NNLO term where the anomalous dimensions are replaced by their large-N behaviour. Finally, we add toR gg the large-N behaviour of V (2) (N ).

Finite order truncation of the resummed expression
In order to provide a valid prediction up to NNLL+NLO, we need to match the resummed results to the fixed-order calculation. This matching procedure guarantees the correct behaviour from perturbative calculations up to a specified order and incorporates the large logarithms from resummation at higher-orders. The matching consists on adding to the resummed results to the perturbative calculations truncated at a given order in α s and subtract to the whole the expansion of the resummed results at the same order. The matched expression of the combined matched cross section is The subscript F.O indicates the order at which the expression is truncated. In practice, the logarithmic accuracy of the resummed expression is determined by expanding the resummed component Σ ab in Eq. (2.9) as a series in α s at the same order as the fixed-order results. Such expansion will give rise to a series of the form The perturbative coefficients Σ (n) are polynomials in the logarithm variable ln(χ). They vanish when ln(χ) = 0 (i.e.b = −N 2 b 2 0 ) and lead to threshold resummation when b = 0. On the other hand, the function H (n) contains all the terms that behave as constants in the large-b limit. The explicit expressions of the two perturbative coefficient functions Σ (1) and Σ (2) which contribute at NNLO are given by (2.40b) The functions Σ (i;j) (N ) are purely functions of N . In order to match Eq. (2.39) to the fixedorder results, we need to Fourier-invert logarithms of the form ln k (χ), such computation requires the evaluation of the following integral at a given order in k (2.41) Thus, according to Eq. (2.15), the computation of I k for a fixed value of k just amounts to the computation of the following limit with M expressed as in Eq. (2.15). This allows us to analytically compute the inverse Fourier of the finite order expansion of the consistent resummed expression up to NNLO, where the K (i) 1 are the i-th derivatives of the Bessel function w.r.t. the argument. One can check that in the small-p T limit-meaning taking N → exp(−γ E )-the above expressions exactly reproduce the standard small-p T expressions given by Eqs. (121)-(124) of ref. [29]. The equivalence between taking p T → 0 and setting N → exp(−γ E ) deserves further comment. On can check that the modified logarithmic variable ln(1 +b 2 /b 2 0 ) of ref. [29] is obtained from our modified ln(χ(b)) by settingN = 1 which in Eqs. (2.43) translates into taking N → exp(−γ E ).

Phenomenological studies
In the previous sections, we have developed a framework that allows us to numerically perform the inverse Fourier-Mellin transform of the resummed expression and to match it to the fixed-order result. We now turn to the main result of this paper, namely a detailed phenomenological study of the impact of consistent and combined resummation on transverse momentum distributions in the case of Higgs and Z-boson production. Our aim is to assess the potential impact of: (i) the improved transverse momentum resummation through the inclusion of the soft contributions which we expect to be more visible in the small-p T region, and (ii) the soft resummation at fixed-p T which we expect to improve the matching to fixed-order in the intermediate and large-p T regions.
All results are produced using the NNPDF31_nnlo_as_0118 set of parton distributions from NNPDF3.1 [42] through the LHAPDF [43] interface in the manner described in refs. [32,40] by expanding the parton luminosity in the basis of the Chebyshev polynomials; We expect the results to be largely independent of the PDF set used. In order to assess perturbative uncertainties and higher order corrections, we perform standard variations of factorization and renormalization scales using the seven-point method. In addition to defining a proper value to the cutoff C entering in the Borel prescription, one also has to define the free parameters in the matching function in Eq. (2.6). The following results are produced by setting k = 3 and m = 2. There is some arbitrariness in the definition of k and m; however, we did explicitly check that results do not change provided that k is chosen from 2 to 5 and m < k. The resummed results are produced from our own resummation code while the fixed-order part is provided by the NNLOJET [44] implementation.

Consistent resummation in the large-b limit
In this section, we start by checking that our results reproduce the standard transverse momentum resummation when the soft contributions are switched off. Besides being a consistency check, this allows us to verify that the Borel prescription agrees with the Minimal Prescription (MP) of ref. [29]. It has been shown analytically in ref. [16] that consistent resummation produces exactly standard small-p T in the large-b limit, here, we show numerically that this is indeed the case using the Borel prescription. Takingb → ∞ means that we adopt the following replacements throughout our expressions: Notice that in the expression of χ,N has been replaced by 1. As suggested in ref. [29], shifting the argument of the logarithm by 1 reduces the effect of unjustified small-p T resummation effect in the large-p T region (equivalentlyb → 0). The replacements in Eq. (3.1) implies that we have to change the function M in the following way, Finally, in order to correctly reproduce the standard small-p T resummation result, one has to implement exactly the same solution to the evolution equation as in ref. [29]. The results are shown in Fig. 4 from which one can see that our formalism is consistent with the standard small-p T results from HqT [29] for Higgs and DYqT [35,45] for DY which both use the Minimal Prescription. This confirms two things: first, that the consistent resummation indeed reproduces the standard small-p T in the correct limit; and second, that the Borel method is a meaningful prescription to perform the inverse Fourier-Mellin.

Higgs boson production at LHC and DY production at Tevatron
We begin by studying the standard model Higgs produced via gluon-gluon fusion at LHC with a center of mass energy of √ s = 13 TeV. The computations have been performed in the large-top-mass (m T ) limit with a Higgs mass m H = 125 GeV GeV. As mentioned previously, the factorization and renormalization scales are varied according to the seven-point method by a factor of two in either direction around m H .
The left plots of Fig. 5 study the impact of the pure consistent (Fig. 5c) and combined (Fig. 5e) resummation on the transverse momentum distribution by comparing the results to the NLO results. The standard small-p T (SSpT, henceforth) resummation is shown in Fig. 5a. The top panels show the different order of resummation (NLL in blue and NNLL in red) along with the NLO results (orange). The ratio w.r.t. the NNLL result is shown in the lower panels for the different types of resummation (SSpT, consistent, and combined). In contrast to the NLO fixed-order results, the resummation leads to a well behaved transverse momentum distribution that has peak at p T ∼ 10 GeV and vanishes for p T → 0. As expected, we observe that the scale uncertainty is reduced with the inclusion of the NNLL terms. This feature is seen for all three types of resummations, but it is more pronounced in the small-p T region. Consistent and standard small-p T resummation agree in the intermediate and large-p T region at NLL, but a sizeable difference appears in the small-p T range (below ∼ 25 GeV). Indeed, consistent resummation displays faster convergence and hence leads to a better agreement with the fixed-order computation. As a consequence, in consistent resummation, the NLL result is much closer to the NNLL result than in SSpT resummation. However, at NNLL, the difference between SSpT and consistent resummation is almost invisible, but the consistent resummation has a moderately smaller uncertainty. Turning now to the combined result, consistent and combined resummation are very close up to scales of at least 40 GeV where the pure threshold resummation starts to contribute. On the other hand, the inclusion of the threshold resummation leads to a surprising agreement between of the NNLL resummed with the fixed-order result, although a noticeable difference persists between the NLL and NLO results for larger values of p T . Perhaps the most interesting feature of the combined resummation is the fact that even in the absence of matching the resummed computation seems to capture the behaviour of the fixed-order in the region where standard small-p T fails to give accurate predictions. As a side note, unlike in standard small-p T resummation, consistent resummation does not require any regulation procedure to get rid of its large-p T behaviour. Instead, it relies on the fact that its b = 0 limit coincides with the threshold resummed inclusive cross section that fixes the integral of the distribution.
Next, in order to fully obtain the most accurate predictions, we match the resummed expressions to the fixed-order results according to Eq. (2.38). The results are presented on the right-hand side of Fig. 5. The upper panels show the comparison between NLL+LO (blue), NNLL+NLO (red), and NLO (orange) for the three types of resummations. The lower panels show the ratio w.r.t. the central scale (i.e. µ R = µ F = m H ) of NNLL+NLO. By comparing Fig. 5b with Fig. 5d, we see faster convergence for the consistent resummation. However, we notice that in the absence of the additional threshold resummed expression, the NLL+LO consistent resummation underestimates corrections from missing higher orders. This can be seen in Fig. 5d where the NLL+LO band shifts away from both the NLO and NNLL+NLO results. At small p T , consistent and combined resummation are exactly similar. The difference only occurs in the medium-range-p T where the contribution from the soft resummation again improves the agreement of the combined with the fixed-order result. In particular, we would like to highlight the good agreement between the resummed and the fixed-order results at large p T . In the context of consistent and combined results, the matched results have broadly similar features to the purely resummed ones. The main new feature is that the matching of the combined resummation changes the result less than the matching of the consistent resummation, and the latter less than the matching of the SSpT resummation.
We now turn to the application of the consistent and combined resummation to a low mass Z-boson production at the Tevatron via DY mechanism. We consider the production of a Z boson at m Z = 91.187 GeV at the Tevatron with a centre of mass energy √ s = 1.96 TeV because soft contributions for DY processes are expected to be more pronounced at colliders with lower center of mass energy. Fig. 6 shows the results in comparison to the standard small-p T resummation. Similar to Higgs results in Fig. 5, the pure resummed results are shown on the left-hand side. In all three types of resummations, one sees similar features as previously, namely the fact that the distributions go to zero as p T gets very small. The distributions then peak at about p T ∼ 2 GeV before vanishing again at large p T . Both at NLL and NNLL, consistent resummation (Fig. 6c) is analogous to standard transverse momentum (Fig. 6a) at low p T ; small but yet noticeable differences occur at large p T values, where consistent resummation displays slightly better convergence. This feature is more accentuated in the combined resummation (Fig. 6e) yielding a good agreement with the NLO (orange) results. Again, these results suggest the validity of the consistent and combined resummation even at large values of p T where the use of only standard p T resummation is not justified.    This behaviour remains when we match the pure resummed calculation to the fixedorder. The results are shown on the right-hand side of Fig. 6. In particular, we can notice the good agreement between the NLL+LO, NNLL+NLO, and NLO results in the combined case (Fig. 6f). Not only is the NNLL+NLO band smaller compared to the NNLL+LO, but the latter is also contained in the former suggesting a good convergence of the resummed perturbative expansion. These results suggest that for the transverse momentum distribution of the Z boson where x is often far from unity, the effect of the threshold resummation is less pronounced.

Conclusion
In this paper we have presented a phenomenological application of the combined threshold and consistent resummation introduced in ref. [16] for transverse momentum distributions. This combined resummation allows for an improvement of the transverse momentum distribution that holds for all values of momenta.
We have presented in Section 2 various analytical frameworks needed for the phenomenology studies. A particular attention was focused on the implementation of the Borel method as an alternative prescription to perform the inverse Fourier and Mellin transform. It was shown that the Borel prescription owes good numerical stability while providing complete control of the power correction terms, introduced by the consistent resummation procedure, allowing us to reproduce standard small-p T results in the limit p T → 0. The results of our phenomenological studies were presented in Section 3 where we matched our resummed expression with fixed-order calculations. We found that the effects of the consistent resummation have a modest significance for a low mass Z boson production via DY mechanism at small and moderate p T while sizeable improvement can be seen in the tail of the distribution. The difference, however, is rather significant in the case of Higgs boson produced via gluon fusion, which highlights the potential relevance of threshold resummation to gluon-induced processes [46][47][48][49]. As a main result we found that while consistent resummation enhances the convergence at small-p T , its combination with threshold resummation improves the agreement with fixed-order results in the intermediate and large-p T regions.
As possible directions for future development of this work, we will extend the result in ref. [16] to N3LL and interface our implementation with Monte Carlo codes in order to produce N3LL+NNLO predictions. In addition, it will be interesting to study further the relation between the consistent and threshold resummation in order to identify the missing soft logarithms with an ultimate goal of developing a combined expression that does not rely on a profile matching function. This, in turn, can be extended to account for the highenergy (or small-x) resummation, which quite recently has been jointly performed with small-p T [50] and threshold resummation [32].