Four-dimensional unsubtraction from the loop-tree duality

We present a new algorithm to construct a purely four dimensional representation of higher-order perturbative corrections to physical cross-sections at next-to-leading order (NLO). The algorithm is based on the loop-tree duality (LTD), and it is implemented by introducing a suitable mapping between the external and loop momenta of the virtual scattering amplitudes, and the external momenta of the real emission corrections. In this way, the sum over degenerate infrared states is performed at integrand level and the cancellation of infrared divergences occurs locally without introducing subtraction counter-terms to deal with soft and final-state collinear singularities. The dual representation of ultraviolet counter-terms is also discussed in detail, in particular for self-energy contributions. The method is first illustrated with the scalar three-point function, before proceeding with the calculation of the physical cross-section for $\gamma^* \to q \bar{q}(g)$, and its generalisation to multi-leg processes. The extension to next-to-next-to-leading order (NNLO) is briefly commented.


Introduction
Computation of higher-order perturbative corrections in quantum field theories is a very active research field in these days in view of the high-quality data provided by the LHC. The standard approach to perform these calculations in perturbative QCD relies in the application of the subtraction formalism. According to Kinoshita-Lee-Nauenberg theorem [1,2] theoretical predictions in theories with massless particles can only be obtained after the definition of infrared-safe physical observables. These involve performing a sum over all degenerate states, which means adding together real and virtual contributions. After ultraviolet (UV) renormalisation of virtual scattering amplitudes, the remaining contributions develop infrared (IR) singularities that cancel when putting all the terms together. This implies that the IR divergent structure of real and virtual corrections are closely related: thus the subtraction method exploits the factorisation properties of QCD to define suitable subtraction IR counter-terms which mimic their singular behaviour.
There are several variants of the subtraction method at NLO and beyond [3,4,5,6,7,8,9,10,11,12,13], which involve treating separately real and virtual contributions. However, from a computational point of view, these methods might not be efficient enough for multi-particle processes. The main reason for that is related to the fact that the final-state phase-space (PS) of the different contributions involves different numbers of particles. For instance, at NLO, virtual corrections with Born kinematics have to be combined with real contributions involving an additional final-state particle. The IR counter-terms have to be local in the real PS, and analytically integrable over the extra-radiation factorised PS to properly cancel the divergent structure present in the virtual corrections. Building these counter-terms represents a challenge and introduces a potential bottleneck to efficiently carry out the IR subtraction for multi-leg multi-loop processes.
With the aim of eluding the introduction of IR counter-terms, we explore an alternative idea based in the application of the loop-tree duality (LTD) [14,15,16,17,18,19,20,21,22,23,24,25,26,27]. The LTD theorem establishes that loop scattering amplitudes can be expressed as a sum of PS integrals (i.e. the so-called dual integrals) with an additional physical particle. Dual integrals and realradiation contributions exhibit a similar structure, and can be combined at integrand level. As shown in Refs. [19,25,26,27], the divergent behaviour of both contributions is matched, and the combined expression is finite. In other words, working in the context of dimensional regularisation (DREG) [28,29,30,31] with d = 4 − 2 the number of space-time dimensions, the mapped real-virtual contributions do not lead to -poles, which implies that the limit → 0 can safely be considered. This fact has a strong implication: the possibility of carrying out purely four-dimensional implementations for any observable at NLO and higher-orders. Thus, the aim of this paper is to implement a novel algorithm for a four-dimensional regularisation of multi-leg physical cross-sections at NLO free of soft and final-state collinear subtractions.
It is worth to mention that the idea of obtaining purely four-dimensional expressions to compute higher-order observables has been previously studied. For instance, it was proposed to apply momentum smearing [32,33,34,35] to combine real and virtual contributions, thus achieving a local cancellation of singularities. Other alternative methods consist in rewriting the standard IR/UV subtraction counterterms in a local form, as discussed in Refs. [36,37], or modifying the structure of the propagators (and the associated Feynman rules) [38,39,40] to regularise the singularities. Besides that, the numerical computation of virtual corrections has received a lot of attention in recent years [41,42,43,44,45,46,47,48,49,50,51]. For these reasons, through the application of LTD, we will tackle both problems simultaneously; we will express virtual amplitudes as phase-space integrals and combine them with the real contributions, working directly at integrand level. Moreover, physically interpretable results will emerge in a natural way.
The outline of this article is the following. In Section 2 we introduce and review the basis of LTD. As starting point, we describe in detail the unsubtracted implementation of NLO corrections with an scalar toy example. First, the IR singular structure of the scalar three-point function and the construction of the corresponding dual integrals is commented in Section 3. Second, the mapping of momenta between real and virtual corrections is defined in Section 4 for the toy example, and our four-dimensional regularisation of soft and collinear singularities is presented with a strong physical motivation. Then, in Section 5 we study the renormalisation of UV divergences at integrand level in the LTD framework. The discussion is focused on the treatment of scalar two-point functions, and we properly rewrite unintegrated dual UV counter-terms in a fully local way. After that, we carefully analyse the implementation of these techniques to the process γ * → qq(g) in Section 6. We put special emphasis in the algorithmic construction of the integrands, and in the numerical implementation of the purely four-dimensional representation. In Section 7, we generalise the unsubtraction algorithm to multi-leg processes, and briefly comment about the extension of the algorithm to NNLO. Finally, we present the conclusions and discuss the future research directions in Section 8.

Review of the loop-tree duality
In this section we review the main ideas behind the LTD method. The LTD theorem [14] establishes a direct connection among loop and phase-space integrals. Explicitly, it demonstrates that loop contributions to scattering amplitudes in any relativistic, local and unitary quantum field theory can be computed through dual integrals, which are build from single cuts of the virtual diagrams. Let's consider a generic N -particle scalar one-loop integral, i.e.
over Feynman propagators G F (q i ) = (q 2 i − m 2 i + ı0) −1 , whose most general topology is shown in Fig. 1. Then, there is a corresponding dual representation consisting of the sum of N dual integrals: where are dual propagators, and i, j ∈ α 1 = {1, 2, . . . N } label the available internal lines. In Eq. (1) and Eq. (2), the masses and momenta of the internal lines are denoted m i and q i,µ = (q i,0 , q i ), respectively, where q i,0 is the energy and q i are the spatial components. In terms of the loop momentum and the outgoing four-momenta of the external particles p i , the internal momenta are defined as together to the constraint k N = 0 imposed by momentum conservation. On the other hand, the d-dimensional loop measure is given by andδ is used in Eq. (2) to set the internal lines on-shell. Moreover, the presence of the Heaviside function restricts the integration domain to the positive energy region (i.e. q i,0 > 0). Since LTD is derived through the application of the Cauchy's residue theorem, the remaining d − 1 dimensional integration is performed over the forward on-shell hyperboloids defined by the solution of G F (q i ) −1 = 0 with q i,0 > 0. Notice that these on-shell hyperboloids degenerate to light-cones when internal particles are massless.
The dual representation shown in Eq. (2) is built by adding all possible single-cuts of the original loop diagram. In this procedure, the propagator associated with the cut line is replaced by Eq. (6) whilst the remaining uncut Feynman propagators are promoted to dual ones. The introduction of dual propagators modifies the ı0-prescription since it depends on the sign of η · k ji , with η a future-like vector, η 2 ≥ 0, with positive definite energy η 0 > 0, and k ji = q j − q i , which is independent of the loop momentum at one-loop. According to the derivation shown in Ref. [14], η is arbitrary so we can chose η µ = (1, 0) to simplify the implementation.
The difference between LTD and the Feynman Tree Theorem (FTT) [52,53], where the loop integral is obtained after summing over all possible m-cuts, is codified in the dual prescription: correlations coming from multiple cuts in FTT are recovered in LTD by considering only single-cuts with the modified ı0-prescription. In other words, having different prescriptions for each cut is a necessary condition for the consistency of the method. As discussed in Ref. [19], the integrand in Eq. (2) becomes singular at the intersection of forward on-shell hyperboloids (FF case), and forward with backward (q j,0 < 0) on-shell hyperboloids (FB intersections). On one hand, the FF singularities cancel each other among different dual contributions; the change of sign in the modified prescription is crucial to enable this behaviour. On the other hand, the singularities associated with FB intersections remain constrained to a compact region of the loop-three momentum space and are easily reinterpreted in terms of causality. From a physical point of view, FB singularities take place when the on-shell virtual particle interacts with another on-shell virtual particle after the emission of outgoing on-shell radiation. The direction of the internal momentum flow establishes a natural causal ordering, and this interpretation is consistent with the Cutkosky rule. In fact, the total energy of the emitted particles, which is equal to q i,0 − |q j,0 |, has to be positive. Together with the positive energy constraint imposed by the delta distribution in Eq. (6), it restricts the possible situations compatible with a sequential decay of on-shell physical particles.

