Electroweak boson production at large transverse momentum in double polarized pp scattering at next-to-leading order accuracy

We compute the next-to-leading order corrections to the cross section for the production of electroweak gauge bosons (Z, W± and γ∗) with large transverse momentum in double (longitudinally) polarized hadronic collisions. The calculation is fully performed in the HVBM scheme within dimensional regularization with a careful treatment of issues arising due to the appearance of γ in polarization projectors and axial couplings. We study the phenomenological consequences of the NLO corrections at the level of both the double polarized cross section and the corresponding asymmetries at RHIC kinematics. ar X iv :1 91 1. 08 23 6v 1 [ he pph ] 1 9 N ov 2 01 9


Introduction
Over the last decades, our understanding of hadron structure has remarkably improved thanks to impressive experimental and theoretical progress. In particular, there is currently much activity aiming at further unravelling the nucleon's spin structure.
As is well known, the total quark and anti-quark spin contribution to the nucleon spin was found to be only about ∼ 25% [1,2,3], information that has been mostly provided by Deep-Inelastic scattering (DIS) data.
One emphasis is on the determination of the spin-dependent gluon distribution, ∆g, of the nucleon, which ultimately would give the gluon contribution to the nucleon spin. Clear evidence for a non-vanishing polarization of gluons was found in the region of momentum fraction and at the scales mostly probed by the BNL Relativistic Heavy Ion Collider (RHIC) data [1]. But from the available results it is not possible yet to set a precise value on the total gluon contribution to the proton spin and to the flavor decomposition of the quarks † .
Spin asymmetries in high-energy pp scattering can be particularly sensitive to ∆g, for processes where gluons in the initial state contribute already at the lowest order of perturbation theory. Furthermore, the Drell-Yan process, the production of vector bosons [6,5], has been shown to be very relevant to help for quark flavor separation, since quarks couple differently to both W ± and Z bosons [7].
Therefore, one particularly interesting process that satisfies both requirements is the singleinclusive production of large transverse-momentum (q T ) vector bosons, pp → V X, that can help to provide information on both ∆g and the flavour decomposition of the quark contribution.
In order to make reliable quantitative predictions for a high-energy process, it is crucial to determine the next-to-leading order (NLO) QCD corrections to the Born approximation. In general, the key issue here is to check the perturbative stability of the process considered, i.e. to examine to what extent the NLO corrections affect the cross sections and spin asymmetries relevant for experimental measurements. Only if the corrections are under control can a process that shows good sensitivity to, say, ∆g at the lowest order be regarded as a genuine probe of the polarized gluon distribution and be reliably used to extract it from future data.
In the particular case of (large transverse momentum) gauge boson production another perturbative issues appear due to the existence of two physical scales, the transverse momentum q T and the 'mass' Q of the boson (or its virtuality, specially in the case of photon production).
In the large-q T region (q T ∼ Q), where the transverse momentum is of the order of the vector boson mass, the QCD perturbative series is controlled by a small expansion parameter, α s (Q), and calculations based on the truncation of the perturbative series at a fixed order in α s are theoretically justified. In this region, the QCD radiative corrections for the unpolarized cross section are known up to the next-to-leading order (NLO) [8,9,10] in an analytical form and next-to-next-to-leading order (NNLO) corrections were recently obtained in numerical implementations in [13]- [19]. In the polarized case, a first attempt to achieve the cross section at NLO accuracy for virtual photon production was presented in [20], where only the nonsinglet contribution was obtained. The full NLO result for the double polarized cross section, including the relevant gluon initiated channels, was missing so far.
Nonetheless the bulk of the vector boson events (particularly at RHIC with a center-ofmass energy of √ S = 510 GeV) is produced in the small-q T region (q T Q), where the convergence of the fixed-order expansion is spoiled by the presence of large logarithmic terms, α n s ln m (Q 2 /q 2 T ). To obtain reliable predictions, these logarithmically-enhanced terms have to be systematically resummed to all perturbative orders (see e.g. [11], [12] and references therein). But even in that case, the fixed-order calculation becomes essential in order to perform the proper matching with the resummed contribution, affecting the full result even at rather small transverse momentum. Furthermore, since the small-q T region is mostly affected by soft gluon emission, which is independent on the polarization of the emitting parton, the resummation typically affects both polarized and unpolarized cross section in a rather similar way, rendering a very small effect at the level of asymmetries. Therefore, counting with the NLO corrections for polarized hadronic collisions becomes fundamental in both kinematical regimes in order to understand the data produced at RHIC and to extract the corresponding information in terms of polarized partonic distributions. In this paper we compute the NLO corrections to the cross section for the production of gauge bosons (Z, W ± and γ * ) with large transverse momentum in double (longitudinally) polarized hadronic collisions.
The paper is organized as follows. In Sect. 2 we briefly review the procedure to compute the NLO corrections with emphasis on the issues arising due to the treatment of γ 5 in dimensional regularization. In Sect. 3 we present the analytical results for the NLO corrections for gauge boson production (some of the lengthly functions are given in the Appendix). In Sect. 4 we study the phenomenological impact of the NLO corrections at RHIC kinematics, including the measurable double longitudinal asymmetries. Finally, in Sect. 5 we summarize our results.

