Dilepton production through timelike Compton scattering (TCS) within the $k_T$-factorization approach

In this work we consider the dilepton production via timelike Compton scattering (TCS) in electron-proton and proton-proton collisions. In particular, the differential cross section in terms of the dilepton invariant mass and rapidity is computed within the $k_T$-factorization approach. Besides, we utilize distinct unintegrated gluon distributions (UGD) in order to compare their impact on the differential cross section of TCS in $pp$ ($ep$) collisions evaluated at the LHC (LHeC), HL-LHC (LHeC), HE-LHC (LHeC) and FCC-hh (eh) center-of-mass energies.


I. INTRODUCTION
Dilepton production can occur through several mechanisms, being the leading one the ordinary Drell-Yan process. The second most important contribution comes from photon fusion, i.e., γγ → + − , which is used for controlling the luminosity at the LHC. In addition, single and double diffractive Drell-Yan also produce dileptons via different interactions, such as Pomeron -Pomeron (IP IP ), Pomeron -Reggeon (IP IR), Reggeon -Reggeon (IRIR), Pomeron -proton (IP p) and Reggeon -proton (IRp) reactions [1,2]. We can still have the reactions γIP and IP γ, where the underlying process is the timelike Compton scattering (TCS). At last, it should be mentioned that the Bethe-Heitler (BH) mechanism contributes at the amplitude level to the physical process of photoproduction of heavy lepton pairs, and it is known that the BH contribution (and its interference with TCS) is large in contrast to timelike Compton scattering itself.
Timelike Compton scattering has been commonly investigated within the formalism of Generalized Parton Distributions (GPDs) [3][4][5] (see also, for example, Refs. [6][7][8]). One of the goals in the study of these distributions is to understand how quarks and gluons assemble themselves to hadrons [9][10][11][12]. Since the cleanest reactions to obtain the GPDs are the DVCS (deeply virtual Compton scattering) and TCS, studying the latter through distinct reactions could be relevant for their determination. Being the "inverse "process of the former, in TCS a quasireal photon interacts with a proton and the final state after the scattering is an outgoing proton and a timelike virtual photon that subsequently decays into a lepton pair. Namely, the process is the following: γp → γ * p. In the context of GPDs, recently the TCS amplitudes and associated observables have been investigated in leadingtwist approximation [13] and a careful analysis is done in order to substantially reduce the model dependence.
The dilepton production from TCS was adressed in ultraperipheral collisions (UPCs) at a fixed-target experiment (AFTER@LHC) using the nucleon and ion beams in Ref. [14]. Moreover, the linearly polarized photon beam has been considered in [15], where new observables were proposed and the impact on determination of polarized GPDs has been studied. Yet, the NLO corrections to the timelike (TCS), spacelike (DVCS) and double deeply virtual Compton scattering (DDVCS) amplitudes have been fully demonstrated in [16].
The process has been also investigated within the color dipole approach. In Ref. [17], the cross section was first computed using a spacelike approximation for ep and eA collisions. It was a straightforward application of a previous work on the diffractive photoproduction of Z 0 done in Ref. [18] and it is complementary to the predictions for nuclear DVCS [19]. The comparison between the inclusive and exclusive dilepton photoproduction was done in Ref. [20]. The wave function for an outgoing photon with timelike q 2 > 0 was derived in [21]. It was shown that the cross section calculation involves a strong oscillatory integrand, which was solved by taking analytic continuation to complex transverse dipole size, r, with a suitable integration contour. This difficulty does not appear if the transverse momentum space is considered. The TCS process is deeply connected with the DDVCS process, γ * p → γ * p → ( + − )p, which was treated in the context of dipole framework in [22] considering the scattering on nucleons and nuclei.
It is well known that, at asymptotically high energies, BFKL (Balitsky-Fadin-Kuraev-Lipatov) [23,24] evolution describes the gluon dynamics. The corresponding evolution equation describes the x behavior of the unintegrated distribution. The results coming from DGLAP (Dokshitzer-Gribov-Lipatov-Altarelli-Parisi) [25][26][27] evolution coincide with those from BFKL in the double logarithmic limit. Both approaches predict strong rise on F 2 at small Bjorken variable x, as measured at HERA. However, BFKL evolution predicts a strong power-like rise. Similarly to the collinear factorization, one can factorize an observable into a convolution of processdependent hard matrix elements with universal parton distributions, but here the virtuality and transverse mo-mentum of the propagating gluon are no longer ordered and then the matrix elements need to be taken off-shell, which leads to the fact that the convolution also occurs over transverse momentum through unintegrated parton distributions (uPDF). This formalism is the so-called k Tfactorization approach [28]. Hereafter, in most cases we will refer to uPDFs simply as UGDs (unintegrated gluon distributions) since at the small-x regime the dominant partons are the gluons.
At small-x (high energies), since ∆y ∼ log(1/x), the evolution of parton distributions go along over a large region of rapidity. In this context, the effects of finite transverse momenta of partons can become progressively important. Thus, the cross sections may be factorized into a k T partonic cross section and an unintegrated parton distribution, φ(z, k 2 T ), i.e., [28] In this work, we focus on TCS and corresponding dilepton production in both ep and pp reactions within the k Tfactorization approach. In Ref. [29], the process was first calculated in such a formalism and the differential cross sections for dilepton production as a function of invariant mass and energy have been analysed for electron-proton scattering process. There, the authors utilized an unintegrated gluon distribution proposed in Ref. [30].
There are not many studies regarding TCS in literature so far. Therefore, this subject has not been substantially explored yet and, in that sense, further investigations could be highly relevant. The aim here is to extend the analysis carried out in Ref. [29] by taking into account other UGDs in order to single out the model dependence and perform predictions for the future high energy ep colliders LHeC, HL-LHeC, HE-LHeC and FCC-eh. The main goal is to compare the differential cross sections (from different UGDs) with respect to the dilepton invariant mass and rapidity distributions, as well as the total production cross section. Furthermore, we also do the investigation for proton-proton collisions, in which one of the mechanisms for dilepton production is the photon-Pomeron interaction where TCS is present. This paper is organized as follows. In the next section, we summarize the derivation of the cross section for dilepton production by TCS in the transverse momentum space for electron-proton and proton-proton collisions. In Section III, we introduce the different UGDs that will be utilized in this work and apply them to the calculations, discussing the model dependence in the expressions. This aims to understand the theoretical uncertainties and propose observables to be measured in future high energy and high luminosity ep and pp machines. Finally, in Section V we outline the paper and expose our main conclusions. Initially, we address the dilepton production through TCS in electron-proton collisions (γp → γ * p). We will adopt the formalism proposed in Ref. [29], where the imaginary part of the TCS amplitude is calculated in terms of the unintegrated gluon function within the k Tfactorization approach. The underlying process is the color dipole qq interaction with the proton producing an exclusive final state where a QCD Pomeron is exchanged in t-channel. The amplitude is given below: where k 2 ⊥ stands for the transverse momentum squared of the gluons, while κ 2 represents the transverse momentum squared of the quark. In the expression above, z is the longitudinal momentum fraction carried by quark and φ(x, k 2 ⊥ ) is the UGD. The running coupling, α s (µ 2 ), is being taken at µ 2 = max(κ 2 + m 2 f , k 2 ⊥ ). Adopting the prescription of Ref. [30], if α s exceeds 0.82, it is frozen at this value in order to assure a perturbative scale. The functions C 0f , D 0f , C 1f and D 1f are defined as follows: where m f is the quark mass of flavor f and α and β are given by: Having Im M f , we shall define the spectral distribution, which is regarded to the diffractive amplitude for the γp → qqp (the virtual photon is being taken as a quark-antiquark pair) transition as: where κ 2 max = (0.25M 2 qq − m 2 f ). In Eq. (8), W 2 and M 2 qq represent the photon-proton center-of-mass energy squared and the dipole invariant mass squared, respectively. The latter is related to κ 2 and z by M 2 qq = (κ 2 + m 2 f )/z(1 − z). The UGD is taken at x = δ M 2 qq W 2 in order to correct the skewedness effects, where δ = 0.41. Here we are following closely Ref. [29], but different prescriptions for the skewedness corrections could be employed. For instance, a skewedness factor R g using the Shuvaev et al. expression for gluons [31] can be multiplied to the amplitude or the prescription of Harland-Lang [32] can be considered, where the skewed gluon density is simply related to the gluon GPD. The offdiagonal correction is one of the theoretical uncertainties in the calculations.
The TCS scattering amplitude is computed as: The integration in the domains of M 2 + − > 4m 2 f and M 2 + − < 4m 2 f is written in the following way: where PV stands for the Cauchy Principal Value and the auxiliary function ξ has been defined as: At this point, some considerations are in order. The upper bound of the integral in Eq. (9) leads to contributions x 1 or so (see the definition of x). These contributions are suppressed by the 1/M 2 qq factor in the spectral distribution, Eq. (8), and by the large x threshold factor (1 − x) n present in the phenomenological UGDs.
Regarding Eq. (10), Θ denotes the Heaviside function and e 2 f is the squared quark charge of flavor f . Analogously, the real part of the amplitude is evaluated by In the previous expression, the definition of The function ReM f is obtained by using the dispersion relation, . Taking into consideration the definition of the variable x alongside Eq. (3), the derivative in the above expression is easily obtained. For simplicity, the diffraction cone approximation will be used, which enables one to embed a t dependence by means of the following factorization: where B is the slope parameter. In this work, we will adopt B = 4 GeV −2 . By utilizing the optical theorem, evaluating the total cross section for the γp → γ * p process is straightforward, where ImA T CS = ImA T CS f and ReA T CS = ReA T CS f Finally, we can express the differential cross section in terms of the dilepton invariant mass distribution, i.e, Here, we will only consider the charm and the light quarks, whose values are taken according to the corresponding UGD model applied. Lastly, one may integrate Eq. (16) in order to get the cross section integrated over the dilepton invariant mass, M 2 + − , i.e., where (M 2 + − ) min is the (cut) minimum invariant mass of the lepton pair.