Singularities of the scalar three-point function
In this section, we show a detailed derivation of the results presented in Ref. [25] concerning the scalar three-point function with massless internal particles. This discussion is useful to analyse and understand the application of LTD to the realistic case presented in Section 6, and the posterior generalisation to multi-leg processes in Section 7. We consider final-state massless and on-shell momenta labeled as p 1 , p 2 , and the incoming momenta p 3 = p 1 + p 2 ≡ p 12 , by momentum conservation, with virtuality p 2 3 = s 12 > 0. The internal momenta are q 1 = + p 1 , q 2 = + p 12 and q 3 = , where is the loop momentum. The scalar three-point function at one-loop is given by the well-known result [54,55] where c Γ is the usual loop volume factor (see Eq. (126) in Appendix A). The LTD representation of the scalar integral in Eq. (7) consists of three contributions with (−2q 2 · p 2 + ı0) (−2q 2 · p 12 + s 12 + ı0) , In order to simplify the computation of these integrals, we work in the centre-of-mass frame of p 1 and p 2 , and parametrise the momenta as with ξ i,0 ∈ [0, ∞) and v i ∈ [0, 1] the integration variables describing the energy and polar angle of the loop momenta, respectively. Integration of the loop momentum in the transverse plane, which is described by the unit vectors e i,⊥ , is trivial in this case. The scalar products of internal with external momenta are given by and the dual integrals in Eq. (9) are rewritten as with the integration measure in DREG given by the direct product of It is possible to perform these integrals analytically; the result is where the phase-space volume factor c Γ is defined in Eq. (126). As expected, the sum of the three dual integrals agrees with the well-known result from Eq. (7). It is worth noticing here that the integrand of I 1 is both IR and UV divergent. However, the application of DREG leads to equal and opposite -poles, which justifies the result shown in Eq. (14). We will discuss this fact more carefully in Sec. 5.1, since it plays a crucial role in the whole implementation of the LTD approach for physical processes.
Notice that in Eq. (12) the dual +ı0 prescription is crucial for computing I 2 , because 1 − ξ 2,0 changes sign inside the integration region, leading to a threshold singularity. For later use, it will be necessary an explicit expression for the imaginary part of I 2 at the integrand level, which is determined by setting the integrand of I 2 on the negative-energy on-shell mode of G F (q 3 ), i.e.
We remark that G D (q 2 ; q 1 ) = G F (q 1 ) because η · k 12 = −p 2,0 < 0 with η µ = (1, 0), and hence Eq. (15) is consistent with the Cutkosky's rule, which is determined by setting two propagators on-shell and outgoing, namely, by reversing the momentum flow of one of the on-shell internal lines [14]. This is the causality connection mentioned in Section 2, and it becomes relevant in our computation because the -expansion of Eq. (15) reveals the presence of a purely imaginary single-pole in I 2 that will not be cancelled by real corrections. At integrand level, this means that the real part of I 2 presents an integrable singularity in the neighbourhood of ξ 2,0 = 1, but there is also a non-integrable singularity that must be cancelled by removing its imaginary component before performing a four-dimensional numerical implementation. Thus, the real part of I 2 is defined as and, by virtue of the Sokhotski-Plemelj theorem, where we make use of the Cauchy's principal value (PV) to get rid of the +ı0 prescription and the imaginary pole. From the formal point of view, we could perform this computation by simply working with the real part of the integrand (and neglecting the prescription). However, numerical instabilities arise and the application of PV prescription leads to a more efficient implementation.
x + ξ 2 y ; two dimensions (left) and three-dimensions (middle). LTD is equivalent to integrate along the forward-light cones (solid lines in the left plot). Backward light-cones are represented by dashed lines. The intersection of light-cones leads to soft, collinear, and threshold singularities in the loop three-momentum space (right plot), or to the cancellation of singularities among dual contributions.
LTD can further be exploited to have a deeper and detailed understanding of the origin of the singularities of the loop integral under consideration. As commented before, the origin of the singularities can be underlined by analysing the relative position and intersections of the on-shell hyperboloids or lightcones of the propagators entering the loop integrand [19]. In Fig. 2, we plot the light-cones that support each of the dual integrals of the dual representation of the three-point function in Eq. (8). Although the scalar three-point function is UV finite, the individual dual integrals in Eq. (9) diverge in the UV because dual propagators are linear in the loop momentum. However, their sum has the same UV singularities present in the original integral, as expected according to the LTD theorem. Thus, we can focus on its IR behaviour; renormalisation of UV divergences will be considered later in Section 5.
Collinear divergences are associated with an extended region supported in the intersection of lightcones, as shown in Fig. 2. At large loop momentum, the intersections occur among forward light-cones and the collinear singularities cancel in the sum of dual integrals. However, there are still collinear divergences originated in the compact region defined by the intersection of forward and backward lightcones. Soft singularities arise at q 1,0 = 0. Finally, a threshold singularity appears at the dual integral I 2 through the intersection of the backward light-cone of G F (q 3 ) with the forward light-cone of G F (q 2 ). The imaginary part of I 2 is singular but this singularity can be removed by using Eq. (15). The singularity of the real part of I 2 is integrable but will lead to instabilities by a direct numerical computation. In that case, a contour deformation is employed to achieve an stable numerical implementation [21,22,23,24]. Fig. 2 and in order to isolate the IR divergences of the scalar three-point function, we define the soft and collinear components of the dual integrals in well-defined compact regions of the loop three-momentum, i.e.

