1 Introduction

The current precision-physics program at the Large Hadron Collider (LHC) requires Standard Model (SM) theoretical predictions at the highest accuracy. Data belonging to “benchmark” processes, which are measured with the utmost precision at the LHC, need to be tested against theoretical results at the same level of accuracy. This is not only important for the extraction of SM parameters per se, but also for searches of signals of new physics, that can appear as small deviations in kinematic distributions with respect to the SM predictions. Reaching the highest possible level of precision is then the main goal and the calculation of perturbative QCD corrections plays a dominant role in this context.

Until a few years ago, the standard for such calculations was next-to-leading order (NLO) accuracy. In recent years, a continuously-growing number of next-to-next-to-leading order (NNLO) results for many important processes has appeared in the literature, giving birth to the so called “NNLO revolution”. For several “standard candles” processes, the first steps towards the calculation of differential cross sections at \(\hbox {N}^3\hbox {LO}\) have also been taken (see e.g. [1, 2]).

The computation of higher-order terms in the perturbative series becomes more involved due to the technical difficulties arising in the evaluation of virtual contributions and to the increasing complexity of the infrared (IR) structure of the real contributions. In order to expose the cancellation of the IR divergences between real and virtual contributions, the knowledge of the behaviour of the scattering amplitudes at the boundaries of the phase space is then a crucial ingredient and it is indeed what is used by the subtraction methods in order to work. These methods can be roughly divided into local and slicing. Among the first, the most extensively used at NLO were proposed in Refs. [3, 4]. As far as the NNLO subtraction methods are concerned, the past few years have witnessed a great activity in their development: the transverse-momentum (\(q_{\mathrm{T}}\)) subtraction method [5,6,7,8], the N-jettiness subtraction [9, 10], the projection-to-Born [11], the residue subtraction [12, 13] and the antenna subtraction method [14,15,16] have all been successfully applied to LHC phenomenology. The first application of the \(q_{\mathrm{T}}\)-subtraction method to differential cross sections at \(\hbox {N}^3\)LO was recently proposed in Ref. [1], in the calculation of the rapidity distribution of the Higgs boson.

While a local subtraction is independent of any regularising parameter, slicing methods require the use of a cutoff to separate the different IR regions. Such separation of the phase space introduces instabilities in the numerical evaluation of cross sections and differential distributions [17,18,19,20], and some care has to be taken in order to obtain stable and reliable results.

The knowledge of logarithmic and power-correction terms in the cutoff plays a relevant role in the identification of universal structures, in the development of regularisation prescriptions and in resummation programs [21,22,23,24,25,26,27,28,29,30,31]. According to their behaviour in the zero limit, the cutoff-dependent terms can be classified into logarithmically-divergent or finite contributions, and power-correction terms, that vanish in that limit. In particular, the terms that are singular in the small-cutoff limit are universal and are cancelled by the application of the subtraction methods, while finite and vanishing terms are, in general, process dependent. However, after the subtraction procedure, a residual dependence on the cutoff remains as power corrections. While these terms formally vanish in the null cutoff limit, they give a non-zero numerical contribution for any finite choice of the cutoff.

From a theoretical point of view, the knowledge of the power corrections greatly increases our understanding of the perturbative behaviour of the QCD cross sections, since more non-trivial (universal and non-universal) terms appear. The origin of these terms can be traced back both to the scattering amplitudes, evaluated at phase-space boundaries, and to the phase space itself. Thus, several papers have tackled the study of power corrections in the soft and collinear limits [32,33,34], while studies in the general framework of fixed-order and threshold-resummed computations have also been performed [35,36,37,38,39,40,41,42].

From a practical point of view, the knowledge of the power corrections makes the numerical implementation of a subtraction method more robust, since the power terms weaken the dependence of the final result on the arbitrary cutoff. This is not only valid when the subtraction method is applied to NLO computations, but it is numerically more relevant when applied to higher-order calculation, as pointed out, for example, in the evaluation of NNLO cross sections in Refs. [19, 20].

Power corrections at NLO have been extensively studied in Refs. [43,44,45,46,47,48,49,50,51,52] in the context of the N-jettiness subtraction method, and in Refs. [53,54,55,56,57,58] within SCET-based subtraction methods. A numerical extraction of power corrections in the context of NNLL’+NNLO calculations was done in N-jettiness [31], and a general discussion in the context of the fixed-order implementation of the N-jettiness subtraction can be found in Ref. [10].

In this paper we consider the production of a colourless system at next-to-leading order in the strong coupling constant \(\alpha _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{S}}\). In particular, we discuss Drell–Yan (DY) V production and Higgs boson production in gluon fusion at NLO, in the infinite top-mass limit. We impose a transverse-momentum cutoff, \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\), on the colourless system and we compute the power corrections in the cutoff, up to order \((q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}})^4\), for the inclusive cross sections. The knowledge of these terms will shed light upon the non-trivial behaviour of cross sections at the boundaries of the phase space, and upon the resummation structure at subleading orders. In addition, it allows to have a better control on the cutoff-dependent terms, when the \(q_{\mathrm{T}}\)-subtraction method of Ref. [5] is applied to the numerical calculation of cross sections, allowing for a use of larger values of the cutoff. We also describe a process-independent procedure that can be used to compute the all-order power corrections in the cutoff.

The outline of this paper is as follows. In Sect. 2 we introduce our notation, and we briefly summarize the expressions of the partonic and hadronic cross sections, in a form that is suitable for what follows. In Sect. 3 we outline the calculation we have done and in Sect. 4 we present and discuss our analytic results for V and H production, along with a study of their numerical impact. We draw our conclusions in Sect. 5. We leave to the appendixes all the technical details of our calculation.

2 Kinematics and notation

We briefly introduce the notation used in our theoretical framework and we recall some kinematic details of the calculations presented in this paper.

2.1 Hadronic cross sections

We consider the production of a colourless system F of squared invariant mass \(Q^2\) plus a coloured system X at a hadron collider

$$\begin{aligned} h_1+h_2 \rightarrow F+X . \end{aligned}$$
(2.1)

We call S the hadronic squared center-of-mass energy and we write the hadronic differential cross section for this process as

$$\begin{aligned} d\sigma = \sum _{a,b} \int _\tau ^1 dx_1 \int _\frac{\tau }{x_1}^1 dx_2\, f_a\!\left( x_1\right) f_b\!\left( x_2\right) d\hat{\sigma }_{ab}, \end{aligned}$$
(2.2)

where

$$\begin{aligned} \tau = \frac{Q^2}{S}, \end{aligned}$$
(2.3)

\(f_{a/b}\) are the parton densities of the partons a and b, in the hadron \(h_1\) and \(h_2\) respectively, and \(d\hat{\sigma }_{ab}\) is the partonic cross section for the process \(a+b \rightarrow F + X\). The dependence on the renormalisation and factorisation scales and on the other kinematic invariants of the process are implicitly assumed.

In Appendix A we have collected all the formulae for the calculation of the partonic cross section. Using Eqs. (A.19) and (A.20), we write the hadronic cross section as

$$\begin{aligned} \sigma= & {} \sum _{a,b} \int _\tau ^1 dx_1 \int _\frac{\tau }{x_1}^1 dx_2\, f_a\!\left( x_1\right) f_b\!\left( x_2\right) \nonumber \\&\times \int dq_{\mathrm{T}}^2\, dz \,\frac{d\hat{\sigma }_{ab}\!\left( q_{\mathrm{T}},z\right) }{dq_{\mathrm{T}}^2}\, \delta \!\left( \! z- \frac{Q^2}{s}\right) , \end{aligned}$$
(2.4)

where s is the partonic center-of-mass energy, equal to

$$\begin{aligned} s = S \, x_1\, x_2. \end{aligned}$$
(2.5)

We have also made explicit the dependence on z, the ratio between the squared invariant mass of the system F and the partonic center-of-mass energy, and on \(q_{\mathrm{T}}\), the transverse momentum of the system F with respect to the hadronic beams. Using Eqs. (2.5) and (2.3) and integrating over \(x_2\) we obtain

$$\begin{aligned} \sigma= & {} \sum _{a,b} \tau \int _\tau ^1 \frac{dz}{z} \int _\frac{\tau }{z}^1 \frac{dx_1}{x_1}\, f_a\!\left( x_1\right) f_b\!\left( \frac{\tau }{z\,x_1}\right) \nonumber \\&\times \frac{1}{z} \int \! dq_{\mathrm{T}}^2 \, \frac{d\hat{\sigma }_{ab}\!\left( q_{\mathrm{T}},z\right) }{dq_{\mathrm{T}}^2}. \end{aligned}$$
(2.6)

We then introduce the parton luminosity \(\mathcal{L}_{ab}(y)\) defined by

$$\begin{aligned} \mathcal{L}_{ab}(y) \equiv \int _y^1 \frac{dx}{x}\, f_a\!\left( x\right) f_b\! \left( \frac{y}{x}\right) , \end{aligned}$$
(2.7)

so that we can finally write

$$\begin{aligned} \sigma = \sum _{a,b} \tau \int _\tau ^1 \frac{dz}{z} \, \mathcal{L}_{ab}\!\left( \frac{\tau }{z}\right) \frac{1}{z} \int \! dq_{\mathrm{T}}^2 \, \frac{d\hat{\sigma }_{ab}\!\left( q_{\mathrm{T}},z\right) }{dq_{\mathrm{T}}^2}. \end{aligned}$$
(2.8)