Inclusive cross section
We calculate the inclusive cross section of electroweak boson production in the framework of perturbative QCD, considering the case of both polarized and unpolarized initial hadrons. The process is described as where h i , i = 1, 2 are the polarized/unpolarized initial hadrons with momenta P i , and V is either a W ± , Z 0 or virtual photon γ * with Q 2 = M 2 V , energy E Q = Q 0 and large transverse momentum q T with respect to the collision axis. The unpolarized σ and polarized ∆σ cross sections are obtained from the sum or difference of helicity-dependent cross sections given by for the unpolarized case, and for the polarized one. Here the superindices +, − denote the helicities of the two incoming hadrons. The inclusive (polarized) unpolarized cross section can be expressed as Here (∆)f h a (x, µ 2 F ) is the (polarized) unpolarized parton distribution function (PDF) of parton a with momentum fraction x in hadron h, probed at the scale µ 2 F , which are given by the combinations In this case, the superindex +, − denotes de helicity orientation of the parton a with respect to the helicity of the parent hadron. The perturbative (polarized) unpolarized cross section (∆)σ a,b (p a , p b , µ 2 F ) corresponds to the hard-scattering partonic process where p i represents the parton momentum. Here collinear singularities due to the radiation of massless partons are factorized out at the scale µ 2 F and included in the scale dependent (polarized) unpolarized PDFs (∆)f h a (x, µ 2 F ). The partonic momenta can be described in terms of the hadronic ones using the relation p i = x i P i . As it is customary, we introduce the Mandelstam variables for both the hadronic and partonic levels S and s are the hadronic and partonic invariant center-of-mass energies squared of the colliding system, respectively, while s 23 is the invariant mass squared of the system recoiling against the boson V . The partonic cross sectionσ a,b has a singular behaviour in the limit s 23 → 0 related to the cancellation of singularities due to soft gluon emission in the recoiling system and virtual gluon infrared singularities. However, the cross section is actually integrable and the 1/s 23 terms can be dealt with the following change of variables in the momentum fraction integration where the new integration limits are given by

Perturbative calculation techniques
The analytical procedures to compute the perturbative cross section are already well established. At the lowest order (LO) only two channels contribute, as depicted in Fig.1: the annihilation process qq → V g (a) and the Compton process qg → V q (b). At the following order (NLO) virtual and real contributions must be considered. In this case, there are four contributing channels: the process initiated by a quark and an antiquark qq → V (qq, q q , gg), the one initiated by two quarks qq → V qq, the process initiated by a single gluon qg → V qg and the two gluon process gg → V qq.
In all processes we evaluate the Feynman diagrams, performing the Dirac traces with the TRACER package [21], and we integrate over the virtual loop and real emission momenta. Ghost graphs are taken into account to simplify the polarization sums of external gluons. Real contributions are computed using traditional partial fractioning and master integrals, while virtual loop contributions are calculated using the FeynCalc [22] and Feynhelpers [23] packages. Dimensional regularization is used to regulate infrared and ultraviolet singularities, working in a space-time of d = 4 − 2 dimensions. As it will be addressed in the following subsections, dimensional regularization requires a special treatment for the Levi-Civita tensor µνρσ and the Dirac matrix γ 5 , which is present in both the electroweak couplings and quirality projectors. The factorization of the mass singularities in the polarized case is also affected by the scheme used to deal with γ 5 .

