Hard scattering processes at high-energy colliders, such as the Large Hadron Collider (LHC), characterized by large scales of energy (M) transferred, allows us to probe the dynamics of fundamental interactions at short distances. In this regime, theoretical predictions for cross sections can be evaluated with perturbative techniques. In particular accurate results require the inclusion of the dominant effects from strong interactions through the calculation of the higher-order terms in Quantum Chromodynamics (QCD) as a series expansion in the coupling \(\alpha _S(M)\). In order to match the experimental kinematical cuts on the measured final states, it is essential to obtain predictions for fiducial cross sections and corresponding differential distributions.

The computation of higher-order QCD corrections at fully-differential level is complicated by the presence of infrared singularities at intermediate stage of the calculation which prevents a direct implementation of numerical techniques and enforce the use of an hybrid analytic and numerical approach. At the next-to-leading order (NLO) general subtraction algorithms, which exploit the universality properties of soft and collinear emissions in QCD, have been developed [1,2,3,4]. These methods have been successfully implemented in general purpose Monte Carlo programs which satisfy the needs for the analysis of experimental data. Beyond the NLO, a widely used extension of the subtraction method is the so called transverse-momentum (\(q_T\)) subtraction formalism originally proposed in Ref. [5]. In fact, thanks to its relative simplicity and generality, the method has been successfully applied to fully differential QCD calculations for several hard-scattering processes at the next-to-next-to-leading order (NNLO) (see Ref. [6] and references therein) and, more recently, also at the next-to-next-to-next-to-leading order N\(^3\)LO [7,8,9].

In the case of the \(q_T\) subtraction formalism, a source of numerical uncertainty which is particularly difficult to quantify in a robust way is due to the unphysical power corrections of the type \({\mathcal {O}}(q_T^{\mathrm{cut}}/M)\), where \(q_T^{\mathrm{cut}}\) is the technical parameter necessary to separate the \(q_T\) resolved and unresolved parton emissions. Power corrections ambiguities are particularly severe in the case of fiducial selection cuts which yield an acceptance that has a linear dependence on \(q_T^{\mathrm{cut}}\) [10,11,12]. In principle the effect of these perturbative power corrections can be reduced setting the value of the technical parameter \(q_T^{\mathrm{cut}}\) sufficiently small. However, very small values of \(q_T^{\mathrm{cut}}\) require an extremely precise numerical control of cross sections in the infrared region \(q_T\rightarrow 0\) which are typically very challenging and time consuming.

In the case of the production of colourless high mass systems, sufficiently inclusive cross sections (such as total cross section in absence of fiducial selection cuts) computed within the \(q_T\) subtraction method have a residual dependence on \(q_T^{\mathrm{cut}}\) of order \({\mathcal {O}}((q_T^{\mathrm{cut}}/M)^2)\) originated from the integration of the corresponding scattering amplitudes over the final state kinematics in the small \(q_T\) region [13, 14]. In the case of the production of coloured systems the residual dependence on \(q_T^{\mathrm{cut}}\) is linear even in absence of fiducial cuts [15,16,17]. In this case, in order to remove the linear power corrections in \(q_T^{\mathrm{cut}}\), it is necessary to take care of both the dependence related to the fiducial selection cuts and the dependence which is present at the inclusive level.

Linear fiducial power corrections have been connected with alternating sign factorial growth of the perturbative expansion in Ref. [18] and modifications of the selection cuts typically used in experimental analysis have been proposed in order to eliminate the linear dependence of the acceptance. In Ref. [19], an experimental procedure has been proposed to remove from cross section measurements the effect of selection cuts which are responsible for linear fiducial power corrections. In Ref. [11] it has been shown that linear fiducial power corrections can be consistently removed through the \(q_T\) resummation formalism, if the \(q_T\) recoil due to multi-parton emission is correctly taken into account. Nevertheless fixed-order calculations have a great relevance in precision physics at colliders: they can be computed (in principle) in a definite and unambiguous way, and they are an essential ingredient for all-order resummed predictions. Therefore the goal of having at disposal fixed-order calculations at high numerical accuracy is very relevant, regardless of the effective physical precision of such predictions.