B. TCS in pp collisions
In case of pp collisions, the production of lepton pairs via TCS is carried out by photon-Pomeron (γIP ) and Pomeron-photon (IP γ) mechanisms. Following Ref. [1], the γIP + IP γ contribution for the amplitude of the pp → pp + − process may be obtained within the equivalent photon approximation (EPA) and reads as where k is the photon energy, dn(k)/dk is the photon flux and y pair is the dilepton rapidity. The subscripts + and − are related to the γIP and IP γ subprocesses, respectively. The flux expression will be extracted from Ref. [33], in which it is written as The quantity √ s is the center-of-mass energy of the pp system and the parameter χ is defined as Given the definitions of rapidity and s (Mandelstam variable), one can express the following relations: The expression above relates the photon-proton centerof-mass energy to the proton-proton one. From the experimental point of view, the production of exclusive dilepton events is relatively understood. For instance, ATLAS collaboration has recently performed measurement at 13 TeV for a dimuon invariant mass of 12 < M + − < 70 GeV [34]. Also, CMS collaboration [35] has measured proton-tagged events at the same energy for exclusive dilepton produced at midrapidity with M + − > 110 GeV and one of the two scattered protons is measured in precision proton spectrometer (CT-PPS). ATLAS has reported similar measurement of forward proton scattering in association with dileptons produced via γγ fusion with a significance higher than 5σ [36]. On the phenomenological side, in the new Super-Chic 4 Monte Carlo [37] photon-initiated production in proton-proton collisions has been implemented. The code takes into account the different contributing channels, including proton dissociation.