2.2 Partonic differential cross sections

In this section we recall the formulae for the first-order real corrections to the Drell–Yan production of a weak boson V (W or Z) and to the Higgs boson production in gluon fusion, in the infinite top-mass limit.

The partonic cross sections \(d\hat{\sigma }_{ab}(q_{\mathrm{T}},z)/dq_{\mathrm{T}}^2\) in Eq. (2.8) are computable in perturbative QCD as power series in \(\alpha _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{S}}\)

$$\begin{aligned} \frac{d\hat{\sigma }_{ab}\!\left( q_{\mathrm{T}},z\right) }{dq_{\mathrm{T}}^2}= \frac{d\hat{\sigma }^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2} + \frac{\alpha _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{S}}}{2\pi } \, \frac{d\hat{\sigma }_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2} + \cdots \nonumber \\ \end{aligned}$$
(2.9)

The Born contribution \({d\hat{\sigma }^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)}(q_{\mathrm{T}},z)}/{dq_{\mathrm{T}}^2}\) and the virtual contributions to \({d\hat{\sigma }_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(q_{\mathrm{T}},z)}/{dq_{\mathrm{T}}^2}\) are proportional to \(\delta (q_{\mathrm{T}})\).

Applying the formulae detailed in Appendix A, with a little abuse of notation,Footnote 1 we can write the partonic differential cross sections for the real corrections to V and H production as

2.2.1 \({\varvec{V}}\) production

  • \(q(\bar{q}) + g \rightarrow V + q(\bar{q})\)

    $$\begin{aligned} \frac{d\hat{\sigma }_{qg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2}= & {} \sigma _{qq}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \,T_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{R}}\, z\,\nonumber \\&\times \frac{z\,(1+3z)\,{\displaystyle \frac{q_{\mathrm{T}}^2}{Q^2}}+(1-z)\, p_{qg}(z)}{\sqrt{(1-z)^2-4z\,{\displaystyle \frac{q_{\mathrm{T}}^2}{Q^2}}}}\,\frac{1}{q_{\mathrm{T}}^2},\nonumber \\ \end{aligned}$$
    (2.10)
  • \(q + \bar{q} \rightarrow V + g\)

    $$\begin{aligned} \frac{d\hat{\sigma }_{q{\bar{q}}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2}= & {} \sigma _{qq}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \, C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{F}}\, z\, \nonumber \\&\times \frac{-\,4z\,{\displaystyle \frac{q_{\mathrm{T}}^2}{Q^2}}+2\,(1-z)\,\hat{p}_{qq}(z)}{\sqrt{(1-z)^2-4z\,{\displaystyle \frac{q_{\mathrm{T}}^2}{Q^2}}}}\,\frac{1}{q_{\mathrm{T}}^2},\nonumber \\ \end{aligned}$$
    (2.11)

where

$$\begin{aligned} \sigma _{q q}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} = \frac{\pi }{N_\mathrm{c}}\,\frac{g^2 \left( g_v^2 + g_a^2\right) }{c_\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } W}^2} \frac{1}{Q^2} \end{aligned}$$
(2.12)

is the Born-level cross section for the process \(q \bar{q} \rightarrow V\), with \(c_\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } W}\) the cosine of the weak angle and g, \(g_v\), \(g_a\) the weak, the vector and the axial coupling, respectively. With \(N_\mathrm{c}\) we denote the number of colours, \(C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{F}}=\left( N_\mathrm{c}^2-1\right) /2N_\mathrm{c}=4/3\) and \(T_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{R}}=1/2\). For W production, the flavours of the quarks in Eqs. (2.10)–(2.12) are different, and the corresponding Cabibbo–Kobayashi–Maskawa matrix element has to be included in Eq. (2.12). The expression of the Altarelli–Parisi splitting functions \(p_{qg}(z)\) and \(\hat{p}_{qq}(z)\) are given in Appendix D.

2.2.2 \({\varvec{H}}\) production

  • \(g + q({\bar{q}}) \rightarrow H + q({\bar{q}})\)

    $$\begin{aligned} \frac{d\hat{\sigma }_{gq}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2}= & {} \sigma _{gg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \,C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{F}}\, z\,\nonumber \\&\times \frac{-\,3\,(1-z)\,{\displaystyle \frac{q_{\mathrm{T}}^2}{Q^2}}+(1-z)\,p_{gq}(z)}{\sqrt{(1-z)^2-4z\,{\displaystyle \frac{q_{\mathrm{T}}^2}{Q^2}}}} \,\frac{1}{q_{\mathrm{T}}^2},\nonumber \\ \end{aligned}$$
    (2.13)
  • \(g + g \rightarrow H + g\)

    $$\begin{aligned} \frac{d\hat{\sigma }_{gg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2}= & {} \sigma _{gg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \,C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{A}}\, z\,\nonumber \\&\times \frac{4z \left( {\displaystyle \frac{q_{\mathrm{T}}^2}{Q^2}}\right) ^{\!\!2} - 8\,(1-z)^2\,{\displaystyle \frac{q_{\mathrm{T}}^2}{Q^2}}+2\,(1-z)\,\hat{p}_{gg}(z)}{\sqrt{(1-z)^2 - 4 z\,{\displaystyle \frac{q_{\mathrm{T}}^2}{Q^2}}}}\,\frac{1}{q_{\mathrm{T}}^2},\nonumber \\ \end{aligned}$$
    (2.14)

where

$$\begin{aligned} \sigma _{gg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} = \frac{\alpha _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{S}}^2}{72\pi } \,\frac{1}{ N_\mathrm{c}^2-1}\, \frac{1}{v^2} \end{aligned}$$
(2.15)

is the Born-level cross section for the process \(gg \rightarrow H\), with v the Higgs vacuum expectation value and \(C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{A}}=N_\mathrm{c}=3\). The expression of the Altarelli–Parisi splitting functions \(p_{gq}(z)\) and \(\hat{p}_{gg}(z)\) are given in Appendix D.

We notice that the terms proportional to the Altarelli–Parisi splitting functions in Eqs. (2.10)–(2.14) embody in a single expression the whole infrared behaviour of the amplitudes, i.e. their soft and collinear limits. The structure of these terms was derived in a completely general form, from the universal behaviour of the scattering amplitudes in those limits, in Ref. [59].

3 Description of the calculation

In the small-\(q_{\mathrm{T}}\) region, i.e. \(q_{\mathrm{T}}\ll Q\), the real contribution to the perturbative cross sections of Eqs. (2.10)–(2.14) contains well-known logarithmically-enhanced terms that are singular in the \(q_{\mathrm{T}}\rightarrow 0\) limit [21,22,23,24,25,26,27,28,29,30]. In the context of inclusive NLO fixed-order calculations, the logarithmic terms are cancelled when using the subtraction prescriptions. For more exclusive quantities, such as the transverse-momentum distribution of the colourless system, the same logarithmic terms need to be resummed at all orders in the strong coupling constant to produce reliable results. Although our studies are of value in the context of the transverse-momentum resummation, here we limit ourselves to the case of inclusive fixed-order predictions at NLO, leaving the resummation program to future investigations. In this paper we compute power-correction terms to the cross section that, although vanishing in the small-\(q_{\mathrm{T}}\) limit, may give a sizable numerical contribution when using a slicing subtraction method.

To explicitly present the perturbative structure of these terms at small \(q_{\mathrm{T}}\), it is customary in the literature [47, 59] to compute the following cumulative partonic cross section, integrating the differential cross section in the range \(0 \le q_{\mathrm{T}}\le q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\),

$$\begin{aligned} \hat{\sigma }_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } <}(z) \equiv \int _0^{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\right) ^2}\!\! dq_{\mathrm{T}}^2 \,\frac{d\hat{\sigma }_{ab}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2}. \end{aligned}$$
(3.1)

The cross section in Eq. (3.1) receives contributions from the Born and the virtual terms, both proportional to \(\delta (q_{\mathrm{T}})\), and from the part of the real amplitude that describes the production of the F system with transverse momentum less than \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\). The virtual and real contributions are separately divergent and are typically regularised in dimensional regularisation. Since the total partonic cross section is finite and analytically known for the processes under study, following what was done in Refs. [60, 61], we compute the above integral as

$$\begin{aligned} \hat{\sigma }_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } <}(z) = \hat{\sigma }^\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } tot}_{ab}(z) - \hat{\sigma }_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } >}(z), \end{aligned}$$
(3.2)

with

$$\begin{aligned} \hat{\sigma }^\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } tot}_{ab}(z)= & {} \int _0^{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{max}}\right) ^2}\!\! dq_{\mathrm{T}}^2 \,\frac{d\hat{\sigma }_{ab}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2}, \end{aligned}$$
(3.3)
$$\begin{aligned} \hat{\sigma }_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } >}(z)= & {} \int _{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\right) ^2}^{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{max}}\right) ^2}\!\! dq_{\mathrm{T}}^2\,\frac{d\hat{\sigma }_{ab}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2} , \end{aligned}$$
(3.4)