In this paper we consider standard fixed-order calculations and we discuss a method to remove linear fiducial power corrections (FPC) within the \(q_T\) subtraction formalism. Our method, which is equivalent to the one proposed in Ref. [11] based on a Lorentz decomposition for hadronic and leptonic tensors, introduces an additional subtraction exploiting the recoil procedure of Ref. [20] and it allows us to obtain fiducial fixed-order predictions with a residual uncertainty of the type \({\mathcal {O}}((q_T^{\mathrm{cut}}/M)^2)\) which can be brought down at sub-permille level. We stress that in this paper we are dealing with the linear power corrections to the \(q_T\) subtraction formalism [5] by using the subtraction counterterms related to the specific \(q_T\) resummation formalism of Refs. [21, 22] and to the use of the class of recoil prescriptions which were introduced in Ref. [20]. Therefore the strict validity of our results and the novelty with respect to the previous results in the literature [8, 11] are related to the particular implementation of the subtraction procedure within the formalism of Refs. [20,21,22].

Our results turn out to be crucial in the case of the N\(^3\)LO extension of the \(q_T\) subtraction formalism [7, 9] where it is particularly challenging to obtain precise perturbative predictions for very small values of \(q_T^{\mathrm{cut}}\).

We consider explicitly the case of fiducial cross sections for Drell–Yan lepton pair production at the Large Hadron Collider up to the N\(^3\)LO in QCD. We have implemented our method within the DYTurbo [23] numerical program and we have obtained perturbative predictions which are in agreement at permille level with those obtained with local subtraction formalism at NLO and NNLO.

We consider the hard-scattering process

$$\begin{aligned} h_1(p_1) + h_2(p_2) \rightarrow \sum _i F_i(q_i)+X, \end{aligned}$$
(1)

where \(F_i\) denotes the (colourless) final states with momenta \(q_i\) produced by the colliding hadrons \(h_1\) and \(h_2\) which we collectively identify as the system F(q), with momentum \(q=\sum _i q_i\), invariant mass \(M=\sqrt{q^2}\) and transverse momentum \(q_T\).

We start from the master formula of the \(q_T\) subtraction formalism for the hadronic cross section [5]

$$\begin{aligned} d\sigma ^F=d\sigma ^{F}_{LO}\otimes {\mathcal {H}}^F +\left[ d\sigma ^{{\mathrm{F}}+{\mathrm{jets}}}-d\sigma ^{\mathrm{CT}}\right] , \end{aligned}$$
(2)

where \(\sigma ^{F}_{LO}\) is the Born level cross section, \({\mathcal {H}}^F(\alpha _S)\) is the process-dependent hard-collinear function [24, 25] with the following perturbative expansion

$$\begin{aligned} {\mathcal {H}}^F= & {} 1+ \frac{\alpha _S}{\pi } \,{\mathcal {H}}^{F\,(1)} + \left( \frac{\alpha _S}{\pi }\right) ^2 \,{\mathcal {H}}^{F\,(2)} + \left( \frac{\alpha _S}{\pi }\right) ^3 \,{\mathcal {H}}^{F\,(3)}\nonumber \\&+\sum _{n=4}^{\infty } \left( \frac{\alpha _S}{\pi }\right) ^n\,{\mathcal {H}}^{F\,(n)}, \end{aligned}$$
(3)

\(d\sigma ^{\mathrm{CT}}\) is the subtraction counter-term [21]

$$\begin{aligned} d\sigma ^{\mathrm{CT}} = d\sigma ^{F}_{LO} \otimes \Sigma ^{F}(q_T/M) d^2\mathbf{q_T} \end{aligned}$$
(4)

and the symbol \(\otimes \) stands for convolutions over momentum fractions and sum over flavour indices of the partons. The second term in the r.h.s. of Eq. (2), \(d\sigma ^{{\mathrm{F}}+{\mathrm{jets}}}\), is the differential cross section for the production of F(q) in association with jets and it has to be evaluated at the previous perturbative order with respect to \(d\sigma ^{\mathrm{F}}\). The subtraction counter-term \(d\sigma ^{\mathrm{CT}}\) has the same singular behaviour of \(d\sigma ^{{\mathrm{F}}+{\mathrm{jets}}}\) in the limit \(q_T \rightarrow 0\) which functional form is known from the \(q_T\) resummation formalism [21, 22].