Treatment of the matrix γ 5 and axial couplings
While dimensional regularization involves working in d-dimensional space-time, the γ 5 matrix and the µνρσ tensor are only well defined in the four-dimensional space-time. A consistent way to treat γ 5 and µνρσ in d-dimensions is the HVBM scheme [24,25], which splits the ddimensional Minkowsky space into the usual four-dimensional one and a (d − 4)-dimensional subspace where, for instance, the (d − 4)-dimensional part of the γ µ matrices, represented aŝ γ µ , commutes with the strictly four-dimensional γ 5 . Calculations in the HVBM scheme are algebraically more involved because for unobserved momenta p and γ µ matrices one needs to take into account both their four-dimensional part (p,γ µ ) and their (d − 4)-dimensional one (p, γ µ ), which follow different algebraic treatment. Integration over 'hat' momentap in the real contributions is performed as established in Ref. [26].
One of the first issues related to the HVBM scheme is the definition of the electroweak vertices. The standard definition of the boson vertices, taking into account the flavours f 1 and f 2 of the involved quarks, can be expressed as where the left and right-handed coupling are given by the following expressions, depending on the type of boson Here θ W is the electroweak mixing angle, τ ± = (τ 1 ± iτ 2 )/2 and τ 3 are the weak isospin Pauli matrices, U is the CKM mixing matrix, and e f is the electric charge of the corresponding quark (e f = 2 3 for u, c, t and e f = − 1 3 for d, s, b). However, since the commutation properties of the combination γ µ γ 5 depends on whether µ ≤ 3 or µ > 3, a symmetric definition of the vertices is needed, which corresponds to the replacement [27,28]. After this substitution, the vertex expression of Eq.(10) in d-dimensions can be rewritten as where the (d − 4)-dimensional part of γ µ cancels within the left and right-handed terms.
In the HVBM scheme the presence of γ 5 in both the symmetrized vertex and the quirality projectors over initial partons lead to anomalous terms of order O(d − 4), which in combination with divergent terms can result in spurious finite contributions. Infrared anomalous terms cancel out between virtual and real contributions [29], but in the case of the renormalization of ultraviolet divergencies and mass factorization appropriate additional counterterms are required.

Renormalization and the axial vertex anomaly
For the renormalization of the virtual contributions to the NLO we adopt de usual M S scheme. In this scheme the running strong coupling constant α s (µ 2 ) at a scale µ 2 satisfies the corresponding renormalization-group equation with the β i factors given by where the sum f runs over all the flavours below the mass threshold. However, within the HVBM scheme, since γ 5 no longer anticommutes with γ µ an additional finite renormalization of the axial quark current needs to be taken into account at NLO [30]. The renormalization constant related to the axial current Z 5 is no longer equal to one. The corresponding one loop expression is a b c Figure 2: Schematic of the parton emission processes.
This term gives rise to the following additional finite counterterms that have to be added to the qq and qg channels in order cancel the spurious axial terms appearing in the virtual contributions where d(∆)σ axial ij /d 3 Q is the axial part of the (polarized) unpolarized Born cross section corresponding to either the qq → V g process or the qg → V q Compton process. This axial cross section is obtained by considering only the terms depending on the axial coupling to the boson V .