where \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{max}}\) is the maximum transverse momentum allowed by the kinematics, \(\hat{\sigma }^\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } tot }_{ab}(z)\) is the total partonic cross section and \(\hat{\sigma }_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } >}(z)\) is the partonic cross section integrated above \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\). The advantage of using Eq. (3.2) is that the partonic cross section integrated in the range \(0\le q_{\mathrm{T}}\le q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) is obtained as difference of the total cross section (formally free from any dependence on \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\)) and the partonic cross section integrated in the range above \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) of Eq. (3.4). Since \(q_{\mathrm{T}}>q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}> 0\), the last integration can be performed in four space-time dimensions, with no further use of dimensional regularisation. In Refs. [60, 61] the computation of the cumulative cross section was performed in the limit \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\ll Q\), neglecting terms of \({\mathcal {O}}\!\left( (q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}})^2\right) \) on the right-hand side of Eq. (3.2). In this paper, we compute these terms up to \({\mathcal {O}}\!\left( (q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}})^4\right) \) included.

3.1 \(\varvec{q_{\mathrm{T}}}\)-integrated partonic cross sections

In this section we present the results for the partonic cross section in Eq. (3.4), integrated in \(q_{\mathrm{T}}\), from an arbitrary value \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) up to the maximum transverse momentum \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{max}}\) allowed by the kinematics of the event, given by

$$\begin{aligned} \left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{max}}\right) ^2 = Q^2\, \frac{(1-z)^2}{4\,z}, \end{aligned}$$
(3.5)

at a fixed value of z. The integrations are straightforward and do not need any dedicated comment. To lighten up the notation, we introduce the dimensionless quantityFootnote 2

$$\begin{aligned} a \equiv \frac{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\right) ^2}{Q^2}, \end{aligned}$$
(3.6)

that will be our expansion parameter in the rest of the paper, and we define

$$\begin{aligned} \pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2 \equiv \frac{4az}{(1-z)^{2}}, \end{aligned}$$
(3.7)

that will allow us to write the upcoming differential cross sections in a more compact form.

3.1.1 \({\varvec{V}}\) production

  • \(q(\bar{q}) + g \rightarrow V + q(\bar{q})\)

    $$\begin{aligned} \hat{\sigma }_{qg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } > (1)}(z)= & {} \int _{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\right) ^2}^{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{max}}\right) ^2}\!\!dq_{\mathrm{T}}^2 \, \frac{d\hat{\sigma }_{qg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2} \nonumber \\= & {} \sigma _{qq}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \,T_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{R}}\, z\left\{ \frac{1}{2}(1+3z)(1-z)\,\right. \nonumber \\&\times \sqrt{1-\frac{4az}{(1-z)^{2}}} \nonumber \\&+\,p_{qg}(z)\! \left[ -\log {\frac{az}{(1-z)^{2}}}\right. \nonumber \\&\left. \left. +\,2\log \frac{1}{2} \left( \sqrt{1-\frac{4az}{(1-z)^2}}+1 \right) \right] \right\} \nonumber \\= & {} \sigma _{qq}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \,T_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{R}}\, z\left\{ \frac{1}{2}(1+3z)(1-z)\,\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2} \right. \nonumber \\&\left. + p_{qg}(z) \log \frac{1+\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}}{1-\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}} \right\} \!, \end{aligned}$$
    (3.8)
  • \(q + \bar{q} \rightarrow V + g\)

    $$\begin{aligned} \hat{\sigma }_{q{\bar{q}}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } >(1)}(z)= & {} \int _{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\right) ^2}^{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{max}}\right) ^2}\!\! \frac{d\hat{\sigma }_{q\bar{q}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2} \nonumber \\= & {} \sigma _{qq}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \, C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{F}}\, z \left\{ - 2\, (1-z)\,\right. \nonumber \\&\times \sqrt{1-\frac{4az}{(1-z)^{2}}} \nonumber \\&+ 2 \,\hat{p}_{qq}(z)\! \left[ - \log \frac{az}{(1-z)^{2}} \right. \nonumber \\&\left. \left. + 2 \log \frac{1}{2} \left( \sqrt{1-\frac{4az}{(1-z)^2}}+1 \right) \right] \right\} \nonumber \\= & {} \sigma _{qq}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \, C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{F}}\, z \left\{ - 2\, (1-z)\,\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2} \right. \nonumber \\&\left. + 2 \,\hat{p}_{qq}(z) \log \frac{1+\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}}{1-\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}} \right\} \! . \end{aligned}$$
    (3.9)

3.1.2 \({\varvec{H}}\) production

  • \(g + q({\bar{q}}) \rightarrow H + q({\bar{q}})\)

    $$\begin{aligned} \hat{\sigma }_{gq}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } >(1)}(z)= & {} \int _{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\right) ^2}^{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{max}}\right) ^2}\!\! \frac{d\hat{\sigma }_{gq}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2} \nonumber \\= & {} \sigma _{gg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \,C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{F}}\, z \left\{ -\frac{3(1-z)^2}{2z}\,\right. \nonumber \\&\times \sqrt{1-\frac{4az}{(1-z)^2}} \nonumber \\&+ p_{gq}(z) \! \left[ - \log \frac{az}{(1-z)^2} \right. \nonumber \\&\left. \left. + 2 \log \frac{1}{2} \left( \sqrt{1-\frac{4az}{(1-z)^2}}+1 \right) \right] \right\} \nonumber \\= & {} \sigma _{gg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \,C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{F}}\, z \left\{ -\frac{3(1-z)^2}{2z}\,\sqrt{1- \pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2} \right. \nonumber \\&\left. + p_{gq}(z) \log \frac{1+\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}}{1-\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}} \right\} \!, \end{aligned}$$
    (3.10)
  • \(g + g \rightarrow H + g\)

    $$\begin{aligned} \hat{\sigma }_{gg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } >(1)}(z)= & {} \int _{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\right) ^2}^{\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{max}}\right) ^2}\!\! \frac{d\hat{\sigma }_{gg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(q_{\mathrm{T}},z)}{dq_{\mathrm{T}}^2} \nonumber \\= & {} \sigma _{gg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \,C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{A}}\, z \left\{ -\frac{4(1-z)^{3}}{z}\,\right. \nonumber \\&\times \sqrt{1-\frac{4az}{(1-z)^{2}}} \nonumber \\&+ 4 \,z \left[ \frac{1-z}{2z}\,\frac{(1-z)^{2}}{6z}\, \sqrt{1-\frac{4az}{(1-z)^{2}}}\,\right. \nonumber \\&\times \left. \left( 1+\frac{2az}{(1-z)^{2}}\right) \right] \nonumber \\&+\, 2 \, \hat{p}_{gg}(z) \left[ -\log {\frac{az}{(1-z)^{2}}}\right. \nonumber \\&\left. \left. +\,2\log {\frac{1}{2} \left( \sqrt{1-\frac{4az}{(1-z)^{2}}}+1 \right) } \right] \right\} \nonumber \\= & {} \sigma _{gg}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \,C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{A}}\, z\! \left\{ -\frac{11}{3} \frac{(1-z)^3}{z}\right. \nonumber \\&\times \left( 1 - \frac{\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}{22}\right) \sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2} \nonumber \\&\left. + 2 \, \hat{p}_{gg}(z) \log \frac{1+\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}}{1-\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}} \right\} . \end{aligned}$$
    (3.11)

We do not consider the process \(q\bar{q} \rightarrow H g\) since it is not singular in the limit \(q_{\mathrm{T}}\rightarrow 0\) and the corresponding analytic/numeric integration in the transverse momentum can be performed setting explicitly \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}= 0\).

A couple of further comments about the above expressions are also in order. In first place, the part of the cross sections proportional to the Altarelli–Parisi splitting functions in Eqs. (3.8)–(3.11) has a universal origin, due to the factorisation of the collinear singularities on the underlying Born. The rest of the above cross sections is, in general, not universal. In addition, for Higgs boson production, the NLO cumulative cross sections that we have computed coincide exactly with the jet-vetoed cross sections \(\sigma ^\mathrm{veto}(p_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^\mathrm{veto})\) of Ref. [62], provided we identify \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}=p_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^\mathrm{veto}\).

3.2 Extending the integration in \({\varvec{z}}\)

According to Eq. (2.8), in order to compute the hadronic cross section we need to integrate the partonic cross sections convoluted with the corresponding luminosities. In the calculation of the total cross sections, the upper limit in the z integration is unrestricted and is equal to 1. When a cut on the transverse momentum \(q_{\mathrm{T}}\) is applied, the reality of Eqs. (3.8)–(3.11) imposes the non negativity of the argument of the square roots, i.e.

$$\begin{aligned} 1- \pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2\ge 0, \end{aligned}$$
(3.12)

that in turn gives

$$\begin{aligned} z \le z^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{max}}\equiv & {} 1-f(a), \nonumber \\ f(a)\equiv & {} 2 \sqrt{a}\left( \sqrt{1 + a} - \sqrt{a}\right) . \end{aligned}$$
(3.13)