The terms \(d\sigma ^{{\mathrm{F}}+{\mathrm{jets}}}\) and \(d\sigma ^{\mathrm{CT}}\) in Eq. (2) are separately divergent due to infrared singularities at \(q_T=0\) and a technical parameter \(q_T^{\mathrm{cut}}\) has to be introduced. For \(q_T^{\mathrm{cut}}>0\) the sum of the terms in the square bracket of Eq. (2) is IR finite (or, more precisely, integrable over \(q_T\)) and it should be evaluated in the limit \(q_T^{\mathrm{cut}}\rightarrow 0\) to obtain the “exact” (free form residual \(q_T^{\mathrm{cut}}\) dependence) value of the cross section. However for finite value of \(q_T^{\mathrm{cut}}\) the cross section in Eq. (2) contains power corrections \({\mathcal {O}}((q_T^{\mathrm{cut}}/M)^p)\), with \(p>0\) [10]. The exact value of the exponent p depends by the cuts on the final states which define the fiducial cross section

$$\begin{aligned} \sigma ^F_{\mathrm{fid}} = \int _{\mathrm{cuts}} d\sigma ^F, \end{aligned}$$
(5)

we thus have

$$\begin{aligned} \sigma ^F_{\mathrm{fid}}= & {} \int _{\mathrm{cuts}} d\sigma ^{F}_{LO}\otimes {\mathcal {H}}^F +\int _{\mathrm{cuts}} \left[ d\sigma ^{{\mathrm{F}} +{\mathrm{jets}}}_{q_T>q_T^{\mathrm{cut}}}-d\sigma ^{\mathrm{CT}}_{q_T>q_T^{\mathrm{cut}}}\right] \nonumber \\&+{\mathcal {O}}\left( (q_T^{\mathrm{cut}}/M)^p\right) . \end{aligned}$$
(6)

In Ref. [12] has been shown that, in the case of Drell–Yan lepton pair production, typical cuts on the transverse momenta and rapidities of the final state particles \(F_i(q_i)\) leads to linear power corrections (\(p=1\)), which corresponding systematic uncertainty may spoil the accuracy of fixed-order calculations within the \(q_T\) subtraction formalism.

Clearly the effect of perturbative power corrections \({\mathcal {O}}((q_T^{\mathrm{cut}}/M)^p)\) can be reduced setting the value of the technical parameter \(q_T^{\mathrm{cut}}\) sufficiently small. However, very small values of \(q_T^{\mathrm{cut}}\) lead to large cancellations among the terms in the square bracket of the r.h.s. of Eq. (2), which in turns give rise to larger numerical integration uncertainties. These cancellations are particularly challenging at NNLO and N\(^3\)LO where the precise knowledge of the fully differential calculations of F in association with jets at NLO and NNLO is respectively required. Eventually a trade-off between the systematical and statistical uncertainties of the computation have to be found and, more importantly, a robust systematic uncertainty to the missing perturbative power corrections has to be computed. The systematic uncertainty can be estimated by evaluating the cross section at different values of \(q_T^{\mathrm{cut}}\) and carrying out a \(q_T^{\mathrm{cut}}\rightarrow 0\) extrapolation [6]. This is not a trivial task due to the large numerical uncertainties associated to the \(q_T^{\mathrm{cut}}\rightarrow 0\) limit.

We now discuss the method which enable us to consistently remove the FPC within the \(q_T\) subtraction formalism thus leaving a quadratic residual uncertainty \({\mathcal {O}}(({q_T^{\mathrm{cut}}}/M)^2\)). The starting point is the observation that the FPC have a kinematical origin [10,11,12]. They are generated by the selection cuts on the final state particles and they are indeed absent in fixed-order [10, 14, 26, 27] or \(q_T\) resummed calculations inclusive over the final state F decay products [28] and also in the case of \(q_T\) resummation with fiducial cuts when the \(q_T\) recoil due to multi-parton emission is correctly taken into account [11]. According to the \(q_T\) resummation formalism of Refs. [21, 22] the fiducial cross section in Eq. (6) can be schematically written in the following form:

$$\begin{aligned} \widetilde{\sigma }_{\mathrm{fid}}^{F}= & {} \int _{\mathrm{cuts}} d\widetilde{\sigma }^{F}_{LO}\otimes {\mathcal {H}}^F \times S(q_T,M)\nonumber \\&+\int _{\mathrm{cuts}} \left[ d\sigma ^{{\mathrm{F}}+{\mathrm{jets}}}_{q_T>q_T^{\mathrm{cut}}} -d\widetilde{\sigma }^{\mathrm{CT}}_{ q_T>q_T^{\mathrm{cut}}}\right] +{\mathcal {O}}\left( (q_T^{\mathrm{cut}}/M)^2\right) ,\nonumber \\ \end{aligned}$$
(7)

