Isolated photon-hadron production in high energy $pp$ and $pA$ collisions at RHIC and LHC

We compute the isolated photon production in association with a charged hadron at mid rapidity in $pp$ and $pA$ based on the Color Glass Condensate (CGC) framework of high energy QCD where, for the first time, we incorporate the Sudakov effect of soft gluon emissions. Our results are based on the leading order $q g \to q \gamma$ channel in the CGC framework and confronted with the recent data from RHIC and LHC concerning the angular distributions and out-of-plane transverse momentum distributions. We find that, while the CGC computation alone results in too narrow distributions, with the help of the Sudakov effect, we can get a satisfactory description of the data. With this as a benchmark, we provide predictions for the magnitude of the nuclear effect brought by the phenomena of gluon saturation in the CGC.

We compute the isolated photon production in association with a charged hadron at mid rapidity in pp and pA based on the Color Glass Condensate (CGC) framework of high energy QCD where, for the first time, we incorporate the Sudakov effect of soft gluon emissions. Our results are based on the leading order qg → qγ channel in the CGC framework and confronted with the recent data from RHIC and LHC concerning the angular distributions and out-of-plane transverse momentum distributions. We find that, while the CGC computation alone results in too narrow distributions, with the help of the Sudakov effect, we can get a satisfactory description of the data. With this as a benchmark, we provide predictions for the magnitude of the nuclear effect brought by the phenomena of gluon saturation in the CGC.