Since our aim is to make contact with the transverse-momentum subtraction formulae, that describe the behaviour of the cross sections in the soft and collinear limits, we need to extend the integration range of the z variable up to 1, i.e. the upper integration limit of z in a Born-like kinematics. In fact, only in the \(z \rightarrow 1\) limit we recover the logarithmic structure from the soft region of the emission. In order to obtain explicitly all the logarithmic-enhanced terms in the small-\(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) limit, we have then to expand our results in powers of a. Since both the integrand and the upper limit of the integral depend on a, the naïve approach of expanding only the integrand does not work, due to the appearance of divergent terms in the \(z \rightarrow 1\) limit, that have to be handled with the introduction of plus distributions.

Using the notation of Ref. [60], we first introduce the function \(\hat{R}_{ab}(z)\), defined by

$$\begin{aligned} \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle }<}= & {} \tau \int _\tau ^{1-f(a)} \frac{dz}{z} \, \mathcal{L}_{ab}\!\left( \frac{\tau }{z}\right) \frac{1}{z} \,\hat{\sigma }_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } <}(z) \nonumber \\\equiv & {} \tau \int _\tau ^1 \frac{dz}{z} \, \mathcal{L}_{ab}\!\left( \frac{\tau }{z}\right) \hat{\sigma }^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \hat{R}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } }(z), \end{aligned}$$
(3.14)

where the upper integration limit in z in the last integral is exactly 1 and \(\hat{\sigma }^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)}\) is the partonic Born-level cross section for the production of the colourless system F. The function \(\hat{R}_{ab}(z)\) can be written as a perturbative expansion in \(\alpha _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{S}}\)

$$\begin{aligned} \hat{R}_{ab}(z) = \delta _\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } B}\,\delta (1-z) + \sum _{n=1}^\infty \left( \frac{\alpha _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{S}}}{2\pi }\right) ^n {\hat{R}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (n)}(z), \end{aligned}$$
(3.15)

where the \(\delta (1-z)\) term is the Born-level contribution, and \(\delta _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm B}=1\) when partons a and b are such that \(a+b \rightarrow F\) is a possible Born-like process, otherwise its value is 0.

The coefficient functions \({\hat{R}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (n)}(z)\) can be computed as power series in a. It is in fact well known in the literature [60] that the NLO coefficient \({\hat{R}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\) has the following formFootnote 3

$$\begin{aligned} \hat{R}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)= & {} \log ^2(a) \, {\hat{R}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,2,0)}(z) + \log (a) \, {\hat{R}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,1,0)}(z) \nonumber \\&+ {\hat{R}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,0,0)}(z) + {\mathcal {O}}\!\left( a^\frac{1}{2} \log {a}\right) , \end{aligned}$$
(3.16)

and the aim of this paper is to compute the first unknown terms in Eq. (3.16) that were neglected in Refs. [60] and [61], namely \(R_{ab}^{(1,m,r)}(z)\), for r up to 4 and for any m.

In a way similar to what was done in Eq. (3.14) for \(\hat{R}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } }(z)\), we introduce the function \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } }(z)\) defined by

$$\begin{aligned} \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle }>}= & {} \tau \int _\tau ^{1-f(a)} \frac{dz}{z} \, \mathcal{L}_{ab}\!\left( \frac{\tau }{z}\right) \frac{1}{z} \,\hat{\sigma }_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } >}(z) \nonumber \\\equiv & {} \tau \int _\tau ^1 \frac{dz}{z} \, \mathcal{L}_{ab}\!\left( \frac{\tau }{z}\right) \hat{\sigma }^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \, \hat{G}_{ab}(z). \end{aligned}$$
(3.17)

Since

$$\begin{aligned} \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } <} + \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } >} = \tau \int _\tau ^1 \frac{dz}{z} \, \mathcal{L}_{ab}\!\left( \frac{\tau }{z}\right) \hat{\sigma }^\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } tot}_{ab}(z) \equiv \sigma ^\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } tot}_{ab}, \end{aligned}$$
(3.18)

and \(\sigma ^\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } tot}_{ab}\) is independent of a, the coefficients of the terms that vanish in the small-\(q_{\mathrm{T}}\) limit in the series expansion in a of \(\hat{R}_{ab}(z)\) and \(\hat{G}_{ab}(z)\) are equal but with opposite sign, at any order in \(\alpha _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{S}}\). We recall that \(\hat{R}_{ab}(z)\) contains terms of the form \(\delta (1-z)\), coming from the Born and the virtual contributions, that are independent of a and are obviously absent in \(\hat{G}_{ab}(z)\).

In the rest of the paper we compute the first terms of the expansion in a of \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\), that will be obtained from the following identity

$$\begin{aligned} \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle }> (1)}= & {} \tau \int _\tau ^{1-f(a)} \frac{dz}{z} \, \mathcal{L}_{ab}\!\left( \frac{\tau }{z}\right) \frac{1}{z} \,\hat{\sigma }_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } >(1)}(z)\nonumber \\= & {} \tau \int _\tau ^1 \frac{dz}{z} \, \mathcal{L}_{ab}\!\left( \frac{\tau }{z}\right) \hat{\sigma }^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \, \hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z) . \end{aligned}$$
(3.19)

We have elaborated a process-independent formula to transform an integral of the form of the first one in Eq. (3.19) into the form of the second one, producing the series expansion of \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\) in a. The application of our formula reorganizes the divergent terms in the \(z \rightarrow 1\) limit into terms that are integrable up to \(z=1\) and logarithmic terms in a. Since this is a very technical procedure, we have collected all the details in Appendix B, and we refer the interested reader to that Appendix for the description of the method.

4 Results

In this section we summarize our findings. We present in Sect. 4.1 the analytic results for the \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\) functions we have computed. In the calculation of these functions, we kept trace of all the terms originating from the manipulation of the contributions proportional to the Altarelli–Parisi splitting functions, in the partonic cross sections of Eqs. (2.10)–(2.14). These terms constitute what we call the “universal part” of our results, as detailed in Sects. 2.2 and 3.1. We will indicate these terms with the superscript “U”, while the remaining terms will have a superscript “R”. We stress here that the distinction between universal and non-universal part is purely formal, and it does not have a physical implication. The reason of this separation is to have hints on the general structure of the \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) dependence of inclusive cross sections for the production of arbitrary colorless systems. We comment on the results that we have obtained in Sect. 4.2.

In Sect. 4.3 we study the numerical significance of the power-correction terms we have computed, discussing first their impact on the different production channels for Drell–Yan Z boson and Higgs boson production in gluon fusion. Then we present their overall effect, normalising the results with respect to the total NLO cross section, in order to have a better grasp on the size of these contributions.

4.1 Results for the \(\varvec{\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)}\) functions

We indicate with \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm U (1)}_{ab}(z)\) the universal part of the \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\) functions, and with \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{ab}(z)\) the remaining part, stripped off of a common colour factor. Our expressions for \(\hat{G}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{ab}(z)\) contain derivatives of the Dirac \(\delta \) function, \(\delta ^{(n)}(z)\), up to \(n=5\), and plus distributions up to order 5. We report here the definition of a plus distribution of order n

$$\begin{aligned}&\int _0^1 dz\, l(z) \left[ g(z) \right] _{n+} \nonumber \\&\quad \equiv \int _0^1 dz \left\{ l(z) - \sum _{i=0}^{n-1} \frac{1}{i!} \, l^{(i)}(1)\, (z-1)^i \right\} g(z) , \end{aligned}$$
(4.1)

where g(z) has a pole of order n for \(z=1\), and l(z) is a continuous function in \(z=1\), together with all its derivatives up to order \((n-1)\). For completeness, we collect in Appendix E more details on the plus distributions, and the identities we have used to simplify our results.