Motivated by
where 0 < w < 1 is a cut in the energy of the internal on-shell particle. The collinear singularity of the dual integral I 2 appears at v 2 = 1 with ξ 2,0 ∈ [0, 1], but I 2 develops also a threshold singularity at ξ 2,0 = 1. For that reason we have imposed a finite w-cut to include the threshold region in the definition of I 1 are due to collinear singularities only. There is some arbitrariness in the definition of the integration regions of these integrals. Provided that we include the soft and collinear singularities, different definitions will differ only in the finite contributions. Indeed, we will redefine them later with a better motivated physical target. These are, however, the simplest choice for the current illustrative purpose. Analytic integration gives Using Pfaff and shift identities, the hypergeometric function in I (c) 2 can be written in the physical region with w > 0 and s 12 > 0. This leads to As expected, the soft integral in Eq.
Outside the region that contains the IR poles, the sum of the dual integrals is finite, although they are separately divergent. A suitable combination is required to obtain finite results. So, we define the forward and the backward regions as those contained in v i ≤ 1/2 and v i ≥ 1/2, respectively. This separation does not have any physical meaning; it is just convenient for the analytic computation. Explicitly, we define where we have performed a trivial change of the integration variables, i.e.
and we have taken the limit → 0 at integrand level. Notice that each dual integrand is still individually singular. For instance, I 1 and I 3 are divergent at v 1 = 0 = v 3 but their sum is finite in this limit, although UV divergences survive in the sum. The divergent high-energy behaviour of I 1 + I 3 is cancelled once we add the contribution of I 2 . These cross-cancellations of singularities allow to perform the integral of the forward contribution with = 0. The logarithmic terms in Eq. (23) are originated from the fact that we are using different coordinate systems for each dual integral. This produces a mismatch of the integration measures that is of O( ). Since the integral behaves as O( −1 ) in the collinear limit, a non-vanishing finite contribution arises from the collinear region. Explicitly, the expansion of v −1− by using Eq. (136) leads to where we removed the C-distribution. These logarithms are avoidable by a proper reparametrisation of the integration variables; we will return to this point in Section 3.1. From Eq. (23), we obtain In an analogous way, we can compute the finite contribution originated in the backward region (v i ≥ 1/2), which is given by where we throw the +ı0 prescription because w > 0 excludes the threshold singularity from the integration region. Also, to obtain Eq. (27) starting from Eq. (12), we used the following change of variables and then we took the limit → 0 at integrand level. Similarly to the forward integral, there is a cancellation of collinear singularities among I 1 and I 2 , which takes place at v 1 = 1 = v 2 in this case; thus I 1 + I 2 is IR-finite but it is still UV-divergent. To regularise the high-energy behaviour, we need to add the contribution of I 3 . Again, notice that we have introduced some logarithmic terms in Eq. (27). These terms are due to the mismatch in the collinear behaviour of I 1 and I 2 at O( ). As we did for the forward case, we derived the logarithmic corrections by expanding the collinear factor (1 − v) −1− and the integration measure after the implementation of Eq. (28); i.e.
and, again, the C-distribution can be removed because the integrand is regular for v = 1. So, the integral over the sum of the three dual integrands can be performed with = 0, and we obtain The sum of Eq.
which is independent of w. This was expected because w is a non-physical cut. It is worth noting that only I IR contains the -poles and the remaining contributions have been computed directly with = 0. Moreover, through the application of Eqs. (25) and (29), the integrand can easily be expressed as the → 0 limit of the original DREG expression plus some logarithmic corrections which lead to the right result.

Unifying the coordinate system
As a final remark, we point out that it is possible to avoid the introduction of the extra logarithmic terms in I (f) and I (b) . They are originated from the fact that each dual integrand has been expressed in a different coordinate system and they approach the collinear limit in a slightly different way at O( ). The solution consists in using the same coordinate system for all the dual integrals, i.e. we have to map exactly the loop three-momenta q i . Although for analytic computations this leads to more complex integration limits, for numerical applications it is indeed the natural choice for the implementation of this approach. For instance, for the forward integral I (f) the integration variables (ξ 1,0 , v 1 ) must be written in terms of (ξ 3,0 , v 3 ). Notice that q 3 = is set on-shell in I 3 , but not in I 1 where q 2 1 = 0. From Eq. (10), the spatial components of the internal loop momenta are parametrised according to with q (+) 1,0 = (q 3 + p 1 ) 2 − ı0 when q 1 is on-shell, which leads to With this change of variables, we find the following representation of the forward integral with = 0, which is free of the logarithmic contributions that appear in Eq. (23), and leads to the same result as in Eq. (26). Notice that in Eq. (34) we used v 2 = v 3 = v and ξ 2,0 = ξ 3,0 = ξ 0 . A similar representation is available for the backward integral, where we must combine I 1 and I 2 , by expressing (ξ 1,0 , v 1 ) in terms of (ξ 2,0 , v 2 ). Explicitly, we use the same change of variables as in Eq. (33), which leads to where we also applied v 2 = v 3 = v and ξ 2,0 = ξ 3,0 = ξ 0 . The integration limits in Eq. (35), codified through the Heaviside theta functions, are more cumbersome than in Eq. (27), but the result of both expressions is the same, and is given by Eq. (30).
In summary, the IR singularities of loop integrals are restricted to a compact area of the integration domain, and the finite remnants are expressible in terms of pure four-dimensional functions, which implies that DREG could be avoided. Still, we need to keep d = 4 to deal with I IR ; in the following we will show how to overcome this issue to achieve a fully four-dimensional implementation.