I. INTRODUCTION
Photon-hadron production in high energy pp and pA collisions has been put forward [1,2] as a promising probe of the small-x hadron wavefunction [3][4][5]. The distinguished character of this process lies in the fact that a well isolated photon does not participate in strong interactions and therefore the uncertainties related to hadronization of the final state are reduced compared to the more abundant di-hadron (hh) production. The isolated photon component is defined to be separated by a suitable isolation algorithm from the so-called decay photons, which are radiated from the π 0 → γγ decays. This, in turn, suppresses photons fragmenting off jets in the course of hadronization. Recently, the PHENIX collaboration at RHIC [6,7] and the ALICE collaboration at the LHC [8] reported on the measurements of the cross section of an isolated photon in association with an unidentified charged hadron in pp and pA that is within the kinematic region potentially sensitive to small-x physics. In this work we are motivated to pursue the phenomenological implications of those measurements for the first time.
A long standing hypothesis in high energy QCD is that the rapid growth of gluon radiation inside the hadron wavefunction at small-x is balanced out by a gluon recombination process, which leads to the phenomena of gluon saturation. The theoretical framework behind this basic picture is called the Color Glass Condensate (CGC) [3][4][5] where the energy, or small-x, enhanced logarithms ∼ α S log(s) ∼ α S log (1/x), can be consistently resummed at each order in the strong coupling α S . The resulting gluon radiation in CGC becomes concentrated around a dynamically generated transverse momentum -the saturation scale Q S , that is build up from multiple scattering of partons on the dense target. In the context of the simultaneous production of two particles, the back-to-back transverse momenta kinematics of the conventional leading order 2 → 2 partonic process gets disrupted in the CGC. This becomes more dramatic for a heavy nuclei, where a scaling Q 2 S ∼ A 1/3 is expected. It is with the above picture in mind that γh [9][10][11][12][13] and γ * h [14][15][16] production in pA (and recently also the related γ − jet production in eA [17]) have been explored as possible pathways to a phenomenological validation of the CGC. Saturation effects become prominent when the imbalance momenta of the photon (k γ⊥ ) and the underlying final state parton (q ⊥ ) is of the order of the saturation scale k ⊥ ≡ |k γ⊥ + q ⊥ | ∼ Q S . In terms of the distribution over the azimuthal angle ∆φ ≡ φ γ − φ h , the general expected feature is a broadening of the away side peak, ∆φ = π, as the kinematics condition k ⊥ → Q S becomes satisfied.
While both γh and hh production are understandable with such a simple physics argument of the broadening of the away side peak, further investigation into the two processes revealed that they are sensitive to two different unintegrated gluon distributions (UGD)s [18]. At leading order, γh production is probing the gluon dipole distribution, while hh production is instead probing a combination of the dipole and the Weiszäcker-Williams distribution [19,20]. These two UGDs have a completely different behavior at small k ⊥ (see e. g. [20] and also Eq. (5) in the following section). Thus, according to CGC, the theoretical prediction for γh production at the away side peak may look nothing like the prediction for hh production. Indeed this is another reason to explore γh production in CGC.
Alternatively, the broadening of the away side peak may be completely unrelated to the initial state, but simply be induced by the soft gluon radiation, the so-called Sudakov effect [21][22][23][24][25]. Soft gluon radiation generates double logs as α S log 2 (k ⊥ /Q), where Q is a hard scale in the process, and therefore becomes enhanced when the underlying individual partonic transverse momenta is hard, k ⊥ Q. For the high energy kinematics case when in addition k ⊥ ∼ Q S holds, saturation effects and Sudakov effects become two contributing mechanisms and so clearly both should be taken into account. A recent development [26,27] showed that it is possible to perform a simultaneous resummation of both the Sudakov logs and the small-x logs in a consistent way. On the other hand, since Q 2 S ∼ A 1/3 for a heavy nuclei, it is expected the saturation effect would become more important in the pA collision, while both saturation and Sudakov effect are needed for a proper description in pp collisions. This is indeed confirmed by recent studies in hh production [28,29] but also in dijet [30], Z-boson [31] and Z-jet [32] productions as well as in related processes in eA collisions [33][34][35], which prompts us to take into account both the small-x and Sudakov effects in order to establish a realistic baseline for a proper extraction of the appropriate nuclear modifications based on a comparison of the obtained results in pp and pA collisions.
The main purpose of this work is to perform a numerical computation of the isolated photon-hadron cross section within the CGC framework, while taking into account Sudakov effects for the first time. In Sec. II we explain our approach with the main formulas given in (13) and (14). Our numerical results, shown in Sec. III, concern the azimuthal angle (∆φ) distributions measured at RHIC at 200 GeV and 510 GeV in pp and the LHC at 5.02 TeV in pp and pPb and also the out-of-plane transverse momentum distributions measured at RHIC. Our main finding is that while the CGC results alone produce too narrow distributions compared to the data a reasonable agreement with the data is possible when the Sudakov effect is taken into account. Additionally we provide predictions for the magnitude of the nuclear effect in pA vs pp also arguing in favour of more symmetric transverse momentum kinematics of the γh system where the effect of gluon saturation would be better resolved. In the final Section IV we summarize our findings while also providing an outlook to further theoretical investigations.