4.2 \({\varvec{V}}\) production

  • \(q(\bar{q}) {+} g \rightarrow V + q(\bar{q})\)

    $$\begin{aligned} \hat{G}^{ {\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{qg}(z){=} T_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{R}}\, \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{qg}(z), \quad \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{qg}(z) = \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm U (1)}_{qg}(z) + \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{qg}(z),\nonumber \\ \end{aligned}$$
    (4.2)

    where

    $$\begin{aligned} \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm U (1)}_{qg}(z)= & {} - p_{qg}(z) \log (a) - p_{qg}(z)\, \log \frac{z}{(1-z)^{2}} \nonumber \\&+\, \left\{ \delta ^{(1)}(1 - z) - 3\,\delta (1 - z) \right\} a \log (a) \nonumber \\&+\, \left\{ \delta (1 - z) - 2 \, z\, p_{qg}(z) \left[ \frac{1}{(1-z)^2}\right] _{2+} \right\} a \nonumber \\&+\, \left\{ - 9\, \delta (1 - z) + \frac{21}{2}\, \delta ^{(1)}(1 - z)\right. \nonumber \\&\left. -\, 3 \,\delta ^{(2)}(1 - z) + \frac{1}{4} \,\delta ^{(3)}(1 - z) \! \right\} \! a^2 \log (a) \nonumber \\&+\, \left\{ + 2 \,\delta (1 - z) + \frac{7}{4} \,\delta ^{(1)}(1 - z)\right. \nonumber \\&\left. -\, \frac{5}{4} \,\delta ^{(2)}(1 - z) + \frac{1}{6}\, \delta ^{(3)}(1 - z) \right. \nonumber \\&\left. -\, 3 \, z^2 p_{qg}(z) \left[ \frac{1}{(1-z)^4} \right] _{4+} \right\} a^2 \nonumber \\&+\, {\mathcal {O}}\!\left( a^\frac{5}{2}\log (a)\right) , \end{aligned}$$
    (4.3)
    $$\begin{aligned} \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{qg}(z)= & {} \frac{1}{2}\, (1+3z) (1-z) - z (1+3z) \left[ \frac{1}{1-z} \right] _+ a\nonumber \\&-z^2 (1+3z) \left[ \frac{1}{(1-z)^3}\right] _{3+} \!\! a^2 \nonumber \\&+\, 2\, \delta (1 - z) \, a \log (a) - 2 \, \delta (1 - z) \, a \nonumber \\&+\, \left\{ 5 \,\delta (1 - z) - \frac{11}{2}\, \delta ^{(1)}(1 - z) \right. \nonumber \\&\left. + \delta ^{(2)}(1 - z) \right\} a^2 \log (a) \nonumber \\&+\, \left\{ - 2\, \delta (1 - z) - \frac{3}{4}\, \delta ^{(1)}(1 - z)\right. \nonumber \\&\left. +\, \frac{1}{2}\, \delta ^{(2)}(1 - z) \right\} a^2 + {\mathcal {O}}\!\left( a^\frac{5}{2}\log (a)\right) ,\nonumber \\ \end{aligned}$$
    (4.4)

    and

    $$\begin{aligned} p_{qg}(z) = 2z^2-2z+1 . \end{aligned}$$
    (4.5)
  • \(q + \bar{q} \rightarrow V + g\)

    $$\begin{aligned} \hat{G}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{q\bar{q}}(z) = C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{F}}\, \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{q\bar{q}}(z), \quad \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{q\bar{q}}(z) = \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm U (1)}_{q\bar{q}}(z)+ \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{q\bar{q}}(z) ,\nonumber \\ \end{aligned}$$
    (4.6)

    where

    $$\begin{aligned} \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm U (1)}_{q\bar{q}}(z)= & {} \delta (1 - z) \log ^2(a) - 2 \, p_{qq}(z) \log (a) \nonumber \\&- \frac{\pi ^2}{3} \delta (1 - z) - 2 \, \hat{p}_{qq}(z) \log (z)\, \nonumber \\&+ 4 \, (1-z) \, \hat{p}_{qq}(z) \left[ \frac{\log (1-z)}{1-z} \right] _+ \nonumber \\&+ \left\{ 6\, \delta (1 - z) - 8\, \delta ^{(1)}(1 - z) \right. \nonumber \\&\left. + 2 \,\delta ^{(2)}(1 - z) \right\} a \log (a) \nonumber \\&+ \left\{ - 6\, \delta (1 - z) + 4 \,\delta ^{(1)}(1 - z) \right. \nonumber \\&\left. -4 \,z \, (1-z) \, \hat{p}_{qq}(z) \left[ \frac{1}{(1-z)^3}\right] _{3+} \right\} a \nonumber \\&+ \left\{ 3 \, \delta (1 - z) - 12 \, \delta ^{(1)}(1 - z) \right. \nonumber \\&\left. + \frac{21}{2}\, \delta ^{(2)}(1 - z) - 3 \, \delta ^{(3)}(1 - z) \right. \nonumber \\&\left. + \frac{1}{4}\, \delta ^{(4)}(1 - z) \!\right\} a^2 \log (a) \nonumber \\&+ \left\{ - 4\, \delta (1 - z) + 6 \,\delta ^{(1)}(1 - z) \right. \nonumber \\&- \frac{1}{2}\,\delta ^{(2)}(1 - z) - \delta ^{(3)}(1 - z) \nonumber \\&+ \frac{1}{6}\, \delta ^{(4)}(1 - z) - 6 z^2 (1-z)\, \nonumber \\&\times \left. \hat{p}_{qq}(z) \left[ \frac{1}{(1-z)^5}\right] _{5+} \right\} a^2 \nonumber \\&+ {\mathcal {O}}\!\left( a^{\frac{5}{2}}\log (a)\right) , \end{aligned}$$
    (4.7)
    $$\begin{aligned} \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{q\bar{q}}(z)= & {} - 2 \, (1-z) + 4\,z \left[ \frac{1}{1-z} \right] _+ a\nonumber \\&+ 4 \, z^2 \, \left[ \frac{1}{(1-z)^3}\right] _{3+}\!\! a^2 \nonumber \\&- 2 \,\delta (1 - z) \, a \log (a) + 2 \,\delta (1 - z) \, a \nonumber \\&+ \left\{ - 2\, \delta (1 - z) + 4 \, \delta ^{(1)}(1 - z) \right. \nonumber \\&\left. - \delta ^{(2)}(1 - z) \right\} a^2 \log (a) \nonumber \\&+ \left\{ 2 \,\delta (1 - z) - \frac{1}{2}\, \delta ^{(2)}(1 - z) \right\} a^2 \nonumber \\&+ {\mathcal {O}}\!\left( a^{\frac{5}{2}}\log (a)\right) , \end{aligned}$$
    (4.8)

    and

    $$\begin{aligned} \hat{p}_{qq}(z) = \frac{1+z^2}{1-z}, \quad p_{qq}(z)= \frac{1+z^2}{(1-z)_+} . \end{aligned}$$
    (4.9)

    In Eq. (4.7) we have written the \((1+z^2)\) terms coming from the numerator of the \( \hat{p}_{qq}(z)\) splitting function as

    $$\begin{aligned} 1+z^2 = (1-z) \, \hat{p}_{qq}(z), \end{aligned}$$
    (4.10)

    in order to keep track of the universal origin of those terms.

4.3 \({\varvec{H}}\) production

  • \(g + q({\bar{q}}) \rightarrow H + q({\bar{q}})\)

    $$\begin{aligned} \hat{G}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{gq}(z) = C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{F}}\, \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{gq}(z), \quad \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{gq}(z) = \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm U (1)}_{gq}(z) + \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{gq}(z),\nonumber \\ \end{aligned}$$
    (4.11)

    where

    $$\begin{aligned} \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm U (1)}_{gq}(z)= & {} - p_{gq}(z) \log (a) - p_{gq}(z)\log \frac{z}{(1-z)^{2}} \nonumber \\&+ \delta ^{(1)}(1 - z)\, a \log (a) + \left\{ \delta (1 - z) \right. \nonumber \\&\left. - 2\, z \, p_{gq}(z) \left[ \frac{1}{(1-z)^2}\right] _{2+}\right\} a \nonumber \\&+ \left\{ - \frac{3}{2} \,\delta (1 - z) + \frac{3}{2}\, \delta ^{(1)}(1 - z) \right. \nonumber \\&\left. - \frac{3}{4}\, \delta ^{(2)}(1 - z) + \frac{1}{4} \,\delta ^{(3)}(1 - z) \right\} a^2 \log (a) \nonumber \\&+ \left\{ \frac{1}{4} \,\delta (1 - z) + \frac{1}{4} \,\delta ^{(1)}(1 - z) \right. \nonumber \\&\left. + \frac{1}{4} \,\delta ^{(2)}(1 - z)+ \frac{1}{6} \,\delta ^{(3)}(1 - z) \right. \nonumber \\&\left. - 3 \, z^2 \, p_{gq}(z) \left[ \frac{1}{(1-z)^4} \right] _{4+} \right\} a^2 \nonumber \\&+ {\mathcal {O}}\!\left( a^\frac{5}{2}\log (a)\right) , \end{aligned}$$
    (4.12)
    $$\begin{aligned} \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{gq}(z)= & {} - \frac{3}{2z} (1-z)^2 + 3\, a + 3 \, z \left[ \frac{1}{(1-z)^2}\right] _{2+} \!\! a^2 \nonumber \\&+ \frac{3}{2} \left\{ \delta (1 - z) - \delta ^{(1)}(1 - z) \right\} a^2 \log (a) \nonumber \\&- \frac{3}{4} \left\{ \delta (1 - z) + \delta ^{(1)}(1 - z) \right\} a^2 \nonumber \\&+ {\mathcal {O}}\!\left( a^\frac{5}{2}\log (a)\right) , \end{aligned}$$
    (4.13)

    and

    $$\begin{aligned} p_{gq}(z) = \frac{z^2-2z+2}{z} . \end{aligned}$$
    (4.14)
  • \(g + g \rightarrow H + g\)

    $$\begin{aligned} \hat{G}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{gg}(z)= & {} C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{A}}\, \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{gg}(z), \nonumber \\ \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{gg}(z)= & {} \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm U (1)}_{gg}(z) + \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{gg}(z), \end{aligned}$$
    (4.15)

    where

    $$\begin{aligned} \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm U (1)}_{gg}(z)= & {} \delta (1 - z) \log ^2(a) -2 \, p_{gg}(z) \log (a) \nonumber \\&- \frac{\pi ^2}{3} \delta (1 - z) - 2 \, \hat{p}_{gg}(z) \log (z)\nonumber \\&+ 4 \, (1-z) \, \hat{p}_{gg}(z) \left[ \frac{\log (1-z)}{1-z} \right] _+ \nonumber \\&+ \Big \{ 12 \,\delta (1 - z) - 8 \, \delta ^{(1)}(1 - z) \nonumber \\&+ 2 \,\delta ^{(2)}(1 - z) \Big \}\, a \log (a) \nonumber \\&+ \left\{ - 6\, \delta (1 - z) + 4\, \delta ^{(1)}(1 - z) \right. \nonumber \\&\left. \,- 4 z \,(1-z) \,\hat{p}_{gg}(z) \left[ \frac{1}{(1-z)^3}\right] _{3+} \right\} a \nonumber \\&+ \left\{ 18 \,\delta (1 - z) - 30\,\delta ^{(1)}(1 - z) \right. \nonumber \\&+ 15\, \delta ^{(2)}(1 - z) - 3 \,\delta ^{(3)}(1 - z) \nonumber \\&\left. + \frac{1}{4}\, \delta ^{(4)}(1 - z) \right\} a^2 \log (a) \nonumber \\&+ \left\{ - \frac{15}{2}\, \delta (1 - z) + 3 \,\delta ^{(1)}(1 - z)\right. \nonumber \\&\left. + \frac{5}{2} \,\delta ^{(2)}(1 - z) - \delta ^{(3)}(1 - z) \right. \nonumber \\&+ \frac{1}{6} \,\delta ^{(4)}(1 - z) - 6\, z^2 \, (1-z) \nonumber \\&\times \left. \hat{p}_{gg}(z) \left[ \frac{1}{(1-z)^5}\right] _{5+} \right\} a^2 \nonumber \\&+ {\mathcal {O}}\!\left( a^\frac{5}{2}\log (a)\right) , \end{aligned}$$
    (4.16)
    $$\begin{aligned} \hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{gg}(z)= & {} - \frac{11}{3 z} \, (1-z)^3 + 8\,(1-z) \, a\nonumber \\&+ 6 \, z \left[ \frac{1}{1-z} \right] _+ \!\! a^2 \nonumber \\&- 3\, \delta (1 - z) \, a^2 \log (a) \nonumber \\&- \frac{5}{2} \, \delta (1 - z) \, a^2 + {\mathcal {O}}\!\left( a^\frac{5}{2}\log (a)\right) , \end{aligned}$$
    (4.17)

    and

    $$\begin{aligned} \hat{p}_{gg}(z)= \frac{2 (z^2-z+1)^2}{z(1-z)}. \end{aligned}$$
    (4.18)

    In Eq. (4.16) we have written the \(2 (z^2-z+1)^2/z\) terms coming from the numerator of the \(\hat{p}_{gg}(z)\) splitting function as

    $$\begin{aligned} \frac{2 (z^2-z+1)^2}{z} = (1-z) \, \hat{p}_{gg}(z), \end{aligned}$$
    (4.19)

    in order to keep track of the universal origin of those terms.