Unsubtraction of soft and collinear divergences
In Section 3, we have illustrated in detail the application of LTD to a scalar one-loop Feynman integral and we have isolated its infrared divergences in the function I IR (see Eq. (22)), which is obtained from a compact region of the loop three-momentum. In the framework of LTD, a suitable mapping of external and loop momenta between virtual and real corrections allows to cancel the IR singularities at the integrand level, such that a full four-dimensional implementation is achieved without the necessity to introduce soft and collinear subtraction terms [25]. We illustrate the method with a simplified toy scalar example before affording a complete calculation in a realistic physical process in Section 6.
We consider the one-loop virtual corrections to the cross-section that are proportional to the scalar three-point function where is the Born cross-section, dΦ 1→2 is the integrated phase-space volume, given in Eq. (133), and g is a generic coupling.
On the other hand, it is necessary to include also the real radiation due to 1 → 3 processes. The momenta configuration is p 3 → p 1 + p 2 + p r , where we keep the same incoming momentum as in the 1 → 2 contributions, where p 3 → p 1 + p 2 with p 2 3 = s 12 , since we aim to obtain a local cancellation of singularities. The real radiation correction to the cross-section is given by with s ir = (p i + p r ) 2 . The real corrections included in Eq. (38) can be understood as the interference of the two scattering amplitudes corresponding to the emission of the real radiation from each of the outgoing particles. We do not take into account for the moment the squares of these amplitudes, which are proportional to 1/s 2 ir . They are topologically related to self-energy diagrams and will be considered explicitly in Section 6 for the physical process γ * → qq(g). For the current illustrative purpose it is enough to consider this interference.
Then, we split the three-body phase-space to isolate the different IR singular regions. This strategy is a common practice in the context of subtraction methods [3,4], because it allows to optimise the local cancellation of the collinear singularities at integrand level. Since there are three particles in the final state and the incoming one is off-shell, it is enough to separate the three-body phase-space into two pieces by making use of the identity which leads to the definitions where y ir = s ir /s 12 are dimensionless variables. Analogously, we define the corresponding dual contributions to the virtual cross-section as with the dual components of the one-loop scattering amplitude according to the decomposition suggested in Eq. (8). So, we claim that with i ∈ {1, 2}, is finite in the limit → 0 and can be expressed using a purely four-dimensional representation. It is worth appreciating that the dual integral I 3 is not necessary to cancel the IR singularities present in the real corrections. This behaviour was expected from the analysis shown in Section 3, explicitly from Eqs. (19) and (22). The dual integral I 3 does not lead to collinear divergences that are not cancelled by the other dual contributions, therefore it is not necessary to define the corresponding σ 3 . In fact, I 3 will solely contribute to the definition of the IR finite virtual remnant, formerly described in terms of the backward and forward integrals.
We now implement a mapping between the final-state momenta of the loop amplitudes {p 1 , p 2 }, the loop three-momentum , and the final-state momenta of the real amplitudes {p 1 , p 2 , p r }. Momentum conservation and on-shell constraints must be fulfilled by p i and p i , simultaneously. Hence, assuming q 1 on-shell and q 3 = q 1 − p 1 off-shell, we propose to perform the evaluation of the dual cross-section in Eq. (43). This mapping has many interesting properties, which deserve to be discussed. In first place, momentum conservation is automatically fulfilled as p 1 + p 2 + p r = p 1 + p 2 , and all the final-state momenta in Eq. (44) are on-shell. Also, it is suitable to describe collinear configurations where p 1 q 1 , which are reached for α 1 → 0; this is crucial to properly combine the divergent regions of the loop and real contributions to the cross-section and achieve a fully local regularisation. Notice that although we are restricted to 1 → 2 and 1 → 3 kinematics, the mapping in Eq. (44) can easily be extended to processes with an arbitrary number of external particles (see Section 7).
The next step consists in using the parametrisation of q i and p i proposed in Eq. (10), together with the mapping in Eq. (44), to rewrite the two-body kinematic invariants y ij in terms of the integration variables (ξ 1,0 , v 1 ). Expressing the scalar products p i · p j with both sets of variables, we find Since this mapping is optimised for the description of the collinear limit p 1 q 1 , it must be used in the regions of the two-body and three-body phase-space where y 1r < y 2r * . The lower limit in the value of y 2r avoids to deal with p 2 q 1 . Thus, a second mapping is necessary to treat the collinear limit with y 2r → 0, that is isolated in the region corresponding to y 2r < y 1r . With q 2 on-shell and q 1 = q 2 − p 2 off-shell, we define and the two-body invariants are given by (47) * To be mathematically rigourous, the transformation proposed in Eq. (45) is a diffeomorphism connecting the physical three-body phase-space and its image in the integration domain of the dual contributions. Thus, in principle, it would not be necessary to define a new mapping, since it covers the whole phase-space in this simple example, although it is not really optimised for the collinear limit p 2 q 1 . Figure 3: The dual integration regions in the loop three-momentum space.
By virtue of Eq. (39), the complete three-body phase space for the real radiation can be parametrised by applying Eq. (44) and Eq. (46); in fact, it is useful to define to parametrise the integration regions for σ 1 and σ 2 , respectively. A graphical representation of the integration regions defined by Eqs. (48) and (49) is shown in Fig. 3.
Once the momenta are properly parametrised, we proceed to evaluate together the real and virtual contributions at integrand level. For y 1r < y 2r , we obtain represents the analogous expression for y 2r < y 1r . The integrand in Eq. (50) has the form and the function f (ξ 1,0 , v 1 ) vanishes in the soft and collinear regions (i.e. ξ 1,0 = 0 and/or v 1 = 0); On the other hand, the integrand of σ The delta function in Eq. (51) cancels the imaginary part of the I 2 dual integral, which is given by Eq. (15). But indeed the condition y 2r < y 1r excludes the threshold singularity of I 2 from the integration region with the exception of the single point at (ξ 2,0 , v 2 ) = (1, 1). This fact allows to calculate Eq. (51) by removing the +ı0 prescription and the delta function. We obtain It is important to notice that this result can be reached following two different paths. The first one consists in using DREG and integrating with d = 4 − 2 ; once an analytic expression is obtained, we verify that no -poles are present and we take the limit → 0. On the other hand, we can consider directly the limit → 0 at integrand level; the expressions obtained are integrable and the results agree in both cases.
After combining virtual and real corrections, we define the virtual remnant σ V as the sum of the three dual integrals excluding the regions of the loop three-momentum already included in Eq. (50) and Eq. (51): This expression is analogous to the sum of the forward and backward contributions defined in Section 3, but does not require any unphysical cut w to deal properly with the threshold singularity. In the previous expression, we have identified all the integration variables, ξ 0 = ξ 2,0 = ξ 3,0 and v = v 2 = v 3 , while (ξ 1,0 , v 1 ) are expressed in terms of (ξ 3,0 , v 3 ) by using the change of variables in Eq. (33) to directly avoid the appearance of logarithmic contributions from the expansion of the integration measure. The integration regions are defined through The explicit expression for R 1 (ξ 0 , v) is, however, rather cumbersome, although this should not be a problem in numerical computations. For the analytic integration, we use a clever expansion of R 1 (ξ 0 , v) that exploits both reference systems; The first two terms in the right-hand side of Eq. (57) contribute at large loop three-momenta rendering the integral defined by Eq. (55) finite in the UV. The next two terms provide a finite contribution. The virtual remnant in Eq. (55) is also IR finite, and therefore it can be calculated with = 0. In particular, we get The sum of all the contributions, Eq. (53), Eq. (54) and Eq. (58), gives a total cross-section of O( ), in agreement with the result that would be obtained from the standard calculations in DREG.