III. RESULTS AND DISCUSSION
As pointed out in the Introduction, we aim to calculate the TCS process in ep and pp collisions within the k T -factorization approach. To do so, an unintegrated gluon distribution is needed as the non-perturbative input of the formalism. In that sense, in this work we will  consider four UGDs containing different physical informations. Initially, we take the KS (Kutak-Sapeta) UGD, specifically its non-linear set [38], which takes into account parton saturation effects. The KS distribution was imported from the TMDlib (Transverse Momentum Dependent parton distributions library) [39], which provides a large number of uPDFs. Due to non-positive definite kernel, the basic formulation of the NLO BFKL equation is unstable. In order to stabilize it, one should resume a subset of higher order corrections. In Ref. [40], the authors took the higher order corrections from consistency constraint on emission of real gluons. In addition, more corrections are performed by running the constant coupling and other contributions come from non-singular pieces of the DGLAP splitting functions. In this context, the authors of Ref. [38] proposed the KS model, whose expression may be seen in the quoted reference.
The next UGD considered has been developed in Ref. [41] (we will name it as MPM hereafter), where the authors utilized the geometric scaling framework to build an expression for the gluon unintegrated function that depends on the variable τ , being τ = T is the transverse momentum of the gluons and Q 2 s is the saturation scale. Along with the variable τ , the MPM has also three other parameters (see the quoted reference for details). Therein, in order to avoid the divergence of jet production in the infrared sector (IR), the saturation scale is taken as an effective regulator of the gluon propagator compatible with a Yukawa potential, φ(k 2 T ) ∼ α s k 2 T /(1 + k 2 T /µ 2 ), leading to a distribution of the following form: where Q 2 s = (x 0 /x) 0.33 and δn = aτ b . The parameters σ 0 , x 0 , a and b were fitted against DIS data for x < 0.01. The model describes simultaneously the DIS data at small-x and the spectra of produced hadrons in pp/pp collisions. The MPM model is based on geometric scaling arguments and Tsallis-like behavior of the measured spectra. Furthermore, it has no dependence on the scale µ 2 and a coupling constant, α s = 0.2, is assumed.
The third UGD is based on the GBW (Golec-Biernat-Wüsthoff) parametrization [42]. Having parameters fitted from DIS data at small-x, the expression is analytical  and given by [42]: Above, the variable Q s is the saturation scale and its value is Q 2 s = (x 0 /x) λ , while σ 0 = 27.32 mb, λ = 0.248 and x 0 = 4.2 × 10 −5 [43]. The model above holds the small-x region and then presents the geometric scaling property with dependence on the ratio k 2 T /Q 2 s . Likewise the MPM model, the GBW parametrization has no dependence on µ 2 and α s = 0.2.
Finally, we also accounted for the UGD proposed in Ref. [30] (we will call it IN hereafter), which is the dis- tribution used in [29] . This parametrization has been separated into two parts, namely the soft and hard ones. The latter is also divided in two underlying contributions. For large Q 2 , the UGD is simply taken as the derivative of the gluon PDF (for instance, GRV, CTEQ, MRS...) with respect to ln Q 2 . On the other hand, for k 2 ⊥ 1/R c (κ 2 is the gluon transverse momentum), the UGD dependence on κ 2 is similar to the Yukawa screened flux of photons in the positron (see Eq. (4) of Ref. [30]). The quantity Q c is inversely proportional to the screening radius, Q c ∼ 1/R c , with R c ∼ 0.27 f m. This variable denotes the propagation/screening of the perturbative color fields (Yukawa-Debye screening). Regarding the soft part of the ansatz, it considers large dipoles in which the dipole cross section does not depend on the energy. In order the verify the explicit expression of this UGD and the corresponding detailed analysis, see Ref. [30] and Eq. (43) therein. The fact is that the IN distribution is not based on saturation physics arguments and the transition between the hard and soft regimes is set by the fixed IR scale Q c . The model is quite successful and describes correctly the structure functions F 2 , F L , F cc 2 and gives reasonable predictions for exclusive vector meson production [44][45][46], as well as the exclusive jet production [47]. Fig. 1, where they are presented as a function of gluon transverse momentum, k T , for fixed values of Bjorken variable: x = 10 −2 , 10 −3 , 10 −4 , 10 −5 . It is noticed that the behavior in k T > 1 GeV is somewhat similar for the MPM, KS and IN models despite their distinct overall normalizations. In these UGDs that region is driven by the DGLAP-like behavior for the integrated gluon distribution. The GBW model presents the already knonw fall-off for large momentum. The transition to the soft region at small k T is different in each model. In the IN parametrization it occurs for a fixed momentum value, k T = Q c , whereas it is dynamical for models based on saturation formalism. The critical line is established by the saturation scale, Q s (x). We anticipate that the predicted cross section will have significant variability in their overall normalizations.