4.4 Comments on \(\varvec{\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)}\)

The leading-logarithmic (LL) and next-to-leading-logarithmic (NLL) coefficients of the \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\) functions that we have computed agree with the ones in the literature, along with the finite term. Their values have been known for a while [26, 29] and are related to the perturbative coefficients of the transverse-momentum subtraction/resummation formulae for V [28] and Higgs boson production [63], as pointed out in Sect. 3.2. The coefficients of the terms of order \(a\log (a)\) and a, and of order \(a^2\log (a)\) and \(a^2\) are instead the new results computed in this paper.

The general form of the \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\) functions we have computed readsFootnote 4

$$\begin{aligned} \hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)= & {} \log ^2(a) \, {\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,2,0)}(z) \nonumber \\&+\, \log (a) \, {\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,1,0)}(z) + {\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,0,0)}(z) \nonumber \\&+ \, a\log (a) \,{\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,1,2)}(z) + a \, {\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,0,2)}(z) \nonumber \\&+ \, a^2 \log (a) \,{\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,1,4)}(z) + a^2 \, {\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,0,4)}(z) \nonumber \\&+\, {\mathcal {O}}\!\left( a^\frac{5}{2}\log (a)\right) , \end{aligned}$$
(4.20)

all the other coefficients being zero.

We will refer to the terms in the first and in the second line of Eq. (4.20) as leading terms (LT). These terms are either logarithmically divergent or finite in the \(a\rightarrow 0\) limit. We name the terms in the sum in the third line of Eq. (4.20) as next-to-leading terms (NLT), and the terms in the fourth line as next-to-next-to-leading terms (\(\hbox {N}^2\)LT), and so forth.

We notice that the NLT and \(\hbox {N}^2\)LT terms are at most linearly dependent on \(\log (a)\), consistently with the fact that the LL contribution is a squared logarithm. In addition, no odd-power corrections of \(\sqrt{a} = q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}/Q\) appear in the NLT and \(\hbox {N}^2\)LT terms. This behaviour is in agreement with what found, for example, in Ref. [48], i.e. that, at NLO, the power expansion of the differential cross section for colour-singlet production is in \((q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}})^2\). We do not expect this to be true in general when cuts are applied to the final state.

4.4.1 Soft behaviour of the universal part

The origin of some of the terms in the diagonal channels, i.e. the \(q\bar{q}\) channel for V production and the gg channel for H production, can be traced back to the behaviour of the Altarelli–Parisi splitting functions in the soft limit, i.e. \(z\rightarrow 1\). In fact, in this limit,

$$\begin{aligned} \hat{P}_{qq}(z) \approx \frac{2C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{F}}}{1-z}, \quad \hat{P}_{gg}(z) \approx \frac{2C_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{A}}}{1-z}, \end{aligned}$$
(4.21)

so that

$$\begin{aligned} \hat{p}_{qq}(z) \approx \hat{p}_{gg}(z) \approx \frac{2}{1-z}\equiv \hat{p}(z). \end{aligned}$$
(4.22)

Inserting \(\hat{p}_{qq}(z)\) and \(\hat{p}_{gg}(z)\) in Eqs. (3.9) and (3.11), respectively, they give rise to a contribution of the form

$$\begin{aligned}&\int _\tau ^{1-f(a)} \frac{dz}{z} \, {\mathcal {L}}\left( \frac{\tau }{z}\right) \, 2\, \hat{p}(z)\, \log \frac{1+\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}}{1-\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}} \nonumber \\&\quad = \int _0^1 \frac{dz}{z} \,\mathcal{L}\left( \frac{\tau }{z}\right) \left\{ \delta (1 - z) \log ^2(a) \right. \nonumber \\&\qquad \left. - 2\, p(z) \log (a) - \frac{\pi ^2}{3} \, \delta (1 - z) + \cdots \right\} \end{aligned}$$
(4.23)

where, following the notation of Appendix D, we have defined

$$\begin{aligned} p(z) = \left[ \frac{2}{1-z}\right] _+ . \end{aligned}$$
(4.24)

The details for the derivation of Eq. (4.23) are collected in Appendix C.5. Inspecting the first three terms of the universal function \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm U (1)}_{q\bar{q}}(z)\) in Eq. (4.7) and \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm U (1)}_{gg}(z)\) in Eq. (4.16), we recognize exactly the three terms on the right-hand side of Eq. (4.23).

4.4.2 The non-universal part

It is also interesting to notice that the non-universal part of the \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\) functions contains terms proportional to \(\log (a)\), multiplied by powers of a. These powers are controlled by the form of the non-universal parts in Eqs. (3.8)–(3.11), and to the way they enter in our generating procedure described in Appendix B. In fact, by inspecting Eq. (B.10), we see that they contribute to \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{ab}\) with terms of the form

$$\begin{aligned}&\int _\tau ^{1-f(a)} dz \, (1-z)^n \sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2} \nonumber \\&\quad = \left\{ \begin{array}{ll} + a^2 \log (a) + a\log (a) + \cdots \qquad &{} n=1\\ -2\, a^2 \log (a) + \cdots \qquad &{} n=2\\ + a^2 \log (a) + \cdots \qquad &{} n=3\\ -6\,a^3 \log (a) + \cdots \qquad &{} n=4\\ +2\,a^3 \log (a) + \cdots \qquad &{} n=5 \end{array} \right. \end{aligned}$$
(4.25)

where the dots stand for power terms in a with no logarithms attached. This also explains why, for V production, \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{qg}(z)\) in Eq. (4.4) and \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{q\bar{q}}(z)\) in Eq. (4.8) contain both terms \(a\log (a)\) and \(a^2\log (a)\): they receive contributions from all the terms in Eq. (4.25) starting from \(n=1\), since Eqs. (3.8) and (3.9) contain a term proportional to \((1-z)\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}\). Instead, \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{gq}(z)\) in Eq. (4.13) and \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm R (1)}_{gg}(z)\) in v (4.17) contain only the term \(a^2\log (a)\), since they receive contributions from the terms in Eq. (4.25) starting from \(n=2\), due to the fact that Eqs. (3.10) and (3.11) contain a term proportional to \((1-z)^2\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}\) and \((1-z)^3\sqrt{1-\pi _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm T}^2}\), respectively.

As far as the finite term in the diagonal channels is concerned, we notice that, in the \(q\bar{q}\) channel of DY production, the first term in Eq. (4.8) happens to correspond to the first-order collinear coefficient function defined in the “hard-resummation scheme”, introduced in Refs. [64, 65] within the \(q_{\mathrm{T}}\)-subtraction formalism. Instead, the first term in the gg channel of H production in Eq. (4.17) has no connection with the first-order collinear coefficient function, that is zero for this production channel. In conclusion, the structure of the terms in the non-universal part depends on the peculiar form of the differential cross sections.