Ultraviolet renormalisation
In the previous section we have shown how to avoid the introduction of subtraction counter-terms to cancel soft and collinear singularities by a suitable mapping of the momenta entering virtual and real corrections. In any practical computation in QFT, UV divergences must also be taken into account.
Another advantage of LTD is to enlighten the physical aspects of renormalisation. In order to explain the proposed approach, we consider first the simplest scalar two-point function, with massless internal lines. There is only one external momenta, p µ , and the two-point function is free of IR singularities if the incoming momentum is not light-like. Due to the fact that the virtuality of the incoming particle is the unique physical scale involved in the problem, the integral vanishes if we set p 2 = 0. So, the non-trivial massless scalar two-point function requires p 2 = 0 and it is only UV divergent. Labelling the internal momenta as q 1 = + p and q 2 = , then we have as shown in the literature † . The LTD representation of the scalar two-point function reads with the dual integrals where for simplicity we consider p 0 > 0 and p 2 > 0 (i.e. the incoming particle has positive energy and we work in the TL region). Following the discussion presented in Section 3, we parametrise the momenta using which is equivalent to settle in the rest frame of the incoming particle. With this choice, the dual integrals are rewritten as where we follow the definition of the d-dimensional integration measure given in Eq. (13); the only difference is that we must set (2p 0 ) 2 = 4p 2 instead of s 12 as the normalisation scale of the system: The integration in ξ i,0 and v i can be performed analytically, resulting in The sum of both contributions gives the standard DREG result in Eq. (59). The imaginary part of the scalar two-point function is calculable as which agrees with the result directly obtained from Eq. (59). This imaginary component is associated with I 1 and it is due to the presence of a threshold singularity. This threshold behaviour can graphically be explained from Fig. 4 (left); the forward light-cone associated with I 1 intersects only once the backward region of the other light-cone, and it is the only dual integral which requires to keep explicitly the +ı0 prescription in Eq. (63) .
In order to build a suitable local UV counter-term of the two point-function, we follow the ideas presented in Ref. [36] and we consider where µ UV is an arbitrary energy scale, and q UV = + k UV with k UV an arbitrary four-momentum. Notice that the counter-term is expressed as a Feynman integral using the customary +ı0 prescription.
We should now construct the corresponding dual representation and combine it with Eqs. (66) and (67). The counter-term in Eq. (69) exhibits a double pole in the complex plane of the loop energy component which requires to apply the extended version of the LTD theorem [17]. There are two possible paths to follow: either to compute the residue at the relevant higher-order pole and use the Cauchy's residue theorem to perform the 0 integral, or to apply integration-by-parts (IBP) identities to rewrite Eq. (69) in terms of a massive tadpole. Figure 4: Light-cones of the two-point function in the loop coordinates µ = p 0 (ξ 0 , ξ ⊥ , ξ z ) for a timelike, p 2 > 0 (left), and a light-like, p 2 = 0 (right), configuration. The internal momenta are q 1 = +p and q 2 = . In the time-like case, a threshold singularity appears when the backward light-cone (dashed) of G F (q 2 ) intersects with the forward light-cone (solid) of G F (q 1 ). In the light-like case, the IR singularities are restricted to the compact region ξ 1,0 ≤ 1.
Let's start with the first approach. The residue of a function f (z) at a multiple pole is given by with n the multiplicity of the pole at z 0 . The location of the double pole of the UV propagator in Eq. (69) is obtained from the on-shell condition, i.e.
where we just keep the solution q (+) UV,0 because it lies in the lower part of the complex plane and describes a positive-energy particle. The calculation of the residue gives and is a possible dual representation of the UV counter-term.
On the other hand, we can apply IBP to lower the power of the propagator, and then use LTD in its usual form. It is straightforward to obtain Both Eq. (73) and Eq. (74) are equivalent dual representations of the UV counter-term in Eq. (69). If we perform the change of variables ξ UV = q 2 UV /p 0 , m UV = µ UV /p 0 , the UV counter-term gets either of the two following dual forms: which corresponds to the explicit expressions for Eq. (73) and Eq. (74), respectively. We have dropped the +ı0 prescription because it is not necessary. After explicit computation, both representations in Eq. (75) lead to the same result: which successfully reproduce the single -pole present in the scalar two-point function. The prefactor ‡ is defined as S = (4π) Γ(1 + ).
The next step consists in combining a dual representation of the UV counter-term with the dual integrals I 1 and I 2 . So, we define the renormalised scalar two-point function as with ξ = ξ i,0 = ξ UV and v = v i = v UV , that verifies which is free of -poles. The integrand of Eq. (77) is integrable in the limit → 0, and the computation can be fully performed with = 0. So, we succeeded in finding a purely four-dimensional representation of the renormalised two-point function. Consistently, the integration and the limit → 0 commute.

Scaleless two-point function
As mentioned before, if we consider p 2 = 0 then L (1) (p, −p) = 0 since it does not contain any scale. From a physical point of view, this implies that self-energy corrections to on-shell massless particles are zero. However, a vanishing integral in DREG does not imply a vanishing integrand. This is particularly an issue in LTD where the aim is to cancel singularities locally, and therefore it is relevant to properly characterise the IR and the UV behaviour separately. So, let's consider the LTD representation of the ‡ We distinguish S from the usual MS scheme factor S MS = (4π) exp(− γ E ) or S = (4π) /Γ(1 − ) as used in Ref. [10]. At NLO all these definitions lead to the same expressions. At NNLO, they lead to slightly different bookkeeping of the IR and UV poles at intermediate steps, but physical cross-sections of infrared safe observables are the same. massless two-point function in the light-like case, i.e. p 2 = 0. We start from Eq. (61), and use the parametrisation for the incoming particle, whilst q µ i is given by Eq. (62). Then, the dual contributions become whose sum is equal to zero because they have opposite sign. In the context of DREG, the scaleless scalar two-point function develops both IR and UV divergences that cancel each other because the parameter IR regularising the IR singularities is identified with UV , which regulates the UV divergences. The important fact is that we can exploit LTD to separate them, proceeding analogously as we did in Section 3.
Analysing the integration domain of the dual contributions, and the relative position of the light-cone as shown in Fig. 4 (right), we realise that the IR singularity is associated with I 1 , because its forward light-cone overlaps with the backward light-cone of G F (q 2 ) in the region ξ 1,0 ≤ 1. Thus, we define which contains a single -pole. Again, as we found for the massless three-point function, the IR singularities are confined in a compact region of the loop three-momentum space. Outside this region, the remnant is given by which is renormalised with the UV counter-term defined in Eq. (73). It is worth noting that LTD leads naturally to this separation of IR/UV regions, which is crucial to achieve a local cancellation of singularities in the computation of physical observables at higher-orders.

Renormalisation of scattering amplitudes and physical interpretation
In general, the UV counter-terms of scattering amplitudes are derived by expanding the internal propagators around the UV propagator [36]. For a single propagator: with k i,UV = q i − q UV , and similarly with numerators. In order to improve the convergence in numerical implementations, the authors of Ref. [37] propose to expand even to higher-powers of G F (q UV ) −1 .
Once the desired UV expansion is obtained, we shall derive the corresponding dual representation following the procedure described before. To obtain a suitable LTD representation of the UV counterterms, it is necessary to deal with multiple-poles and non-trivial numerators depending on q UV . The following identities are enough to construct most of the UV counter-terms: These expressions have been obtained from Eq. (70). Alternatively, we can apply IBP to obtain equivalent dual representations. Explicit examples will be presented in Section 6. Only the genuine UV singularities of the original scattering amplitudes need to be subtracted with this procedure. The spurious UV singularities of the individual dual integrals are cancelled in the sum of all the dual contributions.
To conclude this section, let us remind about the physical interpretation of the arbitrary energy scale µ UV introduced in the renormalisation procedure [25]. This arbitrary scale can be interpreted as a renormalisation scale, since the UV counter-term affects the behaviour of the integrand only in the highenergy region. In fact, as seen in Fig. 5, its dual representation only contributes for loop energies larger than k UV,0 + µ UV , although it is unconstrained in the loop three-momentum. So, we should choose µ UV ≥ O(Q), with Q the physical hard scale that determines the size of the compact region in the loop momentum space where the IR and threshold singularities are located. However, the intersection of the integration domains of the dual integrands gives rise to singularities, thus we could improve the UV cancellations if we avoid those intersections. In other terms, we must choose µ UV in such a way that the on-shell hyperboloids of the UV propagator do not intersect with any of the on-shell hyperboloids of the original integral. Since the UV forward and backward on-shell hyperboloids are separated by a distance 2µ UV , it is possible to give a physical motivation for an optimal minimal choice of µ UV and k UV . If we take µ UV = Q/2, and k UV in the centre of the physical compact region, these conditions are fulfilled in a minimal way, i.e. we naturally avoid intersections with the physical on-shell hyperboloids. Figure 6: Momentum configuration and Feynman diagrams associated with the process γ * → qq(g), up to O(α α S ). Notice that we also consider self-energy corrections to the on-shell outgoing particles, even if their total contribution is zero.