where

$$\begin{aligned} S(q_T,M)=\int _0^\infty db \frac{b}{2} J_0(b q_T) \exp ({\mathcal {G}}(\alpha _S)) \end{aligned}$$
(8)

and

$$\begin{aligned} d\widetilde{\sigma }^{\mathrm{CT}} = d\widetilde{\sigma }^{F}_{LO} \otimes \Sigma ^{F}(q_T/M) d^2\mathbf {q_T}. \end{aligned}$$
(9)

The first term on the r.h.s. of Eq. (9) is the resummed component of the cross section which collects in the generalized form factor \(\exp {({\mathcal {G}}(\alpha _S))}\) and resums to all orders (in the Fourier–Bessel conjugated impact-parameter b space) the enhanced logarithmic corrections of the type \(\alpha _S^n\ln ^m (M^2/q_T^2)\) which are present in the transverse momentum distribution at small \(q_T\) [21]. The second term on the r.h.s. of Eq. (9) is the fixed-order finite component of the cross section and \(q_T^{\mathrm{cut}}\) represents the minimum value of \(q_T\) used to compute such term. In the resummed formula in Eq. (7) the underlying amplitude of the Born level cross-section \(d{\widetilde{\sigma }}^F_{LO}\), which enters also in the term \(d\widetilde{\sigma }^{\mathrm{CT}}\), differs from the corresponding quantity \(d{\hat{\sigma }}^{(0)}\) in Eqs. (2) and (4) for the fact that it is not evaluated with the leading order (LO) kinematics but following the prescription introduced in Appendix A of Ref. [20], which takes into account the \(q_T\) recoil originated from the multiple radiation of soft and collinear partons in a kinematically consistent way.

Fig. 1
figure 1

Fiducial cross section for the production of \(l^+ l^-\) pairs from \(Z/\gamma ^*\) decay at the LHC (\(\sqrt{s}=13\) TeV). NLO results with the \(q_T\) subtraction method (blue squared points) and \(q_T\) subtraction method without FPC (red circled points) at various values of \(q_T^{\mathrm{cut}}\), and with a local subtraction method (black line). Error bars indicate the statistical uncertainties from Monte Carlo numerical integration

Full size image

Exploiting the resummation formula in Eq. (7) we are thus able to construct the following modified \(q_T\) subtraction formula which is free from linear fiducial power corrections [8, 11]:

$$\begin{aligned} \sigma ^F_{\mathrm{fid}}= & {} \int _{\mathrm{cuts}} d\sigma ^{F}_{LO}\otimes {\mathcal {H}}^F +\int _{\mathrm{cuts}} \left[ d\sigma ^{{\mathrm{F}} +{\mathrm{jets}}}_{q_T>q_T^{\mathrm{cut}}} -d\sigma ^{\mathrm{CT}}_{ q_T>q_T^{\mathrm{cut}}}\right] \nonumber \\&+ \int _{\mathrm{cuts}} d\sigma ^{\mathrm{FPC}} +{\mathcal {O}}\left( (q_T^{\mathrm{cut}}/M)^2\right) , \end{aligned}$$
(10)

where

$$\begin{aligned} d\sigma ^{\mathrm{FPC}}=\left[ d\widetilde{\sigma }^{{CT}}_{q_T<q_T^{\mathrm{cut}}}-d\sigma ^{\mathrm{CT}}_{ q_T<q_T^{\mathrm{cut}}}\right] . \end{aligned}$$
(11)

