Scaling properties of direct photon yields in heavy ion collisions

A recent analysis from the PHENIX collaboration of available direct photon measurement results in collisions of various systems such as Au+Au, Cu+Cu, and Pb+Pb, at different beam energies ranging from 39 to 2760 GeV, has shown a universal, within experimental uncertainties, $multiplicity$ scaling, in which direct photon $p_{T}$-spectra for transverse momenta up to 2 GeV/$c$ are scaled with charged hadron pseudorapidity density at midrapidity raised to power $\alpha=1.25$. On the other hand, those direct photon $p_{T}$-spectra also exhibit $geometrical$ scaling in the similar $p_{T}$ range. Assuming power-law dependence of the scaled photon spectra for both scaling laws, we formulate two independent conditions for the power $\alpha$, which overshoot experimental data by $\sim 10\%$ on average. We discuss possible sources that might improve this estimate.


Introduction
Measurements of direct photons provide unique opportunities in probing and studying the properties and evolution of the matter produced in heavy ion collisions (HIC). These photons are defined to be produced from all the sources except for hadronic decays. Since they hardly interact with the "fireball" of quarks and gluons due to a small interaction cross section with the medium, the information they carry from the time of their production is not washed out by final state interactions. Experimentally measured low momentum direct photon p T -spectra by PHENIX (in Au+Au at  [1,2,3,4,5,6,7] and ALICE (in Pb+Pb at √ s N N = 2760 GeV) [8] collaborations in HIC are enhanced with respect to N coll (number of binary nucleon collisions) scaled reference yield (measured or calculated) in p+p collisions. Low momentum direct photon measurements by STAR collaboration show less enhancement [9]. Earlier low energy WA98 data [10] have mostly upper bounds in the relevant p T ≤ 2 GeV/c region. Direct photons in HIC also show large anisotropy (elliptic flow) [7,11,12].     Fig. 1 Direct photon invariant yield as a function of p T for various colliding nuclei and collision centrality selections at three center-of-mass energies. The plot labels give the references of the shown data. The vertical lines of the data points describe statistical uncertainties, the boxes describe systematic uncertainties. Fig. 1 shows several data sets of direct photon p Tspectra at low-(< 1 GeV/c) and intermediate-(from 1 GeV/c up to ∼ 5 GeV/c) p T regions. There have been many theoretical attempts to reproduce the photon yields shown in Fig. 1 and flow coefficients (that are not discussed here) with, however, mixed success. Hydrodynamical simulations of the fireball evolution [13,14,15,16], calculations in the framework of the ellipticfireball expansion scenario [17,18,19], Parton-Hadron-String Dynamics transport approach [20,21,22,23], as well as the spectral function approach [24,25,26,27], have encountered difficulties in the simultaneous description of the observed large yields and large anisotropies, which have given rise to a challenging problem, commonly referred to as direct (thermal) photon flow puzzle [28] (for recent review see [29]).
The early attempts to describe the photon yields by initial state models can be found in Refs. [30,31,32,33,34]. More recently initial state models have been used both as a part of hydrodynamic evolution [16] and in a bottom-up thermalization scenario [35,36]. In the latter case [36] good fits for the photon yields have been obtained for PHENIX and ALICE data (the anisotropies are not discussed there).
Much work still needs to be done to include the anisotropies in the initial state model calculations, in order to directly address the aforementioned puzzle. Ref. [33] discusses the possibility of late stage elliptic flow effects to be explained in the combined Glasma and thermal photon production scenarios, by which the photon mean emission time is shifted toward a later time scale as compared to the mean emission time obtained only from the thermal ansatz. Nonetheless, as it has been pointed out in [36], a more realistic description of the photon production might come out within the combined framework of the initial and late stage models. In this regard a promising step is already undertaken in [16], where the event-by-event hydrodynamical model uses IP-Glasma initial state conditions. However, the most important argument from [30,31,36], directly applicable to the studies of this paper, indicates that the initial state models potentially allow the photon p T -spectra to exhibit geometrical scaling (GS) 1 . In fact, two types of scalings are observed in direct photon p T distributions in HIC and hadron collisions: Direct photon multiplicity scaling (MS). It is an experimental observation in HIC, which shows that direct photon invariant yields follow one universal curve within experimental statistical and systematic uncertainties for various colliding species and system sizes at 1 We will perform detailed studies of GS of photon p T -spectra in HIC elsewhere [37]. different center-of-mass energies, when scaled by charged hadron multiplicities at midrapidity raised to the power of α = 1.25 [5]: where G is a universal energy-and multiplicity-independent function of p T (see Fig. 2) 2 , and Q 0 ∼ 1 GeV/c. MS holds for small and intermediate p T up to ∼ 2 ÷ 2.5 GeV/c, precisely in the enhancement region discussed above.
[GeV/c]          Geometrical scaling (GS). Here the invariant yields of direct photons (and also of charged hadrons) can be related to a universal energy-independent function of scaling variable τ [31,37,38,39,40,41]: where S T is a parameter characterizing geometrical overlap area of colliding nuclei. The variable τ is given by GS has been for the first time observed in deep inelastic scattering (DIS) in e+p collisions [42]. The quantity Q s (x) in (3) is the saturation scale, i.e., the transverse momentum at which the number of gluons with longitudinal momentum fraction x, called Bjorken variable, cannot grow any more due to the nonlinearity of QCD evolution. This happens at gluon occupancy numbers being of the order of 1/α s (α s -strong coupling constant). The phenomenon of gluon saturation is generated from such dynamics [43,44,45,46,47,48,49,50] and can be described by the Color Glass Condensate (CGC) effective theory of high energy scattering [51, 52,52, ?,55].
One still should remember that GS is more general. Indeed, it has been shown in Ref. [56] that the approximate GS is in fact a general property of solutions of the QCD collinear evolution with general boundary conditions. Furthermore, the results of Ref. [57] indicate that the rapidity evolution preserves GS beyond saturation when starting from the boundary conditions in the saturation region. It means that the scaling properties of an initial state are preserved and even built up by the QCD evolution. This may explain the fact that the geometrical scaling is observed in the hadronic p T -spectra and even more so in the case of photons.
Saturation models predict specific dependence of Q s on the collision energy, transverse momentum and collision centrality class, quantified by the number of participating nucleons with δ ∼ 2/3 and λ ∼ 0.2 ÷ 0.35 [41]. Here we have chosen x 0 = 10 −3 for a typical value of Bjorken x where the saturation effects become important. The direct photon GS is illustrated in Fig. 3.
In the present paper we formulate two conditions that allow one to calculate the power α of MS in terms of the parameters that enter the theoretical parametrization of GS, assuming power-law dependence of the scaled photon spectra. These two conditions result in two predictions for α: 1.34 and 1.40, i.e. ∼ 10% above the experimental value. Given the simplicity of the present analysis and approximate nature of both scaling laws, this result may be considered as a satisfactory one. Nevertheless, we also discuss possible corrections that might improve this theoretical prediction of α.

Basics of Geometrical Scaling
Geometrical Scaling is based upon the assumption that bulk properties of charged particles and photons share scaling properties encoded in the initial state wave function of colliding hadrons/nuclei. While detailed models of direct photon production in the early and intermediate stages of HIC can be found in Refs. [30, 32, 33, 35,  36] as it was noted in the introduction, here we discuss general assumptions and necessary conditions leading to GS. Differential gluon production cross-section can be written in terms of the unintegrated gluon distributions [58]: where C contains color factors and numerical constants.
The Bjorken x's of colliding partons read In the midrapidity region y ≈ 0, hence There exist many models of unintegrated gluon distributions ϕ(x, k 2 T ) (see e.g. [59,60,61,62]) that enter (5). Most of them share two important features: geometrical scaling and dependence on the transverse area parameter A T : where φ is a dimensionless function of the scaling variable τ 2 = k 2 T /Q 2 s (x). The precise meaning of A T is best understood in a picture where the impact parameter is also taken into account [63,64].
Ignoring momentum dependence of the strong coupling constant, on dimensional analysis grounds, we arrive at T A where F is a universal, energy-independent function of the scaling variable τ in (3). If one assumes that A T is an energy-independent constant, which is true in the case of the GBW model [59,60] for DIS, then it is the differential cross-section that should exhibit GS. Indeed, it has been shown in Ref. [65] that for charged particles the differential crosssection in p+p scattering scales better (over larger p T interval) than multiplicity, and the exponent λ ≈ 0.3 is compatible with the DIS scaling [66]. In contrast, for multiplicity scaling λ ≈ 0.2.
In order to obtain the multiplicity distribution, we should divide both sides of (8) by an appropriate inelastic cross-section: This distribution would scale if were energy-independent. This assumption leads to (2). It is the energy dependence of S T that is responsible for different scaling properties of the cross-section and multiplicity of charged hadrons in p+p collisions [66].
In the case of HIC for a fixed centrality, S T has geometrical interpretation as an overlap area of two colliding nuclei [61,68] that scales like N 2/3 part . Throughout this paper we shall assume Q s ∼ S T ∼ N δ part . Possible weak energy dependence of S T would lead to the violation of geometrical scaling. For small systems (like e.g. dA or pA) one has two saturation scales that differently scale with N part and the effective saturation scale has to be taken as Q sat [large]Q sat [small] [31,67]. We address this issue in [37]. Here we concentrate on large systems only.
The functional form of the saturation scale Q s (x) in (4) follows from the nonlinear QCD evolution [49], and in the case of heavy ions, also from the collision geometry [61,68], resulting in the following form of the scaling variable τ : where W = √ s N N × 10 −3 is an effective energy scale shown in (4). Should there be any energy dependence of S T (which enters the definition of Q s ), it is already included in λ, as we fix λ from the fits to data that do not distinguish different sources (QCD nonlinear evolution, S T energy dependence) of the energy increase of Q s (W ). These fits [37], as well as a direct inspection of Fig. 3, indicate that for direct photons λ ≈ 0.2, rather than 0.3. The same conclusion has been reached in an earlier study of Ref. [31].

Charged particle multiplicity
Let us first calculate the charged particle multiplicity for large collision systems from the scaling formula (2). To this end we will need a Jacobian to change the integration variable from p T to τ : which yields where the constant κ includes the integral of F ch , S T and other irrelevant constants. Note that the energy dependence ∼ W 2λ/(2+λ) is compatible with multiplicity growth in p+p as measured at the LHC [69]. For small λ and δ = 2/3 we have the quantity N 3+λ 2+λ δ part ≈ N part , as has been assumed e.g. in Ref [68], though it has been now established that better scaling of charged multiplicities is obtained when N part is replaced by the number of participating quarks (see below).
In order to assess the quality of the scaling formula (13), we plot the scaled charged hadron multiplicity in Fig. 4: for the data from Refs. [70,71,72,73] as a function of N part for δ = 2/3 and for two different choices of λ. We see that, in contrast to the photon GS-scaled yields, the best alignment of all data is achieved for λ ∼ 0.35.
Here Q 0 ∼ 1 GeV/c as before. Some systematic shift is seen in the case of Cu+Cu data. We would expect that the scaled multiplicity in (14) should not depend on N part , while we see that it rises with N part . This is a well-known problem in HIC, and it has been argued that better scaling is obtained with the number of participating quarks [70]. Here it is enough to slightly change the value of δ from 2/3 to 3/4 to ensure N part independence of (14). However, as we shall see in the following, the value of δ does not enter the conditions for the MS exponent α.
The fact that λ for charged particle spectra is different than the one for photons is -to our best knowledge -not understood theoretically. It may be due to the interactions of charged particles in medium. The other    reason might be that for charged particles we extract λ from integrated yields that include the scaling part up to 2 ÷ 3 GeV/c [41] and the non-scaling tail, while in obtaining Eq. (13) we have assumed that GS is present in the whole experimentally measured range. In the following we will accept this difference as a phenomenological observation, and distinguish λ for photons from λ ch for charged particles.

Relating the scaling laws
Let us rewrite the scaling laws (2) and (1) in the following way which can be satisfied only when the left hand side and the right hand side of Eq. (15) have the same p T , W and N part dependence. It is clear that one cannot proceed further without explicit knowledge of the functions F γ and G. To this end we shall assume a power-law dependence Indeed, as shown in Figs. 2 and 3, where we plot F γ and G for Au+Au [1,2,3], Cu+Cu [4] and Pb+Pb [8] direct photon data at different energies and centralities for λ = 0.2 and δ = 2/3, the power-law fall-off [31] with n, m ∼ 4, 4.5 works pretty well at small and intermediate τ (or p T ). Note that m = 4 is a generic prediction for radiation from the CGC [67] at intermediate p T , while -as we have checked by looking at prompt photons up to transverse momenta ∼ 2 GeV/c -fits to the data prefer m = 5. Altogether where we have used S T ∼ N δ part . The functions in (17) are proportional to each other if m = 2 + λ 2 n, 4 + n 4 δ = 3 + λ ch 2 + λ ch αδ, λ 2 n = 2λ ch 2 + λ ch α.
From the first equation we see that the power-like falloff should be faster for MS than for GS, and this prediction works quite well, as can be seen from Figs. 2 and 3. Thereby, from (18) we obtain the following: Both estimates give the power α compatible with the experimental value of 1.25 (given our crude assumptions and approximations). Note that the first estimate follows from the assumptions concerning the geometry of A+A collisions (although it does not depend on the actual value of δ), whereas the second one follows from the energy dependence.

Discussion and conclusions
Let us first observe that the equality of two estimates for α, (19) and (20), within an acceptable range of n, λ and λ ch , cannot be obtained if λ = λ ch . As mentioned in Sec. 3 there is no theoretical understanding of this difference. One possibility might be the apparent violation of GS that has been mentioned in Sec. 2, namely the energy dependence of S T . We do not have any experimental handle on this dependence, however recent Glauber model calculations [74] show slight increase of S T with energy, depicted in Fig. 5, which can be effectively parametrized as a power-law: with λ ∼ 0.06. While (21) may explain a difference in quality of GS for multiplicity vs. cross-section, it is difficult to disentangle this dependence from the genuine energy dependence of Q s on x. This is because we effectively include the energy dependence of S T in the process of tuning λ in search for GS. The naive inclusion of (21) in Eqs. (15) and (13) would only modify the second equation for α in (20) shifting it to 1.29. However, a proper way of including the energy dependence of S T would require a global fit with λ, λ and δ as free parameters, rather than adding λ in (20). There are many other sources that lead to the violation of GS in HIC, which we have ignored. For example, the factorization of the unintegrated gluon densities in (7) into A T and φ, neglecting the momentum depen-dence of α s in (5), as well as simplifying the assumptions about N part dependence and the above-mentioned energy increase of S T . It is obvious from Fig. 5 that these types of dependence do not factorize, moreover, much better scaling is obtained with the power δ = 3/4. The latter, however, is of no importance for our analysis, as δ drops out from the second equation of (18). Finally, simple power-law dependence of the photon spectra is also a simplification that affects the accuracy of our analysis. A choice of more accurate data-driven functions could give an improved value of α, closer to that extracted from experiment. This would require, however, a numerical analysis that we wish to avoid for the clarity of argumentation. One should not forget that both scaling laws observed experimentally are only approximate. This concerns also the charged particle spectra used here do calculate dN ch /dη in (13).
Even though one might argue that ∼ 10% accuracy is quite satisfactory, given a simplicity of the present analysis, our results may indicate that there might be other components of direct photons that exhibit GS but have different functional dependences than the initial stage component. This might resolve the difference between λ and λ ch . Identification of such components would be of great interest.
Our conclusion is that we have linked the multiplicity scaling and geometrical scaling of direct photon p T -spectra, estimating the scaling power α from simple assumptions based on the functional forms of the scaling functions G and F γ , as well as on the energy dependence and N part dependence of the saturation scale. We have obtained two independent, albeit consistent, estimates of α that overshoot experimental value by ∼ 10%. This result relies on the fact that λ = λ ch . We have argued that possible explanation of this difference, apart from trivial facts that charged particles interact with the medium while photons do not and that integrated yields of charged particles include non-scaling tails, might also indicate that other contributions of direct photon production exhibit GS but possibly have different functional forms than the one of the early stage production.

Acknowledgments
We are grateful to Axel Drees and Larry McLerran for many important and fruitful discussions. The research of VK is supported under Department of Energy Contract No. DE-FG02-96-ER40988. The research of MP is supported by the NAWA (Polish National Agency for Academic Exchange) Bekker program and partially by the COST Action CA15213 THOR. MP wishes to