Factorization
Mass singularities appear due to the collinear emission of massless partons from one of the incoming partons. These singular terms are detached at a scale µ 2 F to be incorporated into the (polarized) unpolarized distribution functions (∆)f h a (x, µ 2 F ). For the unpolarized case, the results are calculated in the usual M S scheme. However, in the polarized case we use the conventional variation of the M S scheme which takes into account some helicity-conservation violations that arise in the HVBM scheme, particularly in the quark-antiquark annihilation process.
The factorization subtractions ‡ due to collinear emission from an initial parton a, as illustrated in Fig.2 for gluon emission, are given by the following convolution where µ is the scale of dimensional regularization and µ F represents the factorization scale.
Altarelli-Parisi splitting function of a parton a emitting a parton c carrying momentum fraction z, d(∆)σ cb→V X is the (polarized) unpolarized d-dimensional Born cross section of the partonic process cb → V X, and (∆)f ab (z) represents the freedom on subtracting additional finite terms.
In the unpolarized case we rely on the usual M S scheme, which corresponds to f ij (z) = 0, but in the polarized case when the subtraction involves emission of quarks by other quarks, when working with the ∆P qq (z) kernel, we use ∆f qq (z) = −4 C F (1 − z). This is the particular scheme is used in the definition of the polarized parton distribution sets [31].

NLO corrections
In this section we present the NLO contributions to the polarized cross section of electroweak boson production at large transverse momentum in hadronic collisions. For the sake of completeness we also provide the unpolarized results, which are in agreement with those of Gonsalves et al. [10], with the distinction that we corroborated the calculation fully within the consistent HVBM scheme for the axial coupling. The polarized cross section for the virtual photon in the non-singlet case partially agrees with the previous results of Ref. [20] due to an overlooked sign in one of the interferences.

Quark-antiquark annihilation and scattering
The formulas of the inclusive cross section are presented in the same way as in Ref. [10], with the expression separated in contributions with different axial-vector and/or flavour structure. The quark-antiquark cross section formula is given by Here the (∆)B qq (s, t, u, Q 2 ) function corresponds to the (polarized) unpolarized Born diagrams The expression is analogous in the qg case.
. This is related to helicity conservation for massless quarks. The F LR cd (s, t, u, Q 2 ) function is actually equal to its polarized counterpart. The cross section of the process qq → V X can be obtained from Eq.(18) by substituting the corresponding couplings b a c d Figure 5: Diagrams which contribute to the process qq → V qq. The F functions in the cross section formula correspond to these diagrams.

Quark-gluon Compton scattering
The cross section of the Compton process is given by where in this case the (polarized) unpolarized Born terms function reads The (∆)G qg (s, t, u, Q 2 ) terms come from diagrams obtained by crossing from those in Fig.4, while the (∆)V (i) qg (s, t, u, Q 2 ) are again the δ(s 23 ) contributions originating from virtual diagrams contributions, which can be obtained from Fig.3 by crossing, and from the same real contributions. Their explicit formulas are given in the Appendix.
The qg process presents the greatest difference between the polarized and unpolarized cross sections. As a matter of fact, Eq.(22) reveals that at the Born level they have two distinct behaviours: one is symmetric in the variables s and t, while the other one is completely antisymmetric in them.
The three variations of the Compton process can be obtained from Eq.(21) with the following substitutions

Gluon-gluon fusion
The cross section of the gluon-gluon fusion process, which only contributes at NLO, is simply given by where the (∆)G gg (s, t, u, Q 2 ) functions, corresponding to the diagrams obtained by crossing from the ones in Fig.4, are given in the Appendix. The G gg (s, t, u, Q 2 ) and ∆G gg (s, t, u, Q 2 ) functions are fairly similar apart from a sign and a few additional terms.