The inclusion of the additional term \(d\sigma ^{\mathrm{FPC}}\) for \(q_T<q_T^{\mathrm{cut}}\) allows us to produce the correct behavior of the fiducial cross section up to quadratic power corrections in \(q_T^{\mathrm{cut}}\). The terms \(d\sigma ^{{CT}}\) and \(d\widetilde{\sigma }^{{CT}}\) differ for the fact that they are respectively evaluated with the LO (\(q_T = 0\)) and with the recoil (\(q_T \ne 0\)) kinematics. We note that the term \(d\sigma ^{\mathrm{FPC}}\) is universal (i.e. process independent) and it is IR finite (albeit the two terms on the r.h.s. of Eq. (11) are separately divergent). Furthermore the contribution of \(d\sigma ^{\mathrm{FPC}}\) can be treated as a local subtraction: the difference of the terms in Eq. (11) is evaluated pointwise at integrand level and therefore the integration for \(q_T<q_T^{\mathrm{cut}}\) can be extended up to (virtually) arbitrary small value of \(q_T\). In the current numerical implementation [29] we extended the \(q_T\) integration of the term \(d\sigma ^{\mathrm{FPC}}\) down to \(q_T/M \sim 10^{-6}\) GeV (this value is comparable with the typical technical cuts used in local subtraction methods).

We have encoded the formula in Eq. (10), by using the recoil prescription of Ref. [20],Footnote 1 in the public numerical program DYTurbo [23] which implements the \(q_T\) subtraction formalism for Drell–Yan processes. We stress that within our method we are able to correctly remove the linear fiducial power corrections which originates from a mismatch between the real and the counterterm (Born-like) kinematics in the standard implementation of the \(q_T\) subtraction method. In this way the inclusion of fiducial selection cuts does not introduce additional linear power corrections in \(q_T^{\mathrm{cut}}\) with respect to the ones present at inclusive level. Moreover, within our numerical implementation in DYTurbo, we are able to confirm numerically, up to the N\(^3\)LO, that our method correctly removes the linear power corrections and to quantify the residual systematic uncertainty from the \(q_T^{\mathrm{cut}}\) technical parameter in the Drell–Yan fiducial cross section.

We consider the production of \(l^+ l^-\) pairs from \(Z/\gamma ^*\) decay at the LHC (\(\sqrt{s}=13\) TeV) with the following fiducial cuts [31]: the leptons are required to have transverse momentum \(p_T>25\) GeV, pseudo-rapidity \(|\eta |<2.5\) while the lepton pair system is required to have invariant mass \(66<M<116\) GeV and transverse momentum \(q_T < 100\) GeV. We use parton densities functions (PDFs) from the NNPDF3.1 set [32] at NNLO with \(\alpha _S(m_Z^2)=0.118\), and we have evaluated \(\alpha _S(\mu _R^2)\) at \((n+1)\)-loop order at N\(^n\)LO accuracy. Factorization and renormalization scales have been set to \(\mu _F=\mu _R= \sqrt{M^2+q_T^2}\). We use the so called \(G_\mu \) scheme for EW couplings with input parameters \(G_F = 1.1663787\times 10^{-5}\) GeV\(^{-2}\), \(m_Z = 91.1876\) GeV, \(\Gamma _Z=2.4952\) GeV, \(m_W = 80.379\) GeV [31]. We then computed the fiducial cross section for the Drell–Yan process at the LHC with the original \(q_T\) subtraction formula Eq. (2) and with the improved formula Eq. (10).

Fig. 2
figure 2

Fiducial cross section for the production of \(l^+ l^-\) pairs from \(Z/\gamma ^*\) decay at the LHC (\(\sqrt{s}=13\) TeV). NLO results with the \(q_T\) subtraction method (blue squared points) and the \(q_T\) subtraction without FPC (red circled points) at various values of \(q_T^{\mathrm{cut}}\), and with a local subtraction method (black line). Error bars indicate the statistical uncertainties from Monte Carlo numerical integration

Full size image

The results presented in this paper have been obtained by including the hard-collinear functions (see Eq. (3)) from Refs. [24, 33,34,35,36] while the \(\mathrm{Z+jets}\) results have been obtained from Refs. [31, 37, 38]. The counter-term in Eq. (4) have been obtained by analytic expansion of the resummation formula of Ref. [21] using the coefficients calculated in Refs. [39, 40]