NLO corrections to γ * → qq(g)
In this section, we discuss in detail the computation of NLO QCD corrections to the total cross-section for the process γ * → qq(g) by using the LTD approach. We emphasise that it constitutes the first realistic physical application of this method, as already anticipated in Ref. [27].
The computation is done in the context of QCD+QED with massless quarks, up to O(α α S ). The requested Feynman diagrams are shown in Fig. 6. Starting at the LO, we have for γ * (p 12 ) → q(p 1 ) +q(p 2 ) with p 12 = p 1 + p 2 and p 2 12 = s 12 > 0, where e and e q denote the electromagnetic coupling and the quark electric charge, respectively § . The corresponding Born level total cross-section is given by where the two-body phase-space factor is shown in Eq. (133) (Appendix B).
Let's centre in the NLO contributions. The real corrections from the radiative process γ * (p 12 ) → q(p 1 ) +q(p 2 ) + g(p r ) are where the virtuality of the photon is s 12 , the final-state momenta are denoted primed to distinguish them from the Born level kinematics, and y ij = 2 p i · p j /s 12 . To obtain the expressions in Eq. (87) we use the expansion of the three-body phase-space shown in Eq. (135). If we expand this result, we find § As usual, the squared matrix elements are averaged over the number of spin degrees of freedom of the incoming particles, which is taken to be 2(1 − ) for the photon. In any case, we normalise the NLO results by the LO contribution. which contains both double and single -poles, associated with soft and collinear singularities, respectively. On the other hand, the virtual contribution is generated by the interference of the one-loop vertex correction with the Born amplitude, which is given by with d ν 1 ν 2 (p) the gluon (or photon) polarisation tensor. In order to obtain this expression, the internal momenta are defined as q 1 = + p 1 , q 2 = + p 12 and q 3 = , and are considered as in the scalar case discussed in Section 4. In the Feynman gauge, the expression in Eq. (89) takes the form and leads to the following dual contributions with |M qq,i . The individual dual integrals in Eq. (91) contain up to quadratic UV divergences, although the naive power counting in the original loop integral in Eq. (89) leads to a logarithmic behaviour in the UV. The quadratic UV divergences cancel in the sum of all the dual integrals, while the linear UV divergences disappear after integration over the polar angle. Applying the change of variables from Eq. (33) in the first dual contribution, however, produces the cancellation of the linear UV divergences also at the integrand level. As expected, only the logarithmic UV divergences remain in the sum of the dual contributions. Performing the explicit integration over the loop variables leads to Putting together the three dual terms, we obtain i.e., we recover the virtual contribution to the total cross-section at NLO. Notice that it was unnecessary to introduce any tensor reduction; Gram determinants are naturally avoided in LTD, and therefore also are the spurious singularities that the tensor reduction introduces leading to numerical instabilities in the integration over the phase-space. Finally, if we sum the contributions from Eq. (87) and Eq. (93), we obtain which agrees with the well-known result available in the literature. The -poles cancel between real and virtual contributions (as expected from the KLN theorem), so we can safely take the limit → 0 after integration.
It is the purpose of this section to show that → 0 can be considered also before integration, once a proper combination of real and virtual terms is done. In the context of LTD, we shall also consider carefully the contributions introduced by self-energy diagrams. On-shell massless quarks do not introduce further corrections to the total cross-section in Eq. (94) due to the renormalisation of the wave function because IR and UV divergences are treated equally in DREG. In the on-shell scheme the wave function renormalisation constant contains both IR and UV divergences, but they cancel each other, which justifies the exclusion of the corresponding Feynman diagrams when carrying out the computation within the traditional approach.
In order to build a complete LTD representation of the virtual contributions, it is required to include the renormalised self-energy corrections to the external particles and properly disentangle IR/UV singularities at integrand level. This step is crucial to achieve a dual representation with a fully local cancellation of singularities, so it can be integrated in four-dimensions. The quark and antiquark self-energies at one-loop are given by with Σ(p i ) = Σ 2 p / i − Σ 1 m. In these expressions, we keep the same internal momenta q i that were used to define the vertex corrections in Eq. (89). From the usual renormalisation procedure, the self-energy contribution is related with the renormalisation factor Z 2 = 1 + ∆Z 2 . Applying on-shell renormalisation conditions to the quark and antiquark self-energies in the Feynman gauge, we obtain the following contributions Formally, these contributions vanish in DREG, but they feature a non-trivial IR and UV behaviour at the integrand level that we need to make explicit in order to have a local cancellation of all the singularities with those present in the real corrections. Applying LTD to the loop integrals given in Eqs. (96) and (97), the corresponding dual representations are where we also kept explicitly the integration variables associated with each cut.
The UV divergences of the wave function cancel exactly the UV divergences of the vertex corrections, because conserved currents or partially conserved currents, as the vector and axial ones, do not get renormalised. To achieve a local cancellation of the UV divergences, it is relevant to note that the vertex corrections diverge logarithmically in the UV, while the expressions in Eqs. (96) and (97) behave linearly in the UV. However, the linear UV divergence cancels upon angular integration. Therefore, a subtraction UV counter-term is needed to cancel locally also the linear singularities. Assuming k UV = 0 (namely q UV = ) and following the discussion of Section 5, we define whose dual representation is given by with m UV = 2µ UV / √ s 12 . Notice that the term proportional to (1 − 2v UV ) integrates to zero and cancels the linear UV singularity. The integration of Eq. (101) leads to Consequently, M qq |Σ(p 1 ) − Σ UV (p 1 ) develops IR singularities only. The same counter-term cancels the UV divergence of the antiquark leg, i.e. M Similarly, it is necessary to remove the UV divergences from the vertex correction included in |M (1) qq . In fact, working in the Feynman gauge, we use the following counter-term where the subleading terms have been fixed such that only the pole is subtracted, i.e.
The crucial observation here is that Eq. (102) and Eq. (104) share the same divergent structure, which allows to cancel completely the UV divergences.
With all these ingredients, we are able to define a four-dimensional representation of the total crosssection at NLO, with a local cancellation of the IR divergences of the loop and the real corrections. In first place, the UV renormalised virtual cross-section is given by with M qq,UV , which also includes the self-energy corrections for simplicity . The renormalised virtual cross-section σ (1,R) V contains only IR singularities at the integrand level, and the UV counter-terms involve a non-trivial integrand level cancellation of UV singularities that must be taken into account in order to find a proper four-dimensional representation of the total cross-section. For this reason, we start by splitting σ where we define as the original virtual terms (including self-energies) and the UV counter-terms, respectively. From Eq. (107), we collect all the dual terms arising when either of the internal momenta q 1 or q 2 are set onshell, and restrict the loop integration by the dual mapping conditions defined in Eq. (48) and Eq. (49), respectively. This leads us to define the virtual dual contributions to the cross-section as together with the dual remnants which fulfill σ Finally, the UV cross-section is given by Following the discussion presented in Section 4, we also separate the real three-body phase-space to isolate the different collinear configurations; this motivates us to introduce σ (1) that are analogous to Eq. (40) for the toy-model previously analysed and fulfil i σ which are obtained by applying the mappings defined in Eqs. (44)- (47) to Eq. (87), respectively. Since the total cross-section is given by we put together σ (1) V,i and σ (1) R,i to define, that are finite in the limit → 0. Moreover, through the application of the momenta mapping and the separation of the integration region, there is a local cancellation of singularities which allows us to take the limit → 0 at integrand level, i.e. before integration.
From Eq. (115), we appreciate that there are still some missing contributions to the cross-section. In fact, we must combine all the virtual terms and UV counter-terms that have not been included in the virtual dual cross-sections in Eq. (109). So, we define as the remnant virtual correction. It is worth appreciating that this expression admits a four-dimensional representation, which can be build applying the change of variables shown in Section 3.1. Finally, the calculation of the total integrals leads to whose sum gives σ that is in agreement with the → 0 limit of the result obtained through DREG. The crucial difference is that we get these results after the integration of four-dimensional expressions ¶ , thus avoiding to deal with complicated expansions.
To conclude this section, let's notice that the procedure shown here could be simplified. In order to achieve an analytical computation it was necessary to split the integration region of the virtual corrections into different pieces. In a numerical implementation, this procedure will not be necessary. The loop integration will occur unrestricted, and the real corrections will be switch on to cancel the collinear and soft divergences in the region of the loop three-momentum where the respective momentum mapping conditions are fulfilled.