Quark-quark scattering
Finally, the last NLO contribution is the quark-quark initiated, whose cross section is given by where the (∆)H(s, t, u, Q 2 ) functions correspond to the contributions from the diagrams in Fig.6. The global 1/2 factor is a statistical factor needed to avoid double counting events due to the full phase space integration and quantum number summation over the two quarks in the final state. Most of the H functions are related to the F ones from the quark-antiquark process. The H ac (s, t, u, Q 2 ), H bd (s, t, u, Q 2 ), H ad (s, t, u, Q 2 ), H bc (s, t, u, Q 2 ), H LR ab (s, t, u, Q 2 ) and H LR cd (s, t, u, Q 2 ) functions are equal to their polarized counterparts, while the H LL ab (s, t, u, Q 2 ) and H LL cd (s, t, u, Q 2 ) differ only by a sign. The H bb (s, t, u, Q 2 ) and H dd (s, t, u, Q 2 ) functions, as it is the case with F cc (s, t, u, Q 2 ) and F dd (s, t, u, Q 2 ), have significant differences with their polarized counterparts.
The cross section of antiquark-antiquark annihilation is obtained from Eq.(25) by the substitution With the qq and qq cross sections it is possible to calculate the non-singlet cross section, which is given by The non-singlet polarized cross section for virtual photon production has been calculated by Field et al. in Ref. [20]. However, our results do not fully agree with theirs because they assumed b a c d Figure 6: Diagrams corresponding to the qq → V qq process. The H functions in the cross section formula correspond to these diagrams.
that the ∆H LL ab (s, t, u, Q 2 ), ∆H LL cd (s, t, u, Q 2 ), ∆H LR ab (s, t, u, Q 2 ) and ∆H LR cd (s, t, u, Q 2 ) contributions of the qq process cancel against the corresponding ∆F LL cd (s, t, u, Q 2 ) and ∆F LR cd (s, t, u, Q 2 ) functions of the qq process, which is not correct since they actually differ by a sign and add up to an additional contribution. This error was carried over from a mistake in the original unpolarized EMP calculation [8], which it is pointed out in Ref. [10].

Phenomenological Results
In this section we compute the corresponding cross sections to evaluate the impact of the polarized NLO corrections. We apply them to RHIC kinematics, that is, polarized protonproton collisions at a center-of-mass energy of √ S = 510 GeV. We use the MMHT2014nlo [32] and DSSV [31] PDFs sets for the unpolarized and polarized processes, respectively. We begin by analysing the reduction of theoretical uncertainties. These are estimated by studying the dependence of the complete NLO results on the renormalization and factorization scales due to the truncation of the perturbative series expansion. In Fig.7 we present the polarized q T distributions of on-shell Z boson production at both LO and NLO accuracy at y = 0, with their respective confidence bands. Since in this process we have two physical scales given by Q and q T , it is convenient to define a physical scale in between them given by Q 2 s = (Q 2 + q 2 T )/4. The central curves are obtained by fixing the renormalization and factorization scales to µ R = µ F = Q s . The bands are obtained by varying µ R and µ F simultaneously and independently in the range 0.5Q s ≤ µ R , µ F ≤ 2Q s with the constraint 0.5 ≤ µ F /µ R ≤ 2. As expected in perturbative QCD, the second order corrections yield a significant contribution to the Born cross section while reducing uncertainties considerabily.
For a more detailed analysis on the impact of the two scales, in Fig.8 we study the scale dependence of the cross section at a fixed value of q T . The LO curves show stronger dependence on the scale. In this case the µ R dependence corresponds simply to the variations of α s with the scale, while µ F dependence arises exclusively from the PDFs. At NLO these dependences are corrected by logarithmic terms proportional to the β function and the splitting functions, respectively, yielding more stable results. We do not show the results for other particular values of q T since they are fairly similar. Throughout the rest of this section we will set the two scales equal to Q s . The size of the higher order QCD corrections to hadronic processes are usually presented in terms of the 'K-factor', the ratio between the NLO results over the LO ones. However, in order to make this ratio meaningful in the polarized case, the same NLO-evolved parton densities are used to calculate both NLO and LO quantities [33]. This is due to the low constraints to polarized PDFs available from experimental data, especially to ∆g, which can in some cases give very different results when fits are performed at LO or at NLO. For this reason, polarized K-factors can get artificially large or small when the gluons are involved.
In Fig.9 we present the K-factors for the production of on-shell weak bosons at fixed rapidity y = 0 and y = 1, for both the polarized and unpolarized cases. The polarized K factors reveal pp → Z(q T ) S = 510 GeV q T = 25 GeV y = 0  a larger contribution of the NLO corrections for W production, and they show a stronger variation along q T compared to their unpolarized counterparts. The polarized Z production present a remarkably lower K factor at larger q T . At rapidity y = 1 the K-factors are slightly higher, specially at low q T . In this case the bosons W + and W − show similar behaviours. The q T dependence of the polarized K factors is related to the different behaviour of the particular sub-channels contributing to the process. For instance, the changes in the concavity in the curves for y = 0 (barely visible in the Z boson case) are related to the sign change of the qg-channel contribution. This can be observed in Fig.10 for W + production, where we present the sub-channel polarized cross section ratios. In all cases the cross section is governed by the quark initiated channels at low q T , but the qg-channel becomes the most relevant one at high q T . The q T range at which this dominance begins depends on the particular boson studied and its rapidity, taking place at lower q T values in the case y = 1. The gg contribution, that only arises at order O(α 2 s ), is always negligible at the energy of this process. Finally in Fig.11 we show the longitudinal double spin asymmetries for weak boson production at different rapidities, which is defined as The asymmetries grow with q T , reaching sizeable values, and are larger at higher rapidity. The NLO corrections to the asymmetries are generally low, especially at rapidity y = 1. A similar behaviour has been observed in many other processes, where QCD corrections tend to compensate in the ratios of polarized and unpolarized cross sections. However, in the case of Z production the NLO corrections in the polarized case are clearly bigger than the unpolarized ones, accounting for lower values of asymmetry at large q T . This can already be seen in Fig.9 since the ratios of the asymmetries at the two orders are equal to the ratio of the respective K factors, following the relation A N LO = A LO ∆K/K, and the K factors of Z production are the ones that show larger discrepancies between the polarized and unpolarized cases, especially at y = 0.