In Fig. 1 we show the NLO fiducial cross section calculated for different values of the \(q_T^{\mathrm{cut}}\) technical parameter with the original \(q_T\) subtraction method (blue squared points) and with the modified formula in Eq. (10) (labeled as recoil \(q_T\) subtraction, red circled points). As a reference, we also show the result obtained with a local subtraction formalism (black line) which represents the exact (free from significant systematic uncertainties) prediction. The local result is obtained independently with the dipole subtraction formalism [1, 2] as implemented in MCFM [41,42,43]. Error bars in Fig. 1 indicate the statistical uncertainties from Monte Carlo numerical integration which turns out to be completely negligible. From Fig. 1 we observe that the systematic uncertainty (defined as the deviation from the local subtraction result) of the original \(q_T\) subtraction results increase linearly with \(q_T^{\mathrm{cut}}\) and it is around \(0.3\%\) at \(q_T^{\mathrm{cut}}=0.5\) GeV, \(0.6\%\) at \(q_T^{\mathrm{cut}}=1\) GeV, \(1\%\) at \(q_T^{\mathrm{cut}}=2\) GeV and \(2\%\) at \(q_T^{\mathrm{cut}}=4\) GeV. In order to obtain a systematic uncertainty below \(0.1\%\) level a calculation with \(q_T^{\mathrm{cut}}\lesssim 0.1\) GeV is necessary. Conversely the results obtained with the \(q_T\) subtraction without FPC have a systematic uncertainty for \(q_T^{\mathrm{cut}}=1\) GeV which is smaller than the statistical uncertainty of the local-subtraction result, which is \(0.01\%\). In Fig. 1 an interpolation of the \(q_T^{\mathrm{cut}}\) dependence of the modified (original) \(q_T\)-subtraction obtained with quadratic (linear and quadratic) terms is represented by the red dashed (blue solid) line.

In Fig. 2 we show the fiducial cross section at NNLO with the original \(q_T\) subtraction method (blue squared points) and with the modified formula in Eq. (10) (red circled points) together with the result obtained with a local subtraction formalism (black line). The local result is obtained with the sector improved subtraction formalism [44, 45] as implemented in FEWZ [46, 47]. Error bars in Fig. 2 indicate the statistical uncertainties from Monte Carlo numerical integration. Statistical uncertainties are at the level of \(0.1\%\) for the local subtraction results and at the level of \(0.1\%\) or larger (smaller) for the \(q_T\) subtraction results with \(q_T^{\mathrm{cut}}\lesssim 0.1\) GeV (\(q_T^{\mathrm{cut}} > rsim 0.1\) GeV). The \(q_T^{\mathrm{cut}}\) systematic uncertainty of the \(q_T\) subtraction results is around \(0.3\%\) at \(q_T^{\mathrm{cut}}=0.5\) GeV, \(0.6\%\) at \(q_T^{\mathrm{cut}}=1\) GeV and \(0.7\%\) at \(q_T^{\mathrm{cut}}=2\) GeV and \(0.2\%\) at \(q_T^{\mathrm{cut}}=4\) GeV. As in the case of the NLO results, in order to obtain a systematic uncertainty below \(0.1\%\) level a calculation with \(q_T^{\mathrm{cut}}\lesssim 0.1\) GeV is necessary. However this is exactly the IR region where large cancellations give rise to sizable statistical uncertainties due to numerical integration. Conversely the results obtained with the \(q_T\) subtraction without FPC have a systematic uncertainty which is smaller than \(0.04\%\) for \(q_T^{\mathrm{cut}}=1\) GeV. As in Fig. 1 also in Fig. 2 we have shown an interpolation of the \(q_T^{\mathrm{cut}}\) dependence of the results.

Fig. 3
figure 3

Fiducial cross section for the production of \(l^+ l^-\) pairs from \(Z/\gamma ^*\) decay at the LHC (\(\sqrt{s}=13\) TeV). NLO results with the \(q_T\) subtraction method (blue squared points) and the \(q_T\) subtraction without FPC (red circled points) at various values of \(q_T^{\mathrm{cut}}\)

Full size image

Finally, in Fig. 3 we show the fiducial cross section at N\(^3\)LO with the original \(q_T\) subtraction method (blue squared points) and with the modified formula in Eq. (10) (red circled points) for different values of \(q_T^{\mathrm{cut}}\) with the interpolation of the results as in Figs. 1 and 2. No local subtraction results are available at this perturbative order. Moreover in this case we are not able to show results for \(q_T^{\mathrm{cut}}< 4\) GeV. In fact we have checked that our analytic expression for the counter-term \(d\widetilde{\sigma }^{\mathrm{CT}}\) agrees with the small-\(q_T\) limit of the NNLO fixed-order results for the production of a \(Z/\gamma ^*\) boson in association with jets reported in Ref. [31] at permille level down to \(q_T\sim 4\) GeV while below that threshold such agreement deteriorates.