A comparison between the considered UGDs -MPM (solid lines), GBW (dashed lines), KS (dot-dashed lines) and IN (dotted lines) -is done in
As already mentioned, the calculations for ep collisions will be performed at the center-of-mass energies of the proposed future facilities as the High-Luminosity LHeC (HL-LHeC), its high energy update (HE-LHeC) and the Future Circular Collider in the lepton-hadron mode (FCC-eh). Their values of collision energies are outlined in Table I [48]. In Fig. 2, the differential cross sections of dilepton invariant mass, Eq. (16), are presented for the design energies of the projected experiments mentioned. Furthermore, we calculated the integrated cross section for the γp → + − p process over an interval between (M 2 + − ) min and infinity. The values as a function of minimum dilepton invariant mass, (M 2 + − ) min = 0.1 GeV 2 and (M 2 + − ) min = 1 GeV 2 are summarized in Table II, respectively. We also provide the event rates per year for M 2 + − > 1 GeV 2 using the design luminosities in each energy [48]. On average, the behavior for MPM, GBW and KS are quite similar in the dilepton mass interval considered. The IN model produces a higher cross section with sizable   deviation compared to the other UGDs. This fact and the matter of GBW prediction being very close to the remaining UGDs mean that this process is dominated by transverse momentum around critical line, Q s , or so. We call attention once again for the difference on the transition hard-soft in the IN model which occurs at a fixed scale having order of magnitude of a few GeV. At very high energies, the saturation scale is enhanced a lot, Q s (W ) ∼ (W/W 0 ) 0.12 GeV, and therefore in the saturation models (or saturation inspired parametrizations) there is an increasing contribution from transverse momenta in the region k T Q s (W ).
The presented cross section can be directly compared to previous calculations in literature using the color dipole formalism or k T -factorization. Within the color dipole picture in the spacelike approximation [17,20], it was found that the invariant mass distribution is driven by the small mass region and the forward amplitude scales with ∼ Q 2 s (x)/M 2 + − when parton saturation models are considered with the critical line being defined by the saturation scale, Q s (x). This leads to a differential cross section having the behavior d 2 σ/dM 2 The integrated cross section was predicted to be 0.08 nb using the cut M + − ≥ 1.5 GeV [20] for W γp = 1.4 TeV. Having in mind that the spacelike cross section is 3 or 4 times smaller than the timelike one, that calculation is similar to ours for the HL-LHeC energy using the models based on saturation physics. Interestingly, the inclusive dilepton photoproduction has been also estimated in Ref. [20], where it was found to be σ inc = 0.78 nb for the integrated cross section integrated for M + − ≥ 3 GeV. The first calculation using k T -factorization was performed in Ref. [29] using only the IN distribution. The present calculation is fully consistent with that study, with the integrated cross section for HL/HE-LHeC energies being ∼ 0.5 nb for (M 2 + − ) min ≥ 1 GeV. Moving now to pp collisions, Fig. 3 presents the dilepton invariant mass distribution from TCS process at the LHC (Large Hadron Collider) energies of 7 TeV (left panel) and 13 TeV (right panel) at central rapidity (y=0).
In addition, Table III shows the cross sections for the LHC13, HE-LHC (27 TeV) and FCC-hh (100 TeV) [48]. As the predictions for the HL-LHC (14 TeV) are quite similar to the LHC13 ones, we are only displaying the latter. As expected, the behavior follows the pattern verified in ep collisions. Once again, the differential cross section is dominated by dileptons produced with low invariant mass and there is a large dispersion in the absolute value of the cross section when different models are considered. The results can be contrasted with the phenomenological predictions for dilepton production coming from the two-photon channel and hard diffractive channel given by Pomeron-Pomeron interactions. The invariant mass distribution has been predicted for 14 TeV in [1], where it was found to be dσ/dM + − 0.7 nb for M + − = 2 GeV and dσ/dM + − 0.1 pb for M + − = 10 GeV considering the IN UGD. Our results are similar at 13 TeV, where the evaluations using models based on saturation physics give systematically smaller values than the IN parametrization. The same is true for the rapidity distribution, as we will see in what follows. Still, the invariant mass distribution from the twophoton production is around two orders of magnitude higher than the exclusive production through the photon-Pomeron channel [1].
Finally, in Fig. 4 the rapidity distribution for pp collisions at the LHC is displayed for the energy of 13 TeV. At midrapidity, the cross section ranges in dσ/dy(y = 0) 10 − 100 pb and contains huge theoretical uncertainty. This is due to the amplification of deviations coming from different model assumptions for a cross section for an exclusive observable. Namely, the large uncertainty can be traced back to the values of cross sections evaluated in quantities proportional to φ(x, k) squared. In Ref. [1], where only the IN UDG has been utilized, the prediction at 14 TeV is 30 pb for midrapidity. This value is consistent with our calculations in its order of magnitude. The authors in [1] found that the contribution from the process γγ → + − at central rapidity is around 10 nb, whereas the contribution from inclusive single diffractive and central diffractive production of dileptons reaches 1 nb. Similar predictions are also presented in Ref. [2] at 13 TeV, concerning the two-photon and inclusive diffraction channel. In that study, the Forward Physics Monte Carlo (FPMC) has been used and the role played by pair transverse momentum cuts was demonstrated in order to disentangle the exclusive photon-induced production at low-p T from the diffractive sector. Our conclusions about the exclusive dilepton production in pp collisions are the same of those presented in [1], where the cross section for exclusive diffractive production is almost the order of magnitude than that for the central diffractive production mechanism. A comprehensive analysis is needed (combining p T and M + − cuts) in order to demonstrate the feasibility of a experimental measurement.
The same process in pp ultraperipheral collisions has been also investigated within the GPD formalism. The prediction for the integrated cross section using NLO GRVGJR2008 PDFs and the hard scale µ 2 F = 5 GeV 2 is 1.9 pb at 14 TeV [49], which has the order of magnitude similar to our predictions for the saturation models like KS UGD. However, the general trend is the estimates using k T -factorization being higher than those from GPD formalism. The background from the Bethe-Heitler process was estimated to be 2.9 pb [49], which cannot be neglected even in high energy scattering. From the experimental point of view, the dilepton production coming from QED contribution is well understood. For instance, the ATLAS collaboration has measured [50] the cross section at 7 TeV in the electron channel, which is determined to be σ(γγ → e + e − ) = 0.428 ± 0.039 pb, whereas in the muon channel one has σ(γγ → µ + µ − ) = 0.628 ± 0.038 pb (errors summed into quadrature, p T > ∼ 10 GeV and |η| < 2.4). The LHCb collaboration has measured in the dimuon channel the cross section σ(pp → µ + µ − pp) = 67 ± 19 nb [51] (errors summed into quadrature, M > 2.5 GeV and η µ − < 2.4).

IV. SUMMARY AND CONCLUSIONS
In this work we calculated the dilepton invariant mass distribution from the TCS cross section in ep collisions for the LHeC, HL-LHeC, HE-LHeC and FCC-eh energies. Besides, we also evaluated this observable for pp collisions at the LHC7, LHC13, HL-LHC, HE-LHC and FCC-hh, along with the rapidity distribution in the case of the LHC13. It was found that the theoretical uncertainty is quite large when we consider different models for the UGDs including those containing parton saturation effects. We found a deviation around one order of magnitude in the models considered in present study. There are other uncertainties coming from the t behavior of the non-forward amplitude, the ansatz for the skewedness corrections and the prescription for the coupling at very low dipole invariant mass. It is clear that a comprehensive analysis on the kinematics variables should be performed in order to disentangle experimentally the TCS contribution from the similar signal coming from Bethe-Heitler background and also from exclusive dilepton production in two-photon fusion. The exclusive diffraction processes are currently being investigated by CMS and ATLAS collaborations at the LHC [52] and the study presented here is complementary to the usual predictions in two-photon physics and inclusive diffraction (see review in [53]). This interesting subject definitely needs more work and the analysis for nuclear targets is ongoing and will be presented in future contribution.