Summary
In this work, we have presented the first complete calculation at next-to-leading order in perturbative QCD of the cross section for the production of electroweak gauge bosons with large transverse momentum q T in double polarized hadronic collisions. The calculation was done fully in the consistent HVBM scheme for the treatment of γ 5 and the tensor µνρσ in dimensional regularization, carefully dealing with issues arising from the (d − 4)-dimensional part of the momenta and γ µ matrices. We also kept the unpolarized results, easily obtained as a by-product of the polarized calculation, corroborating the previous well-known results [10], but this time fully within the HVBM scheme. The polarized results are also in agreement with a previous partial calculation of the non-singlet polarized production of virtual photons [20], aside from a missed contribution in their work.
Using our results, we studied in some detail the phenomenological impact of the polarized NLO corrections in weak boson production at RHIC, that is, polarized pp collisions with a center-of-mass energy of 510 GeV. We analysed the scale dependence in polarized Z boson production to assess the reliability of the theoretical predictions. We found that the NLO correction have a greatly reduced scale dependence with respect to the Born results. We also studied the polarized K factors, the ratio between the NLO and LO results, comparing them to their unpolarized counterparts. The polarized K factors show larger variations and a more complex dependence in q T , having a slightly higher value than the unpolarized ones except for the Z boson case. The effect of NLO corrections on the double spin asymmetries is generally low for W production, but visibly affects the Z production at higher q T .

Acknowledgements
We are grateful to Werner Vogelsang for helpful discussions. This work has been partially supported by Conicet and ANPCyT.

A Appendix
In this appendix we present the formulae for the functions associated with the NLO contributions in the cross sections of Eqs. (18), (21), (3.3), and (25). They are expressed in terms of the invariants s, t, u, Q 2 and s 23 , but we also define some frequently occurring transcendental functions and denominators. The denominators are: The first group of transcendental functions involves some recurring logarithms: with µ R and µ F the renormalization and factorization scales, respectively. The rest of the transcendental functions include dilogarithms that appear in the virtual contributions: Finally, some expressions include '+' distributions in the variable s 23 . These terms originate from the soft singularities appearing as poles in s 23 . Within dimensional regularization, these poles are explicitly made manifest with the identity 1 s 23 where A is the upper integration limit in Eq. (8). The '+' distributions yield finite contributions and are integrated with the following rules: where f (s 23 ) is any function that depends on s 23 .