We observe, in the case of the \(q_T\) subtraction without FPC, a strong reduction of the dependence from \(q_T^{\mathrm{cut}}\) also at the N\(^3\)LO and we have quantified the residual \(q_T^{\mathrm{cut}}\) dependence by performing an extrapolation fit of the subleading power corrections at NLO, NNLO and N\(^3\)LO. We performed an interpolation fit for \(q_T>q_T^{\mathrm{cut}}\) of the known structure of subleading power corrections [8, 14, 26] and we estimated the uncertainty on the extrapolated cross section for \(q_T^{\mathrm{cut}} \rightarrow 0\). The interpolation curves have been shown in Figs. 1, 2 and 3. The uncertainty on the extrapolated results is below \(0.01\%\) level at the NLO and NNLO (interpolation range \(0.05 \le q_T^{\mathrm{cut}} \le 20 \) GeV) and of \(0.18\%\) at the N\(^3\)LO (interpolation range \(4 \le q_T^{\mathrm{cut}} \le 20 \) GeV).Footnote 2 At the NLO and NNLO, this uncertainty is smaller than statistical errors of our results, which is consistent with the agreement between our predictions and the local results reported in Table 1. In the N\(^3\)LO case the interpolation convergence is less good and we obtain an uncertainty on the extrapolated cross section which is of the same order of the statistical error quoted in Table 2.

Fig. 4
figure 4

Production of \(l^+ l^-\) pairs from \(Z/\gamma ^*\) decay at the LHC (\(\sqrt{s}=13\) TeV). Power correction contributions at \({\mathcal {O}}(\alpha _S)\), \({\mathcal {O}}(\alpha _S^2)\) and \({\mathcal {O}}(\alpha _S^3)\) at various values of \(q_T^{\mathrm{cut}}\)

Full size image

In order to quantify the impact of the calculated fiducial power corrections, we show in Fig. 4 the contribution of the FPC (Eq. (10)) as a function of \(q_T^{\mathrm{cut}}\). First of all we observe that the sign of the FPC contribution changes from \({\mathcal {O}}(\alpha _S)\) to \({\mathcal {O}}(\alpha _S^2)\) and from \({\mathcal {O}}(\alpha _S^2)\) to \({\mathcal {O}}(\alpha _S^3)\). This behaviour is consistent with the observation that linear power corrections in the small \(q_T\) region (produced by the fiducial cuts) results in an alternating-sign factorial growth of the fixed-order perturbative series [18]. The second observation is that the impact of the FPC is not numerically reduced at higher orders and it turn out to be particularly sizable at N\(^3\)LO up to very small value of \(q_T^{\mathrm{cut}}\): for \(q_T^{\mathrm{cut}}= 0.05\) GeV the impact of the N\(^3\)LO FPC is about \(-0.4\%\) and it is the result of a \(+0.3\%\) contribution at \({\mathcal {O}}(\alpha _S^2)\) and a \(-0.7\%\) at \({\mathcal {O}}(\alpha _S^3)\) (the \({\mathcal {O}}(\alpha _S)\) FPC contribution turns out to be negligible at \(q_T^{\mathrm{cut}}= 0.05\) GeV). This means that when standard selection cuts are implemented within the original \(q_T\) subtraction, a permille level systematic accuracy for NNLO and N\(^3\)LO fiducial cross sections cannot be easily reached even with extremely low values of \(q_T^{\mathrm{cut}}\).

Table 1 Fiducial cross sections at the LHC (\(\sqrt{s}=13\) TeV): fixed-order results at NLO and NNLO. The uncertainties on the values of the cross sections refer to the statistical uncertainties from Monte Carlo numerical integration
Full size table

In Table 1 we report the predictions for the cross section in the fiducial region at NLO and NNLO with the \(q_T\) subtraction method for \(q_T^{\mathrm{cut}}= 0.5\) GeV and \(q_T^{\mathrm{cut}}= 1\) GeV, with the recoil \(q_T\) subtraction for \(q_T^{\mathrm{cut}}= 1\) GeV and we compare with the local subtraction results.Footnote 3 Errors indicate the statistical uncertainties from Monte Carlo numerical integration. In the case of the recoil \(q_T\) subtraction the results are nearly independent by \(q_T^{\mathrm{cut}}\) for \(q_T^{\mathrm{cut}}\lesssim {\mathcal {O}}\)(GeV).