4.4.3 Higher-order soft behaviour of the squared amplitudes

In this section, we extend the study performed in Sect. 4.2.1 in order to investigate the origin of the power-suppressed terms \(a \log (a)\) and \(a^2 \log (a)\), present both in the universal and in the non-universal parts. We will show that their origin can be connected to the higher-order soft behaviour of the squared amplitudes. To this aim, we have performed a Laurent expansion in the energy \(k^0\) of the final-state parton of the exact squared amplitudes of Eqs. (A.23), (A.29), (A.34) and (A.39). In the following, we call leading soft (LS) the term proportional to the highest negative power of \(k^0\), next-to-leading soft (\(\hbox {N}^1\)LS) the subsequent term, and so on. All the technical details and expressions of the expansion terms are collected in Appendix A.3.

We have then applied the algorithm described in this paper to each of the terms of the expansions that we have calculated, in order to compute their behaviour as a function of a. This has allowed us to trace the origin of the \(a\log (a)\) and \(a^2 \log (a)\) terms. Our findings are collected in the following:

  • \(q(\bar{q}) + g \rightarrow V + q(\bar{q})\)

    We reproduce the \(a\log (a)\) behaviour of \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{qg}(z)\) in Eq. (4.2) if we consider the soft-expansion of the exact amplitude up to the \(\hbox {N}^1\)LS level, i.e. if we sum Eqs. (A.24) and (A.25), and the \(a^2\log (a)\) behaviour if we consider the soft-expansion up to the \(\hbox {N}^3\)LS level, i.e. if we sum Eqs. (A.24)–(A.27).

  • \(q + \bar{q} \rightarrow V + g\)

    We reproduce the \(a\log (a)\) behaviour of \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{q\bar{q} }(z)\) in Eq. (4.6) if we consider the soft-expansion of the exact amplitude up to the \(\hbox {N}^1\)LS level, i.e. if we sum Eqs. (A.30) and (A.31), and the \(a^2\log (a)\) behaviour if we consider the soft-expansion up to the \(\hbox {N}^2\)LS level, i.e. if we sum Eqs. (A.30)–(A.32).

  • \(g + q({\bar{q}}) \rightarrow H + q({\bar{q}})\)

    We reproduce the \(a\log (a)\) behaviour of \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{gq}(z)\) in Eq. (4.11) if we consider the soft-expansion of the exact amplitude up to the \(\hbox {N}^1\)LS level, i.e. if we sum Eqs. (A.35) and (A.36), and the \(a^2\log (a)\) behaviour if we consider the soft-expansion up to the \(\hbox {N}^2\)LS level, i.e. if we sum Eqs. (A.35)–(A.37).

  • \(g + g \rightarrow H + g\)

    We reproduce the \(a\log (a)\) behaviour of \(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{gg}(z)\) in Eq. (4.15) if we consider the soft-expansion of the exact amplitude up to the \(\hbox {N}^2\)LS level, i.e. if we sum Eqs. (A.40)–(A.42), and the \(a^2\log (a)\) behaviour if we consider the soft-expansion up to the \(\hbox {N}^4\)LS level, i.e. if we sum Eqs. (A.40)–(A.44).

Collecting our result in a table, we have:

 

V

H

\(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{qg}(z)\)

\(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{q\bar{q} }(z)\)

\(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{gq}(z)\)

\(\hat{g}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}_{gg}(z)\)

\(a\log (a)\)

\(\hbox {N}^1\)LS

\(\hbox {N}^1\)LS

\(\hbox {N}^1\)LS

\(\hbox {N}^2\)LS

\(a^2 \log (a)\)

\(\hbox {N}^3\)LS

\(\hbox {N}^2\)LS

\(\hbox {N}^2\)LS

\(\hbox {N}^4\)LS

In summary, the next-to-leading-soft approximation of the exact amplitudes reproduces the \(a\log (a)\) term only for V production and for the non-diagonal channel of H production. For the diagonal channel of H production, only the expansion up to next-to-next-to-leading-soft order reproduces the \(a\log (a)\) term. We would like to point out that, of the three terms contributing to the \(\hbox {N}^2\)LS of Eq. (A.42), only the constant one, i.e. the number 16, is needed to reproduce the \(a\log (a)\) coefficient. The u / t and t / u terms do not give rise to any \(a\log (a)\) contribution.

Fig. 1
figure 1

Difference of the total cross sections \(\left( \sigma ^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } > (1)} - \widetilde{\sigma }\right) \) as a function of \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\), for Z boson production, in the \(qg \rightarrow Z q\) (left pane) and in the \(q\bar{q} \rightarrow Z g\) channel (right pane). The three curves correspond to the three possible choices of \(\widetilde{\sigma }\): results for \(\widetilde{\sigma } = \sigma ^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm LT }\) are displayed in blue, for \(\widetilde{\sigma } = \sigma ^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm LT }+ \sigma ^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm NLT }\) are displayed in black and for \(\widetilde{\sigma } = \sigma ^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm LT }+ \sigma ^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm NLT }+ \sigma ^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm N^2LT }\) are displayed in red. The corresponding values of \(a=\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}/m_\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } Z}\right) ^2\) are displayed on the top of the figure. The statistical errors of the integration are also shown, but they are totally negligible on the scale of the figure

Fig. 2
figure 2

Difference of the total cross sections \(\left( \sigma ^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } > (1)} - \widetilde{\sigma }\right) \) as a function of \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\), for H boson production, in the \(qg \rightarrow H q\) (left pane) and in the \(gg \rightarrow H g\) channel (right pane). Same legend as in Fig. 1. The corresponding values of \(a=\left( q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}/m_\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } H}\right) ^2\) are displayed on the top of the figure

Moreover, only the expansion up to next-to-next-to-leading-soft order in the diagonal channel for V production and in the non-diagonal channel for H production reproduces the \(a^2\log (a)\) coefficient. Higher orders in the expansion in the softness of the final-state parton are needed for the non-diagonal channel of V production and for the diagonal channel of H production.

4.4.4 \(\varvec{q_{\mathrm{T}}}\)-subtraction method

In the original paper on the \(q_{\mathrm{T}}\)-subtraction method [5], the expansion in \(\alpha _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{S}}\) of the transverse-momentum resummation formula generates exactly the three terms in Eq. (3.16), plus extra power-correction terms.

In the formula for \(\hat{R}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\) that we can build from our expression of \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\), by changing the overall sign and adding the \(\delta (1-z)\) contribution from the virtual correction, the power-correction terms are exactly those produced by the expansion of the real amplitudes. If one is interested in using our formula for \(\hat{R}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\) to reduce the dependence on the transverse-momentum cutoff, within the \({q_{\mathrm{T}}}\)-subtraction method, the aforementioned extra terms need then to be subtracted from our expression of \(\hat{R}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\).

Fig. 3
figure 3

Results for \(1 - \left( \sigma ^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } > (1)} - \widetilde{\sigma }\right) /\sigma _\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } NLO}\) as a function of \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\), for Z boson production, in \(pp\rightarrow Zj\). Same legend as in Fig. 1. In the left pane, the low-\(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) region is displayed, while, in the right pane, a larger region in \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) is shown. The total cross section at NLO for Z production, \(\sigma _\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } NLO}\), has been taken equal to 55668.1 pb

Fig. 4
figure 4

Results for \(1 - \left( \sigma ^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } > (1)} - \widetilde{\sigma }\right) /\sigma _\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } NLO}\) as a function of \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\), for H boson production, in \(pp\rightarrow Hj\). Same legend as in Fig. 1. In the left pane, the low-\(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) region is displayed, while, in the right pane, a larger region in \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) is shown. The total cross section at NLO for H production, \(\sigma _\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } NLO}\), has been taken equal to 31.52 pb

4.5 Numerical results

As previously pointed out, NLO (and NNLO) cross sections computed with the \(q_{\mathrm{T}}\)-subtracted formalism exhibit a residual dependence on \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\), i.e. the parameter a we have introduced in Eq. (3.6). This residual dependence is due to power terms which remain after the subtraction of the IR singular contributions, and vanish only in the limit \(a\rightarrow 0\) (limit which is unattainable in a numerical computation). In this section we discuss the residual systematic dependence on \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) due to terms beyond LT, NLT and \(\hbox {N}^2\)LT accuracy.

We present our results for Z and H production at the LHC, at a center-of-mass energy of \(\sqrt{S}=13\) TeV. In our NLO calculations we have set the renormalisation and factorisation scales equal to the mass of the corresponding produced boson, and we have used the MSTW2008nlo parton-distribution function set [66]. The mass of the Z boson \(m_\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } Z}\) and of the Higgs boson \(m_\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } H}\) have been set to the values 91.1876 GeV and 125 GeV, respectively.

As an overall check of our calculation, we compared the results obtained with the analytically \(q_{\mathrm{T}}\)-integrated cross sections in Eqs. (3.8)–(3.11) with the numerically-integrated results computed with both the DYqT-v1.0 [67, 68] and HqT2.0 [6, 69] codes, and found an excellent agreement.