Generalisation to multi-leg processes and NNLO
Assuming there are no initial-state partons (for instance in e + e − annihilation) the generalisation to multileg processes is straightforward from the results presented in the previous sections. The implementation of the method in a Monte Carlo event generator is indeed simpler than presented before. As usual, the ¶ Explicit four-dimensional expressions are shown in Appendix C.
NLO cross-section is constructed from the one-loop virtual correction with m partons in the final state and the exclusive real cross-section with m + 1 partons in the final state where the virtual contribution is obtained from its dual representation In the above equation, M N is the N -leg scattering amplitude at LO, with N > m, and M (1,R) N is the renormalised one-loop scattering amplitude, which also contains the self-energy corrections of the external legs. On the other hand, O N ({p j }) defines a given IR-safe physical observable by constraining the integration domain, for example, a jet function. The delta functionδ (q i ) symbolises the dual contribution with the internal momentum q i set on-shell. By renormalised, we mean that appropriate UV counter-terms have been subtracted locally, according to the discussion presented in Section 5, including UV singularities of degree higher than logarithmic that integrate to zero, and the UV contributions from the wave function of the external particles, in such a way that only IR singularities arise in dσ In Eq. (121), we have also assumed a definite ordering of the external particles that leads to a definite set of internal momenta q i . This means, that the one-loop scattering amplitude M (1,R) N contains not only the contribution from the maximal one-loop N -point function, but also all the terms that can be constructed with the same set of internal momenta. Keeping this ordering is necessary to preserve the partial cancellation of singularities among dual contributions at integrand level. Obviously, all the possible permutations and symmetry factors have to be considered to obtain the physical cross-section.
The real cross-section is given by where the external momenta {p j }, the phase-space and the tree-level scattering amplitude M (0) N +1 are rewritten in terms of the loop three-momentum (equivalently, the internal loop on-shell momenta) and the external momenta p i of the Born process. The momentum mapping links uniquely the soft and collinear states of the real and virtual corrections. Therefore, if the physical observable is IR-safe then O N +1 → O N in all the possible IR-degenerate configurations of the (N +1)-particle process. In this way, it is guaranteed that the simultaneous implementation of the real-emission terms with the corresponding dual contributions leads to an integrand level cancellation of IR singularities.
Analogously to the dipole method [5,6], in order to construct the momentum mapping between the m and m + 1 kinematics, we single out two partons for each contribution. The first parton is the emitter and the second parton is the spectator. The difference with respect to the dipole formalism is that both the emitter and the spectator are initially related to external momenta of the virtual scattering amplitudes, and not to internal or external momenta of the real emission processes. Then, the loop three-momentum and the four-momenta of the emitter and the spectator are used to reconstruct the kinematics of the corresponding real emission cross-section in the region of the real phase-space where the twin of the emitter decays to two partons in a soft or collinear configuration. Explicitly, if the momentum of the final-state emitter is p i , the internal momentum prior to the emitter is q i and it is on-shell, and p j is the momentum of the final-state spectator, then, the multi-leg momentum mapping is given by The incoming initial-state momenta, p a and p b , are left unchanged. Momentum conservation is automatically fulfilled by Eq. (123) because p i + p j + k =i,j p k = p i + p r + p j + k =i,j p k , and all the primed final-state momenta are massless and on-shell if the virtual unprimed momenta are also massless. The momentum mapping in Eq. (123) is motivated by the general factorisation properties in QCD [19,56], and it is graphically explained in Fig. 7.
The emitter p i has the same flavour as p ir , the twin emitter or parent parton (called emitter in the dipole formalism) of the real splitting configuration that is mapped. The spectator p j is used to balance momentum conservation, and has the same flavour in the virtual and real contributions. As for the dipoles, there are alternatives to treat the recoiling momentum, but the option with a single spectator is the most suitable. The radiated particles, p i and p r , might have a different flavour than the emitter p i . If the emitter is a quark or antiquark, the role of p i and p r can be exchanged if the subindex i is used to denote the flavour, as we did in Eq. (46). If the emitter is a gluon, the radiated partons are gluons or a quark-antiquark pair.
The momentum mapping in Eq. (123) is suitable in the region of the loop-momentum space where q i is soft or collinear with p i , and therefore in the region of the real phase-space where p i and p r are produced collinear or one of them is soft. In Eq. (121), we have introduced a complete partition of the real phase-space which is equivalent to divide the phase-space by the minimal two-body invariant y ir . Since the real and virtual kinematics are related, the real phase-space partition defines equivalent regions in the loop three-momentum space. Notice, however, that we have not imposed these constrains in the definition of the virtual cross-section in Eq. (121), as we did for the analytic applications presented in Sections 4 and 6. The actual implementation of the NLO cross-section in a Monte Carlo event generator is a single unconstrained integral in the loop three-momentum, and the phase-space with m final-state particles. By virtue of the momentum mapping, real corrections are switch on in the region of the loop threemomentum where they map the corresponding soft and collinear divergences. This region is a compact region, and it is of the size of the representative hard scale of the scattering process. At large loop threemomentum only the virtual corrections contribute, and their UV singularities are subtracted locally by suitable counter-terms. These are the only counter-terms that are required for the implementation of the method, the IR singularities are unsubtracted because their cancellation is achieved simultaneously. The full calculation is implemented in four-dimensions, more precisely with the DREG parameter = 0. Moreover, there is no need to perform any tensor reduction in the calculation of the virtual contributions, hence avoiding the appearance of Gram determinants leading to spurious numerical instabilities. Inte-We assume that the internal and external momenta are ordered according to Fig. 1. grable threshold singularities of the loop contributions are also restricted to the physical compact region, and are treated numerically by contour deformation [21,22,23,24] in a Monte Carlo implementation.
The case of lepton-hadron and hadron-hadron collisions deserves a comment apart. The crosssection is computed by convoluting the corresponding partonic cross-section with the process-independent parton distributions of the incoming hadrons. Since the initial-state partons carry a well defined momentum, the partonic subprocesses are not collinear safe. By virtue of the universal factorisation properties of QCD for massless incoming partons [57], the initial-state collinear singularities are factorised and reabsorbed in the definition of the non-perturbative parton distributions, and are removed from the partonic cross-section by suitable collinear counter-terms which are proportional to 1/ and the Altarelli-Parisi splitting functions [55,58,59,60,61,62]. The initial-state collinear counter-terms are convoluted with the Born cross-section. However, their standard form is not suitable in the LTD approach because that convolution is a single integral in longitudinal momentum with Born kinematics. Its unintegrated form, which should depend also on the transverse momentum of the real radiation, is necessary, and will be presented in a future publication.
The extension to NNLO, although not straightforward, can instinctively be anticipated and will be treated also in a future publication. The NNLO cross-section consists of three contributions where the double virtual cross-section dσ (2) VV receives contributions from the interference of the two-loop with the Born scattering amplitudes, and the square of the one-loop scattering amplitude with m finalstate particles, the virtual-real cross-section dσ (2) VR includes the contributions from the interference of one-loop and tree-level scattering amplitudes with one extra external particle, and the double real crosssection dσ (2) RR are tree-level contributions with emission of two extra particles. The LTD representation of the two-loop scattering amplitude is obtained by setting two internal lines on-shell [16]. It leads to the two-loop dual components M   N (δ (q j )) . In both cases, we have at our disposal two independent loop three-momenta and m final-state momenta, from where we can reconstruct the kinematics of the one-loop corrections entering dσ (2) VR , and the treelevel corrections in dσ (2)