From the results of Table 1 we observe that the differences between the recoil \(q_T\) subtraction results and the local subtraction results are of \({\mathcal {O}}(0.01\%)\). Thus the modified \(q_T\) subtraction formula in Eq. (10) allows us to obtain accurate permille level predictions for fiducial cross section with values of \(q_T^{\mathrm{cut}}\sim {\mathcal {O}}\)(GeV).

Table 2 Fiducial cross sections at the LHC (\(\sqrt{s}=13\) TeV): fixed-order results at N\(^3\)LO. The uncertainties on the values of the cross sections refer to the statistical uncertainties from Monte Carlo numerical integration
Full size table

In Table 2 we report the predictions for the cross section in the fiducial region at N\(^3\)LO with the \(q_T\) subtraction method and with the recoil \(q_T\) subtraction for \(q_T^{\mathrm{cut}}= 4\) GeV. Local subtraction results are not available at this order.

From the results shown in Fig. 4 we could expect that a value of \(q_T^{\mathrm{cut}}\sim 4\) GeV is associated with a systematic uncertainty due to the FPC of around \(2\%\), which is of the same order of the size of the \(\alpha _S^3\) corrections and thus challenge the reliability of the \(q_T\) subtraction results.Footnote 4 However the cancellation of the alternating sign linear fiducial power corrections shown Fig. 4 makes the impact of the N\(^3\)LO FPC for the particular value of \(q_T^{\mathrm{cut}}= 4\) GeV to be around \(0.2\%\) which is indeed the difference between the \(q_T\) subtraction and recoil \(q_T\) subtraction results reported in Table 2.

We finally note that any numerical implementation of the subtraction method, including the local versions, contains and depends on various technical parameters necessary to avoid the numerical evaluation of singular points. These parameters cannot be arbitrarily large and their actual value has to be chosen in order to make the numerical result independent (within the required numerical accuracy) from their actual value. From this viewpoint the independence by \(q_T^{\mathrm{cut}}\) observed within our method is similar to the one observed in the local version of the subtraction method.

In this paper we have considered higher-order QCD corrections to the production of high-mass systems in hadron collisions within the \(q_T\) subtraction formalism. We have presented a method to consistently remove the linear power corrections in \(q_T\) of the type \({\mathcal {O}}(q_T^{\mathrm{cut}}/M)\), where \(q_T^{\mathrm{cut}}\) is the technical parameter necessary to separate resolved and unresolved parton emission regions, which appears when fiducial kinematical cuts are applied on the final state system. As a first application we have implemented our method within the DYTurbo numerical program and we have considered explicitly the case of fiducial cross sections for \(Z/\gamma ^*\) boson production at the LHC up to N\(^3\)LO in QCD. We have obtained perturbative predictions which are in excellent (permille level) agreement with those obtained with local subtraction formalism at NLO and NNLO and we have computed N\(^3\)LO predictions with a residual \(q_T^{\mathrm{cut}}\) systematic uncertainty at the permille level.

Our results can be helpful in increasing the numerical precision of the existing numerical codes based on the \(q_T\) subtraction formalism and also on improving their time performances. In particular we were able to remove the source of systematic uncertainty at the origin of the discrepancies observed in Ref. [12]. Moreover our method is particularly important in the cases where fully local perturbative calculations for cross section are not available or when large numerical integration uncertainties are associated to the \(q_T\rightarrow 0\) limit such as in the case of N\(^3\)LO predictions in hadron collisions.

Finally, we make some observations about some consequences of our findings on resummed calculations. Our results show that resummed fiducial cross sections calculated by correctly taking into account the \(q_T\) recoil effects [20] (e.g. the resummed cross sections calculated in Ref. [9]) are free from significant numerical systematic uncertainties due to the minimum value of \(q_T\) used to compute the finite component of the cross section in Eq. (9). Moreover our results show that the matching between the resummed and finite (fixed-order) calculations in the small \(q_T\) region has an impact of \({\mathcal {O}}((q_T/M)^2)\) and it is expected to have a very small (negligible) effect for \(q_T/M\lesssim {\mathcal {O}}(10^{-1})\) (\(q_T/M\lesssim {\mathcal {O}}(10^{-2})\)).