A.1 Quark-antiquark annihilation and scattering
In this subsection we present the functions corresponding to the contributions to the quarkantiquark annihilation and scattering processes in both the polarized and unpolarized cases.
We first present the δ(s 23 ) contributions arising from virtual and real contributions Here the ∆B qq (s, t, u, Q 2 ) is the function representing the Born contribution defined in Eq. (19). The function ∆V (2) qq (s, t, u, Q 2 ) arises from the triangle quark loops diagrams, which only contributes to Z production.
The following function is originated from the quark-antiquark annihilation qq → V gg diagrams in Fig.4: In this function we have implicitly included a 1/2 factor to avoid double counting due to the summation over all the final states of the two indistinguishable gluons, and the t ⇔ u indicates that the whole function needs to be repeated switching the variables t and u. The annihilation process also contributes some of the δ(s 23 ) terms which are already included in ∆V qq (s, t, u, Q 2 ). We now proceed to present the functions corresponding to the qq → V qq process diagrams represented in Fig.5: 12s 2 + 3t 2 + 2tu + 3u 2 + 10s (t + u) ds + 2 2s − s 23 + Q 2 + 5 4s .
The diagrams involved in ∆F aa (s, t, u, Q 2 ) function also contribute to the δ(s 23 ) terms in ∆V qq (s, t, u, Q 2 ). The following interference only contributes to Z production: The following functions have also terms proportional to 1/s 23 , but the expression is in fact finite in the limit s 23 → 0: ∆F bc (s, t, u, ∆F ad (s, t, u, Q 2 ) = ∆F ac (s, u, t, Q 2 ), ∆F bd (s, t, u, Q 2 ) = ∆F bc (s, u, t, Q 2 ), ∆F cc (s, t, u, ∆F dd (s, t, u, Q 2 ) = ∆F cc (s, u, t, Q 2 ).
Finally, the interference between the F c and F d diagrams yield different results depending on whether the quirality of the quark and the antiquark are the same (LR) or opposite (LL). This also makes the LR contribution equal to its unpolarized counterpart, while the LL differs by a sign.
∆F LR cd (s, t, u, Q 2 ) = A.1.2 Unpolarized qq → V X Due to helicity conservation, most of the unpolarized results differ from the polarized one only by a sign, e.g. G qq (s, t, u, Q 2 ) = −∆G qq (s, t, u, Q 2 ), so we list only the functions in which this is not true: F LR cd (s, t, u, Q 2 ) = ∆F LR cd (s, t, u, Q 2 ). (51)

A.2 Quark-gluon Compton scattering
In this section we present the functions corresponding to the quark-gluon Compton processes in both the polarized and unpolarized case. The diagrams are obtained from Figs. 3 and 4 by suitable crossing.
We start with the δ(s 23 ) contributions from virtual and real contributions where ∆B qg (s, t, u, Q 2 ) is the Born function defined in Eq. (22). Here again the ∆V (2) qg (s, t, u, Q 2 ) terms arise from the triangle quark loops diagrams, which only contribute to Z production.
Here the 1/s 23 terms without the corresponding '+' distributions are actually finite in the limit s 23 → 0.

A.2.2 Unpolarized qg → V X
The qg → V X process is the one that presents more differences between the expressions of the polarized and unpolarized cases: where B qg (s, t, u, Q 2 ) is the unpolarized Born function defined in Eq. (22). We note that the virtual contributions presented in Eq.(A4) of Ref. [10], in our case included within V qg (s, t, u, Q 2 ), are missing a minus sign in the finite terms not multiplied by the Born function.

A.3 Gluon-gluon fusion
We now present the functions for the gluon-gluon fusion processes for both the polarized and unpolarized cases. The diagrams are obtained from Fig.4 by suitable crossing.
where in this case the t ⇔ u only applies to the additional terms added to ∆G gg (s, t, u, Q 2 ), which itself is already symmetric in t and u.

A.4 Quark-quark scattering
Finally, we present the functions corresponding to the quark-quark scattering processes for both the polarized and unpolarized cases. They are obtained from the diagrams of Fig.6. Most of them are related to the quark-antiquark results.
For the following functions we use the fact that changing an antiquark in the diagrams in Fig.5 with a quark in Fig.6 also changes its coupling to the boson V by the replacement L ↔ −R † , which accounts for a change in sign in the contributions. This relative sign was overlooked in the results presented in Ref. [20].
The last contributions are exclusive to the quark-quark process, and they are equal to their unpolarized counterparts: ∆H ac (s, t, u, Q 2 ) = ∆H bd (s, u, t, ∆H ad (s, t, u, Q 2 ) = ∆H bc (s, u, t, Q 2 ) = C F − C