Conclusions and outlook
In this article, we carefully discussed the implementation of a novel algorithm to compute higher-order corrections to physical observables. This method is based on the LTD theorem, which states that virtual contributions can be expressed as the sum over single-cut (at one-loop) or dual integrals, whose structure closely resembles real-emission amplitudes. We exploit this knowledge to perform an integrand-level combination of real and virtual terms, which leads to a fully local cancellation of singularities and allows to implement the calculation without making use of DREG.
One of the interesting uses of LTD lies in the possibility of exploring the causal structure of virtual contributions and disentangling their singularities. In particular, we applied this technique to prove that the IR divergences are generated in a compact region of the loop-momentum space. This is a crucial fact in order to achieve the real-virtual cancellation of singularities in IR-safe observables.
To illustrate the importance of the compactness of the IR singular regions, we studied the NLO corrections to a 1 → 2 process in the context of a toy scalar model. By using suitable momentum mappings, we generated 1 → 3 on-shell massless kinematics from the 1 → 2 process and the loop three-momentum. These mappings relate exactly the integration regions where the singularities are originated. In this example, we distinguished two regions in the real-emission phase-space and defined suitable mappings to cover them in the dual space of the loop three-momentum. In this way, the combination of dual and real contributions led to expressions that are integrable in four-dimensions; the remainders of the virtual part were also represented by a four-dimensional integral. The algorithm is called unsubtracted because the summation over degenerate soft and collinear final-states is performed thanks to these momentum mappings, then making unnecessary the introduction of IR subtractions.
On the other hand, we investigated the cancellation of UV divergences at integrand level. We started with the simplest example of a massless scalar two-point function. Using the ideas presented in Ref. [36], we obtained the dual representation of the local UV counter-term that exactly cancels the divergences in the high-energy region of the loop momentum, rendering an integrable representation in four-dimensions. We also extended the procedure to deal with arbitrary scattering amplitudes, and provided a subtle physical interpretation of the energy scale entering the UV counter-term as renormalisation scale.
Then, we applied the IR unsubtraction LTD algorithm to the physical process γ * → qq(g) at NLO in QCD. In the first place, we found the dual representation of the virtual contribution and made use of the momentum mappings to perform the real-virtual combination. It is worth appreciating that we had to include also self-energy corrections to the external on-shell legs, while they are usually ignored because their integrated form vanishes in DREG due to the lack of physical scales. Since our approach explicitly splits the IR and UV regions in the dual integrals, it was necessary to disentangle their IR/UV behaviour. In this way, we computed the full NLO correction to γ * → qq(g) by making use of purely four-dimensional expressions also in a physical application.
The generalisation of the algorithm to deal with multi-particle processes is quite straightforward, at least when only final-state singularities take place. Essentially, the real emission phase-space must be split to isolate the different collinear configurations and, then, a proper momentum mapping for each of these configurations is defined. Combining the dual integrands with the real matrix elements in the corresponding regions leads to integrable expressions in four-dimensions. The cancellation of UV divergences is done following the same ideas as in the 1 → 2 case. To conclude, we succinctly sketched the extension of the algorithm to cure initial-state collinear singularities, and discussed briefly its generalisation to NNLO.
In summary, in this work we presented explicitly a well-defined algorithm that allows to override DREG by exploiting LTD. It constitutes a new paradigm in perturbative calculations because it takes advantage of combining directly real and virtual corrections in an integrable four-dimensional representation, while providing an easy physical interpretation of the singularities of the scattering amplitudes and unveiling their hidden nature.
Since in this work we deal with massless internal propagators, then q (+) i,0 = q 2 i − ı0 and the dual integration measure simplifies to which justifies the notation established in Eq. (13).
B Phase-space Following the notation introduced in Eq. (6), the d-dimensional phase-space associated to a N -leg scattering process with m final-state particles is given by the expression with p a the sum of the incoming momenta for either 1 → m or 2 → m processes. The total volume of the 1 → 2 phase-space is given by the well-known formula where s 12 = p 2 12 is the invariant mass of the decaying particle. For the NLO real-radiation correction to 1 → 2 process we need to deal with 1 → 3 decays; the corresponding phase-space is dΦ 1→3 = s 12 2 (4π) 3−2 Γ(2 − 2 ) s 12 µ 2 −2 (y 12 y 1r y 2r ) − dy 1r dy 2r , where the parametrisation used is the one introduced in Eqs. (45) and (47). Because we are interested in the relative factors of dΦ 1→3 compared with dΦ 1→2 , Eq. (134) can be rewritten as dΦ 1→2 (y 12 y 1r y 2r ) − dy 1r dy 2r .
Notice that the kinematical configuration in the TL-region is given by p 3 → p 1 +p 2 (or p 3 → p 1 +p 2 +p r for the real term), with p 3 a massless off-shell particle whose virtuality is p 2 3 = s 12 > 0. The integrand is considered Lorentz invariant, so we integrated out angular degrees of freedom associated with the extra-dimensions.
Finally, let's recall another useful technique for phase-space integration and -expansion of the results, which relies in the application of distributions [4] to rewrite the DREG measure. The basic formula is given by where x C ∈ (0, 1] is an arbitrary cut and f is a generic test function without singularities at x = 0 and without -poles. The C-distributions are defined according to i.e. the cancellation of the integrand is forced in a neighbourhood of the singular point x = 0. Notice that the test function must be an entire function of in order to avoid additional -poles. Besides that, Eqs. (136) and (137) can be adapted to several domains and measures. In particular, the -expansion of the phase-space at large loop momentum becomes expressible as under the assumption that f (q i , ) is an entire function of and it has a vanishing soft limit (i.e. f → 0 as q i,0 → 0). Appreciate the presence of extra-contributions given by the collinear residues at v = 0 or v = 1, as well as some additional terms introduced by the linear dependence of the integrand.
C NLO corrections to γ * → qq(g): explicit 4D formulae In this Appendix, we collect the four-dimensional representation of the integrands associated to the integrals in Eq. (118). Explicitly, we have for the real-virtual combinations, and the dual virtual remnant. In the previous expression, we have identified all the integration variables, ξ = ξ 2,0 = ξ 3,0 = ξ UV and v = v 2 = v 3 = v UV , while (ξ 1,0 , v 1 ) are expressed in terms of (ξ 3,0 , v 3 ) by using the change of variables in Eq. (33). The integration regions are defined through Eq. (57), and we use Eq. (56) to simplify the analytic integration. Notice that the integrand of the dual virtual remnant behaves as d σ in the high-energy limit, and as O(ξ −3 ) after angular integration, thanks to the UV counter-term that subtracts up to linear UV divergences locally.