II. THEORETICAL FRAMEWORK
In this section we give a quick recap of the main formulas that will be used in obtaining our numerical results for p(P p )A(P A ) → h ± (P h )γ(k γ )X production. Throughout this work, we will use p for the projectile proton, A for the target nucleus (or proton in pp collisions). We will also use the shorthand notation h ± ≡ h + +h − to express observables which are obtained by adding the contributions of positive and negative unidentified charged hadrons. We will be using the leading order formulas for γh production -at mid-rapidity the dominant channel is q(p)g(k) → q(q)γ(k γ ), which we firstly address below in the CGC framework, and subsequently incorporate soft gluon resummation in the following subsection. The final formula that we will be using in our numerical computations is given in (13) and (14).
A. Photon-hadron cross section in the CGC framework The following CGC formulas are computed using a dilute-dense framework [3][4][5] where a dilute projectile parton passes through a dense (nuclear) target described by a classical gluon field. In the dilute-dense framework an all-order scattering on the target is taken into account building up a finite transverse momentum in the final state, while the dilute projectile is treated order-by-order in perturbation theory. The leading order pA → hγX inclusive CGC cross section can be straightforwardly obtained from the underlying partonic qg → qγ [36][37][38] channel. In the massless quark limit it simplifies to the following expression, where f q (x p , µ 2 ) is the collinear quark distribution function, D q (z h , µ 2 ) is the relevant collinear fragmentation function of a quark with flavor q and momenta q µ = P µ h /z h to a particular hadron species at a factorization scale µ 2 . In this work we are using the CTEQ6M quark distributions [39] and the DSS fragmentation functions [40].
We have k ⊥ ≡ k γ⊥ + q ⊥ as the imbalance momentum. Furthermore, (πR 2 A ) is the target area. The hard factorσ is given asσ Here, P qγ (z) is the quark-to-photon splitting function with z = k + γ /(k + γ + q + ). The remaining kinematic variables are given as with the projectile and target light-cone momenta in the center of mass frame given as P − A = P + p = s 2 , where s is the center-of-mass collision energy. Here, and in the following, the light-cone variables are defined as p ± = (p 0 ± p 3 )/ √ 2 with the rapidity η = log(p + /p − )/2. η γ and η h are the photon and the hadron rapidities, respectively.
The functionÑ A,Y A (k ⊥ ) is the CGC dipole in the fundamental representatioñ whereŨ (b ⊥ ) is the fundamental light-like Wilson line arising from all order scattering on the dense target. The (1), signifies a broadening of the collinear 2 → 2 away side peak, that would be represented simply by a δ (2) (k ⊥ ) (see (A1)) and that is also shifted from k ⊥ = 0 to k ⊥ ∼ Q S .
As mentioned in the introduction, a well known theoretical distinction of γh with respect to hh correlations is that they probe fundamentally different gluon distributions [19,20]. The leading order γh cross section is proportional to the gluon dipole distribution ϕ DP (Y, k ⊥ ), while the leading order hh cross section is proportional to a combination of the Weiszäcker-Williams gluon distribution ϕ WW (Y, k ⊥ ) and the dipole distribution [19,20]. While these two distributions both display a high-k ⊥ perturbative 1/k 2 ⊥ tail, their behavior is completely different for low k ⊥ where we have This prediction is particularly important at the away-side-peak, ∆φ = π, where small values of k ⊥ would be probed.
In the case of γh correlations, ϕ DP (Y, k ⊥ ) causes a dip in the cross section, see [10] and also [11,[13][14][15][16], since the underlying partonic cross section in (1) is strictly vanishing as the kinematics condition k ⊥ = |k γ⊥ + q ⊥ | becomes satisfied. On the other hand, due to the low k ⊥ behavior of ϕ WW (Y, k ⊥ ) the hh production would rather level off to a constant at the away side peak. In this paper, the Y A dependence of the gluon dipole follows the running coupling Balitsky-Kovchegov evolution equation [41][42][43] (rcBK) which is a good approximation to the more general Jalilian-Marian-Iancu-McLerran-Weigert-Leonidov-Kovner (JIMWLK) evolution [44][45][46][47] of the dipoles. As a general feature of small-x evolution the distribution broadens as x becomes smaller. The initial condition for the rcBK evolution is set at x = x 0 = 0.01 where the dipole is given by the anomalous dimension McLerran-Venugopalan (MV) model, which we will call here the MV γ model. The initial condition for the rcBK evolution of the fundamental dipole is explicitly given bỹ where Y 0 = log 1/x 0 , Q S0 is the initial saturation momentum, γ is the anomalous dimension and Λ IR is the IR cutoff of the model. We use the parameter set [48] γ = 1.119, (Q p S,0 ) 2 = 0.169 GeV 2 , Λ IR = 0.241 GeV. The rcBK equation provides the dipole distribution at x < x 0 while for x > x 0 we use the matching to the collinear gluon PDF as explained in [49]. The matching procedure fixes the proton radius to R p = 0.5257 fm. For the nuclei we use Q 2 S0,A = cA 1/3 Q 2 S0,p where c 0.5 [50]. In our computations we set Q 2 S0,A = 3Q 2 S0,p for heavy nuclei such as Pb and Au.
We will be computing the cross section for an isolated γh production, defined in a standard way through a fixed isolation cut in η − φ space around the photon, defined by (η γ , φ γ ). An isolation cone, R, is introduced to cut out a region where ∆η 2 + ∆φ 2 < R in order to isolate the photon from any soft and collinear hadronic activity within R. Here ∆η = η γ − η h and ∆φ = φ γ − φ h . For the cross section at hand (1) this effectively suppresses its fragmentation component that would appear when γ is being emitted collinearly to the underlying final state parton [9].

B. Implementing the Sudakov effect
The implementation of the Sudakov effect rests on the Collins-Soper-Sterman (CSS) or transverse momentum resummation formalism [24,25]. The inclusion of soft gluon radiation to all orders introduces a transverse momentum dependence into the distribution and fragmentation functions that follows the CSS evolution equation [24,25]. The result of CSS evolution is obtained in b ⊥ -space resulting in a compact formula for the so-called Sudakov factor. Refs. [26,27] further demonstrated that Sudakov resummation can be accomplished on top of the small-x resummation at the one-loop order. In the context of γh production considered here the relevant part of the cross section (1) gets modified through where we multiplied by a prefactor k 2 ⊥ to get the unintegrated gluon distribution 1 . In Eq. (7) the exponential e −S Sud (b ⊥ ,Q) corresponds to the resummed contribution from soft-gluon radiations, where the Sudakov factor S Sud (b ⊥ , Q), has been introduced as the radiation kernel summed over µ b up to the hard scale Q. The Sudakov factor is given by the following generic form [51,52] Here µ b is a factorization scale according to the b * -prescription [24] where b max = 1.5 GeV −1 . With this prescription µ b > 2e −γ E /b max thus preventing the b ⊥ -integral from entering the non-perturbative region. For the hard scale we take Q 2 = x p x A s. The quantities A and B are channel dependent coefficients that can be computed in perturbation theory. For each initial and final state quark (or gluon) in the channel we have [51,52]. According to Refs. [26,27], the pertubative prefactor is absent in the single-log B term for the case of any incoming small-x gluon. Consequently, we do not count the small-x gluon in obtaining this coefficient. For the qg → qγ channel we should use In order to compensate for the missing effect at large b ⊥ it is common to add a non-perturbative Sudakov factor In this work we are using the parametrization [53] where g 1 = 0.212 GeV 2 , g 2 = 0.84, Q 2 0 = 2.4 GeV 2 . Eq. (12) should in principle be used for each quark and gluon in the initial or the final state. The small-x gluon already contains some non-perturbative information through the k ⊥ -dependent distribution. To avoid any possible double counting, the small-x gluon is therefore dropped [28,33] and we have our final formula for the γh cross section reads where we recall that k ⊥ = k γ⊥ + P h⊥ /z h . The explicit computation of Eq. (14) is performed in momentum space through a convolution of the momentum space gluon distribution, k 2 ⊥ÑA,YA (k ⊥ ), with the Fourier transform of W (z h , x p , b ⊥ , Q). The latter is computed numerically with the algorithm in [54]. We also remark that with resumming infinite number of soft gluons exact kinematics relations (3) are lost, and in practice [55,56] one resorts to an approximate relation x p,A = p ⊥ (e ±ηγ + e ±η h )/ √ s, with p ⊥ ≡ max(k γ⊥ , q ⊥ ). According to (13) the Sudakov factor brings an additional dependence of the imbalance k ⊥ -distributions on the hard scale Q. As a typical result of CSS evolution, the k ⊥ distribution gets broadened as a function of the hard scale Q. This can be intuitively understood from the leading log 2 (k 2 ⊥ /Q 2 )-dependence as the increase in the hard scale brings more opportunity for soft gluon radiation. In the CGC framework radiative processes are enhanced as log(1/x). Therefore, since x ∼ Q/ √ s, an increase in the hard scale would lead to a narrower distribution. It is thus theoretically and phenomenologically interesting to check the interplay between the small-x resummation and the Sudakov resummation.

III. NUMERICAL RESULTS AND DISCUSSION
In this Section we perform a numerical computation of the cross section for the isolated γh cross section based on (14). In what follows, we abbreviate this result as as CGC w(ith) Sud(akov). To isolate the effects of both small-x and Sudakov gluons and to gain further insight we will also show results from the CGC framework alone (Eq. (1)), and from the cross-section in the collinear pQCD framework with only the Sudakov effect taken on account. For the latter we refer to (A3) in Appendix A. We abbreviate these results as CGC w/o Sud and Sud, respectively.
To keep the discussion within phenomenological reach, we will mostly study systems and kinematic windows already probed by different experimental collaborations. Firstly, the PHENIX experiment measured isolated γh ± production in pp collisions at √ s = 200 GeV [7], covering the kinematic range where an isolation cut R = 0.3 has been applied to the photon. For √ s = 510 GeV collisions [6] we have the following kinematics |η γ | < 0.35 , 7 GeV < k γ⊥ < 15 GeV , |η h | < 0.35 , 0.7 GeV < P h⊥ < 10 GeV , and R = 0.4. The kinematics (15) and (16) has been further separated in k γ⊥ × P h⊥ bins as indicated in Figs. 1 and 2, respectively. The ALICE experiment measured isolated γh ± production in pp and pPb collisions at √ s = 5.02 TeV [8] with the Isolated γh ± angular distribution for mid-rapidity pp collisions at 510 GeV given in k γ⊥ × P h⊥ bins vs data from PHENIX [6]. The experimental data has been shifted for the underlying event.
following kinematics and R = 0.4. The transverse hadron momenta P h⊥ is further distributed into bins as indicated in Fig. 3. Before considering the results we make a general remark about our computation based on (14). In the kinematics region where the transverse momentum of the final state k γ⊥ + P h⊥ is large, and where also (k γ⊥ + P h⊥ ) 2 Q 2 , the perturbative Sudakov factor, (8), would dominate the overall two-particle momentum imbalance in the cross section. While the non-perturbative Sudakov factor is necessary to carry out the b ⊥ -integral, it is irrelevant for the k ⊥ -spectrum, as demonstrated in [55]. However, in our computations this is not completely the case, since for RHIC kinematics we have Q ∼ 7 − 21 GeV, where the non-perturbative Sudakov factor should also play a role. At the LHC Q ∼ 17 − 55 GeV, but also k ⊥ is possibly larger due to the more asymmetric γh momenta configuration.
We will first show our results for the angular distributions and compare them with the data from RHIC and the LHC. In our compuations we focus only on the shape of the γh yield and thus normalize both our theoretical curves and the experimental data to unity. Similar procedure has been employed e. g. in [28,56,57]. We compute as well the normalized out-of-plane transverse momentum distributions and extract their Gaussian widths and compare both with the data from RHIC. We also provide the predicted nuclear modifications by showing the CGC w Sud computation at RHIC and the LHC.

A. Angular distributions
In Fig. 1 we compare our CGC w Sud result with the 200 GeV pp data at RHIC [7] and find a fair agreement in most of the bins. Here the experimental data have been shifted vertically for the underlying event, that is, by the background of uncorrelated γh ± pairs [7]. The CGC w Sud results are also compared with a CGC w/o Sud computation. The away side peak from the CGC w/o Sud results is clearly too narrow to be able to describe the data. As argued previously, the dip at ∆φ = π, and the resulting double peak structure around it, for the CGC w/o Sud computation is due to the low k ⊥ behavior of the dipole gluon distribution ϕ DP (Y, k ⊥ ) ∼ k 2 ⊥ /Q 2 S . For this kinematics the double peak is strongly focused in a narrow region around ∆φ = π. Let us stress again here that the presence of a double peak in general is a robust prediction [10,11,[13][14][15][16] of the leading order γh production in CGC. From the PHENIX data alone it is difficult to find support for this feature, though it might be simply missed by the  FIG. 3. Isolated γh ± angular distribution for mid-rapidity pp and pPb collisions at 5.02 TeV for 12 GeV < k γ⊥ < 40 GeV for three bins in P h⊥ vs data from ALICE [8]. The magenta points (lines) correspond to pp collisions while the teal squares (lines) correspond to a pPb collision. . Isolated γh ± angular distribution for mid-rapidity pp and pA collisions at 200 GeV. The higher curve (and narrower peak) correspond to the 5 < k γ⊥ < 7 GeV bin, while the lower (wider) peak corresponds to a harder photon momentum bin, 12 < k γ⊥ < 15 GeV. experimental resolution. In any case, our prediction is that, for the kinematics considered here, the double peak is completely washed away by including the Sudakov effect. It is instructive to also compare to a Sud only result based on the leading order qg → qγ collinear formula. While the Sud result gets closer to the data than CGC w/o Sud in general, best results are obtained when both the CGC and the Sudakov effects are taken into account.
In Fig. 2 we compare the CGC w Sud computation with the 510 GeV pp data from RHIC [6] where again we find overall good agreement with the data. In Fig. 3 we compare with the 5.02 TeV pp and pPb data from LHC [8] and find good agreement with both the pp and the pPb data. Due to asymmetric kinematics in Fig. 3 our results show only a small nuclear effect. Below (see Fig. 5) we look at more symmetric configurations.
In Fig. 4 we make predictions for pA collisions at 200 GeV by considering a subset of bins from Fig. 1. We pick up the smallest and the largest k γ⊥ bin from Fig. 1 and distribute the results along P h⊥ bins. For comparison, we also repeat the pp results from Fig. 1. It is useful at this point to discuss the systematics across bins. Lets compare two pp (or pA) curves corresponding to the smallest and the largest k γ⊥ , for a fixed P h⊥ . We see that as k γ⊥ (the transverse momentum of the trigger particle) is increased, the away side peak gets narrower. We can understand this in a intuitive way as follows. At large trigger k γ⊥ also the momentum imbalance k ⊥ eventually increases (P h⊥ is held fixed) and so the Sudakov logs become less prominent. Note that for high k γ⊥ -kinematics the details of the non-perturbative Sudakov factor should be negligible. At the same time, nonlinear effects from the CGC also play less of a role as with high trigger k γ⊥ we are probing the perturbative tail of the gluon transverse momentum distribution.   FIG. 5. Isolated γh ± angular distribution for mid-rapidity pp (solid) and pA (dashed) collisions at 5.02 TeV. Each pp and pA curve pair corresponds to a different binning of photon momentum. In descending order, the pair correspond to bins of 5 < k γ⊥ < 7 GeV, 9 < k γ⊥ < 12 GeV , and 12 < k γ⊥ < 15 GeV. Therefore, we expect the probability for a high-k γ⊥ trigger photon to scatter at an angle ∆φ = π to be strongly suppressed, explaining the narrower shape. Likewise, for smaller k γ⊥ the away side peak will get broadened. Our second point concerns the nuclear effect that is visible by a comparison of the full (pp) and the dashed (pA) curves in Fig. 4. Due to a larger saturation scale in the nuclei than in the proton target we observe a suppression of the away side peak in pA in comparison to pp for all kinematic bins considered.
In Fig. 5 we present projections for mid-rapidity pp and pA collisions at the LHC energy of 5.02 TeV. We keep the associated hadron bins from Fig. 3 and lower the trigger photon k γ⊥ in order to reach a more symmetric configuration. The imbalance momentum of the γh ± system thus approaches the saturation scale at the away side peak leading to a more substantial nuclear effect as can be seen from Fig. 5.

B. Out-of-plane transverse momentum distributions
In addition to the angular distributions, PHENIX also measured [6,7] the so-called out-of-plane p out distributions where p out ≡ P h⊥ sin(∆φ), binned as a function of x E ≡ −k γ⊥ ·P h⊥ /k 2 γ⊥ = −P h⊥ cos(∆φ)/k γ⊥ . A quick computation shows that p 2 out = z 2 h k 2 ⊥ − k 2 γ⊥ (1 − x E /z h ) 2 at the leading order. Close to the away side peak ∆φ = π we have x E z h and so p out z h k ⊥ . By binning the result in x E , we can appreciate that the p out -distributions serve as a proxy for the intrinsic k ⊥ distributions [7,58].
In Fig. 6 we compare our results with the PHENIX 200 GeV data [7] in x E bins. The CGC w Sud results show a good agreement with the data in the small to moderate p out region (up to ∼ 1−2 GeV), and also in the large p out region for the x E bin 0.1 < x E < 0.15. For the remaining two bins, the small to moderate p out region (p out ∼ 1 − 2 GeV) is also nicely described by our result, while as p out increases our results tend to overshoot the data. For comparison, the CGC w/o Sud computation is also shown where similar conclusions hold as for the angular distributions: the CGC w/o Sud computation produces a too narrow distribution which cannot be accommodated within the experimental data. In addition, the CGC w/o Sud computation predicts a double peak in the distribution concentrated in a narrow region around p out = 0. Again, the present PHENIX data sets alone do not allow us to judge clearly whether this feature is supported or not. We can only underline the importance of a more complete computation, which includes the Sudakov resummation, where a double peak structure is not present for these kinematics. For completeness we also show the results of a Sud computation which seems to do a good job in the moderate p out region, but tends to undershoot the data for large p out .
To have a closer look at the small p out region, PHENIX extracted the Gaussian widths of the p out -distributions across 0.1 < x E < 0.5 from Fig. 6 by assuming a Gaussian-like shape in the range −1.1 GeV < p out < 1.1 GeV. The systematic error is estimated by varying this range by ±0.2 GeV [7]. Using the same procedure we compute the Gaussian widths of the p out distributions from Fig. 6. Our results are shown in Fig. 7 as a function of the hard scale k γ⊥ in comparison to the PHENIX 200 GeV and 510 GeV pp data. The theoretical bands correspond to varying the p out Gaussian fit by ±0.2 GeV as in [7]. We see that the best description of both the 200 GeV and the 510 GeV data is obtained with the CGC w Sud computation. In Fig. 7 we also plot Sud results for the Gaussian widths and find that they are below the PHENIX data. The CGC w/o Sud widths are not shown -due to the double peak structure the behavior is clearly not Gaussian-like for small p out . But even if we choose to ignore this issue, based on the fact that the double peak is rather narrow for this particular kinematics, it is visible already by the naked eye from Fig. 6 that the distribution are too narrow in comparison to the data.
In Fig. 8, as a measure of the nuclear effect, we show our predictions for the difference in the Gaussian widths squared computed in mid-rapidity pA and in pp at RHIC (left) and at the LHC (right) kinematics. For convenience we denote this quantity as p 2 out pA − p 2 out pp , though the meaning of the "averaging" procedure is in the sense of fitting the small p out region with a Gaussian -same as we used for Fig. 7. Our choice of kinematics matches the one used in the preliminary RHIC pPb vs. pp results in [58][59][60] with the hard scale 5 GeV < k γ⊥ < 9 GeV and we plot our results as a function of x E . We are tempted to compare the 200 GeV results from Fig. 8 with the preliminary result in Refs. [58][59][60] that shows up to about 0.1 GeV 2 broader widths in pA that in pp, for x E < 0.4, albeit for π 0 h ± production. In our γh ± computation we find 0.15 GeV 2 for x E < 0.4. Increasing the energy to 5.02 TeV we find that the difference in the pA vs. pp Gaussian widths grows, becoming more pronounced at large x E .

IV. CONCLUSIONS
In this paper we have numerically computed the isolated photon-hadron production cross section based on the leading order CGC formula with the Sudakov resummation explicitly taken into account for the first time. We have demonstrated that both the CGC and the Sudakov effects are important to obtain a reasonable description of the data at RHIC and the LHC. We have provided predictions for the nuclear suppression around the away side peak ∆φ = π and the nuclear broadening of the intrinsic transverse momentum distributions at RHIC and at the LHC. While the PHENIX collaboration published only pp data, the ALICE collaboration has results also for pp and pPb.
Unfortunately, there is no clear evidence of a nuclear effect from the data alone. This might be due to asymmetric k γ⊥ × P h⊥ binning by ALICE and we argue that in more symmetric configurations the nuclear effect would be more apparent. In addition to γh ± production the PHENIX experiment measured also π 0 h ± production [6,7]. The data seems to indicate a suppression in the γh ± vs. π 0 h ± production in the bins that are close to the away side peak. It would be interesting to make a detailed side-by-side comparison of γh ± vs. π 0 h ± production in CGC.
As a future work we plan to take into account next-to-leading order corrections which bring 2 → 3 partonic processes into play. Firstly, already in the collinear framework, the presence of an additional (unobserved) parton present in the final state naturally disrupts the γh back-to-back kinematics. As the Sudakov effect is associated with soft gluon branching it is most important in the region close to the away side peak with the momentum imbalance k ⊥ = |q ⊥ + k γ⊥ | such that k ⊥ Q. By contrast, the 2 → 3 processes are genuine hard branchings, that are able to support large momentum imbalances k ⊥ ∼ Q. Therefore, next-to-leading order corrections and Sudakov effects are important to get a more complete phase-space picture of γh correlations. This has recently been applied in describing hard-p ⊥ γ-jet data in pp, see e. g. [61,62] and also [57]. We also mention here recent works on transverse momentum resummation in γ-jet production that also takes into account a preferred direction set by the final state jet [63,64].
In the small-x kinematics region with k ⊥ ∼ Q S we are interested, a complete next-to-leading order treatment would include considering the gg → qqγ [65,66] and qg → qgγ channels [67,68] (see also [69,70] for related higher order computations in eA collisions). With the current photon detectors in the mid-rapidity region for both RHIC and LHC, as a first step, the inclusion of the gg → qqγ channel might be enough. Considering planned forward upgrades [71,72], where the isolated photon signal would be more favourably extracted from the π 0 background, taking into account also the qg → qgγ channel becomes important.
where for the qg → qγ channel [73] we have 2σ with the conventional Mandelstam variablesŝ = (p + k) 2 ,t = (k − k γ ) 2 andû = (p − k γ ) 2 . The parton momenta fractions are fixed as x p = (k + γ + q ) /P + p and x A = (k − γ + q − )/P − A . Taking into account the Sudakov resummation yields the following result The Sudakov factor S Sud (b ⊥ , Q) + S non−pert (b ⊥ , Q) contains perturbative and non-perturbative pieces as in (13). The perturbative piece has the same form as in (8), however, because of the initial stage gluon, in the collinear limit one has B = 2B q + B g (the double-log A-coefficient remains the same as in (10)). In the non-perturbative piece we have now S non−pert (b ⊥ , Q) = 2S q non−pert (b ⊥ , Q) + S g non−pert (b ⊥ , Q). It is useful to explicitly confirm that we can recover (A1) from (1) in the leading twist limit. First lets rewrite (A2) in a more convenient form. Usingû = −2p · k γ = −k 2 γ⊥ /z andŝ = 2q · k γ = k 2 γ⊥ /z(1 − z) we get Recall now that the transverse momentum dependent gluon distribution ϕ(Y, k ⊥ ) function is related to the adjoint dipole N Y (k ⊥ ) as [74] ϕ(Y, k ⊥ ) = (πR 2 A ) Integrating ϕ(Y, k ⊥ ) over k ⊥ returns the gluon distribution xf g (x) = 1 π 2 which may be formally inverted as ϕ(Y, k ⊥ ) = xf g (x)π 2 (2π) 2 δ (2) (k ⊥ ). In the large N c limit, the adjoint dipole and the fundamental dipole in coordinate space are related as N Y (b ⊥ ) =Ñ , see e. g. [75]. In the leading twist approximation this effectively becomes N Y (b ⊥ ) 2Ñ Y (b ⊥ ) and so we can write (in momentum space) Then, Eq. (1) yields Eq. (A1).