Then, in order to study the residual \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) dependence of the NLO cross sections for all the partonic subprocesses, we insert the expansion in Eq. (4.20) into Eq. (3.19), and we introduce the following definitions

$$\begin{aligned} \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm LT }\equiv & {} \tau \int _\tau ^1 \frac{dz}{z} \, \mathcal{L}_{ab}\!\left( \frac{\tau }{z}\right) \hat{\sigma }^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \nonumber \\&\times \left[ \log ^2(a) \, {\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,2,0)}(z) \right. \nonumber \\&\left. + \log (a) \, {\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,1,0)}(z) + {\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,0,0)}(z) \right] , \end{aligned}$$
(4.26)
$$\begin{aligned} \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm NLT }\equiv & {} \tau \int _\tau ^1 \frac{dz}{z} \, \mathcal{L}_{ab}\!\left( \frac{\tau }{z}\right) \hat{\sigma }^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \nonumber \\&\times \left[ a\log (a) \,{\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,1,2)}(z) + a \, {\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,0,2)}(z) \right] , \end{aligned}$$
(4.27)
$$\begin{aligned} \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm N^2LT }\equiv & {} \tau \int _\tau ^1 \frac{dz}{z} \, \mathcal{L}_{ab}\!\left( \frac{\tau }{z}\right) \hat{\sigma }^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (0)} \nonumber \\&\times \left[ a^2 \log (a) \,{\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,1,4)}(z) + a^2 \, {\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,0,4)}(z) \right] , \end{aligned}$$
(4.28)

where we have dropped the \(^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } >}\) and \(^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}\) superscripts for ease of notation, since there is no possibility of misunderstanding in this section, because we present only the NLO results we have computed for the \({\hat{G}}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1)}(z)\) functions.

Fig. 5
figure 5

Same as Fig. 3, but using only the universal part of \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,n,m)}(z)\) in computing the cross sections of Eqs. (4.26)–(4.28)

The \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,n,m)}(z)\) functions in Eqs. (4.26)–(4.28) contain plus distributions up to order 5 and to compute these cross sections we have first built interpolations of the luminosity functions \(\mathcal{L}_{ab}(y)\), defined in Eq. (2.7), for the channels that contribute to Z and H production at NLO. We have expanded the luminosity functions on the basis of the Chebyshev polynomials up to order 30. In this way, the computation of the derivatives of the luminosity functions can be performed in a fast and sound way.

In the forthcoming figures, we plot the following quantities as a function of \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) (the corresponding value of a is given on top of each figure):

  1. 1.

    \((\sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } > (1)} - \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm LT })\),

  2. 2.

    \((\sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } > (1)} - \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm LT }- \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm NLT })\),

  3. 3.

    \((\sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } > (1)} - \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm LT }- \sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm NLT }-\sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm N^2LT })\),

where \(\sigma _{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } > (1)}\) is the cumulative cross section defined on the left-hand side of Eq. (3.19), obtained by integrating the exact differential cross sections of Eqs. (3.8)–(3.11). We expect that, by adding higher-power terms in a, these differences tend to zero more and more quickly when \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\rightarrow 0\). And in fact, the results shown in the following figures confirm this behaviour.

We first present our findings separated according to the partonic production channels. In all the figures presented in this paper, the statistical errors of the integration procedure are also displayed, but they are always totally negligible on the scales of the figures.

In Fig. 1 we collect the results for the aforementioned cross-section differences, as a function of \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\), for the \(qg \rightarrow Zq\) (left) and \(q \bar{q} \rightarrow Z g\) (right) channels, and in Fig. 2 we collect similar results for the \(gq \rightarrow Hq\) (left) and \(gg\rightarrow Hg\) (right) channels. As expected, NLT and \(\hbox {N}^2\)LT contributions increase the accuracy of the expanded cross section, with respect to the exact one.

To give a more quantitative estimation of the power-suppressed corrections, we present results for the total hadronic cross section, normalised with respect to the corresponding exact NLO cross section \(\sigma _\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } NLO}\) (i.e. including also the virtual contributions). The results are shown in Figs. 3 and 4, where we have used \(\sigma _\mathrm{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } NLO}=55668.1\) pb for Z production and 31.52 pb for H production. On the left panes we plot results in a smaller \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) region, while, on the right panes, we extend the \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) interval to higher values.

These plots show exactly how the residual cutoff dependence of the cross sections changes when the \(q_{\mathrm{T}}\)-subtraction counterterm is corrected by the NLT and \(\hbox {N}^2\)LT power terms. For example, for Z production and for \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}=10\) GeV, corresponding to \(a=0.012\), the LT cross section gives an estimate of the exact cross section within the 5‰, that reduces to below the 1‰ when the NLT contribution is added and becomes less than 0.01‰ when also the \(\hbox {N}^2\)LT is present. For Higgs boson production, the residual cutoff dependence is even more pronounced: in fact, at \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}=10\) GeV, corresponding to \(a=0.0064\), the LT is precise within the 1% level. When the NLT is added, the precision reaches the 0.2‰, and is below 0.001‰ with the addition of the \(\hbox {N}^2\)LT.

Fig. 6
figure 6

Same as Fig. 5 but for Higgs boson production

An interesting question is to estimate the impact of the universal parts of the \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,n,m)}(z)\) functions, with respect to the non-universal ones. We have then computed the cross sections in Eqs. (4.26)–(4.28), taking into account only the universal parts of the \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,n,m)}(z)\) functions. Our results are displayed in Fig. 5, for Z production, and in Fig. 6, for H production.

Comparing these figures with the corresponding ones with the full \(\hat{G}_{ab}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } (1,n,m)}(z)\) functions, i.e. Figs. 3 and 4 , we see that the non-universal contributions play a crucial role for Z production, while their role is minor in Higgs boson production. This is due to the fact that the non-universal part in Eqs. (3.10) and (3.11) is suppressed by higher powers of \((1-z)\), with respect to the corresponding expression for Z production, in Eqs. (3.8) and (3.9), confirming the conclusions drawn in Sect. 4.2.2.

5 Conclusions

In this paper we considered the production of a colourless system at next-to-leading order in the strong coupling constant \(\alpha _{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{S}}\). We imposed a transverse-momentum cutoff, \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\), on the colour-singlet final state and we computed the power corrections for the inclusive cross section in the cutoff, up to the fourth power. Although we studied Drell–Yan vector boson production and Higgs boson production in gluon fusion, the procedure we followed is general and can be applied to other similar cases, up to any order in the powers of \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\).

We presented analytic results, reproducing the known logarithmic terms from collinear and soft regions of the phase space, along with the finite contribution, and adding new terms as power corrections in \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\). We found that the logarithmic terms in \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) show up at most linearly in the power-correction contributions, consistently with the fact that the LL contribution is a squared logarithm. In addition, no odd-power corrections in \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) appeared in our calculation, in agreement with known results in the literature for the NLO differential cross section in colour-singlet production. We do not expect this to be true in general when cuts are applied to the final state.

Along the calculation we kept track of the origin of the newly-computed terms, so that we were able to separate them into a universal part, and a part that depends on the process at stake. In particular we derived and identified the contribution to the universal part coming from soft radiation, present in the diagonal partonic channels for Z and H production. We could also explain some features about the presence of power-suppressed logarithmic terms, appearing in the non-universal part of the power corrections. Furthermore we showed that the knowledge of the squared amplitudes at the next-to-leading-soft approximation is not enough to predict the \((q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}})^2\log q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) behaviour of the power corrections. The same conclusion can be drawn for the knowledge of the next-to-next-to-leading-soft approximation in predicting the \((q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}})^4\log q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\) power correction.

We also studied the numerical impact of the power terms in the hadronic cross sections for Z and H production at the LHC at 13 TeV, both by keeping track of the different partonic production channels and by summing over all of them. We plotted the behaviours of the cross sections while adding more and more orders of the power-correction terms, as a function of \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}\), and comparing them with the exact cross sections. For example, in Z production and for a value of \(q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}}=10\) GeV, the sensitivity on the cutoff can be reduced from 5 to 0.01‰, when adding the \((q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}})^4\) contributions to the \((q_{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{T}}^{{\displaystyle } {\scriptstyle } {\scriptscriptstyle } {\scriptscriptstyle } \mathrm{cut}})^2\) ones. Higgs boson production suffers from a larger sensitivity on the cutoff, and the dependence goes from 1% to 0.2‰, when all the power corrections we computed are added. By performing the same numerical comparisons for just the universal part of the power corrections, we showed that the non-universal contributions play a crucial role for Z production, while their role is minor in Higgs boson production.

The knowledge of the power terms is crucial for understanding both the non-trivial behaviour of cross sections at the boundaries of the phase space, and the resummation structure at subleading orders. Within the \(q_{\mathrm{T}}\)-subtraction method, the knowledge of the power terms helps in reducing the cutoff dependence of the cross sections. While the application of the \(q_{\mathrm{T}}\)-subtraction method in NLO calculations is superseded by well-known local subtraction methods, at NNLO it still plays a major role, also in view of the fact that, as shown in Refs. [1, 19], the sensitivity to the numerical value of the cutoff increases at higher orders.