Energy versus centrality dependence of the jet quenching parameter $\hat q$ at RHIC and LHC: a new puzzle?

The central goal of jet quenching studies in high-energy nuclear collisions is the characterization of those QCD medium properties that are accessible by these probes. Most of the discussion in the last years has been focused on the determination of the jet quenching parameter, $\hat q$. We present here an extraction of this parameter using data of inclusive particle suppression at RHIC and LHC energies for different centralities. Our approach consists of fitting a $K$ factor that quantifies the departure of this parameter from an ideal estimate, $K\equiv \hat q/(2\epsilon^{3/4})$, where $\hat q$ is determined by the local medium quantities as provided by hydrodynamical calculations. We find that this $K$ factor is larger at RHIC than at the LHC, as obtained already in previous analyses, but, surprisingly, it is almost independent of the centrality of the collision. Taken at face value, the $K$ factor would not depend on the local properties of the medium as energy density or temperature, but on global collision quantities such as the center of mass energy. This is a very intriguing, unexpected possibility for which we cannot yet provide a clear interpretation. We also comment on the limitations of the formalism that may affect this conclusion.


I. INTRODUCTION
Jet quenching -the suppression of high-energy particles and jets in nucleus-nucleus collisions relative to the expectation from a superposition of nucleon-nucleon ones -is one of the best available tools to characterize the properties of the medium created in collisions of heavy nuclei at high energies. See the recent reviews [1][2][3][4]. A huge number of experimental data is available and has been phenomenologically studied in the last 15 years, starting with the suppression of inclusive particle production at high transverse momentum at the Relativistic Heavy-Ion Collider (RHIC) and now at the Large Hadron Collider (LHC). At present, a large experimental effort is devoted to the study of reconstructed jets, where more differential analyses are possible for a more precise characterization of medium properties.
For reference to the experimental works, see e.g. the reviews [5,6].
The final goal of the jet quenching studies is to extract medium parameters which characterize the QCD matter formed in high-energy nuclear collisions. In this paper we present an analysis of RHIC and LHC data on the nuclear modification factor R AA for inclusive particle production at high transverse momentum using the formalism of the quenching weights proposed in [45][46][47]. The main result of this paper is an extraction of the value for the jet quenching parameterq in a method which has been well tested, is easy to implement and interplay with different hydrodynamical models for the medium (as done previously by some of the authors of the present paper in Ref. [48], see also [15,49,50]), and provides a good description of the experimental data, as we will show. Despite the limitations of the formalism whose applicability is restricted to describing leading particle production in a jet, we argue that the main conclusions of the analysis are solid and somehow unexpected.
Our approach here is to define a local transport coefficient,q, which is solely determined by the local (in position and time) medium quantities as extracted from hydrodynamical models, in particular the energy density. So, we define the jet quenching parameter aŝ q = K2 3/4 , motivated by the ideal estimateq ideal ∼ 2 3/4 [51], and fit the values of K. We do not impose any particular dependence of this K-factor on energy, centrality, temperature, etc. On the contrary, K is the only free parameter in the fit of R AA for each centrality at RHIC and LHC energies. Our main findings are that this K-factor is ∼ 2 − 3 times larger for RHIC than for the LHC (larger values at RHIC than at the LHC have been found before [52]) and, unexpectedly, this K-factor does not seem to depend on the medium parameters, e.g., the temperature, but instead on the center of mass energy of the collision. Indeed, we find these K-factors to be basically independent of centrality both for RHIC and the LHC.
Were the K-factor determined, say, by temperature, then the most central RHIC collisions should present a value similar to semi-peripheral LHC data. This is not the case. We have performed the study with several, quite different, hydrodynamical profiles and these conclusions do not depend on the profile we use, although the values of K do present some dependences and, interestingly, they dramatically change for some different assumptions for the dynamics at initial times before the starting of hydrodynamical evolution.
We have no clear interpretation of this finding. We comment on the different theoretical limitations of the procedure which could affect this result. A more detailed study of these limitations is not easy with present theoretical tools but we expect that it will be possible in the near future. Taken at face value, if our result is not due to a limitation of the technique, it would indicate that the properties of the QCD media produced at RHIC and LHC are different in what concerns the jet quenching process. Some possibilities could be related with the initial stages of the collision, the presence of quantities related with the total energy of the collision (relative differences between local and global bulk/thermodynamical properties, the presence of a magnetic field, etc.) or others. At this point, however, we do not speculate with these possibilities and present our findings as they are obtained. It would be important to check our conclusion with other jet quenching implementations.
The rest of the paper is organized as follows: we present the formalism of energy loss in Section II; in Section III we describe how the energy loss is interfaced with the profiles taken from hydrodynamical models; in Section IV we briefly comment on the details of the three hydrodynamical models used in our analysis; in Section V we present the results of the fits for R AA for both RHIC and LHC; finally, we present our conclusions.

II. ENERGY LOSS FORMALISM
The building block of the energy loss model is the one-gluon inclusive energy spectrum of gluons with energy ω emitted off a highly energetic colored particle traversing a QCD medium, ωdN/dω [47,53]. In order to compute the energy loss, ∆E, the distribution P (∆E) including all possible contributions, not only from single gluon emission but also from two, three, etc., is needed. The way in which the spectrum, dN/dω, is related to the distribution of lost energy P (∆E) is not exactly known: as usual in perturbation theory, an all-order computation is not possible. Starting with the proposal of independent gluon emission approximation [45], several groups have just iterated the one-gluon inclusive in an independent manner [47,54,55], which has been the most standard procedure to deal with the problem in the last years. An implementation of these ideas, the quenching weights, was worked out in Ref. [47] and will be used in this paper. Note that these quenching weights can be obtained from an iterative solution of DGLAP evolution for medium-modified fragmentations functions in the soft limit [56].
An attempt to include the single-inclusive emission in the form of a kernel in the more sophisticated manner of a rate equation was proposed in [57,58] and incorporated into what is known as the AMY framework. This rate equation still assumed with no proof that subsequent emissions are independent. It is with a more developed study of the role of coherence in jet quenching [35,42] that one can prove that in the limit of gluon formation times much smaller than the medium size, τ form L, a resummation is possible [25,44] recovering, as a particular case, the rate equations by the AMY group. The main limitation of these resummations and the corresponding implementation in the AMY approach, is that they apply to arbitrarily large medium length L, leading to an enhanced energy loss, while finite length effects are relevant in phenomenological implementations, in particular to avoid over-representation of the soft part of the spectrum. The main advantage of the quenching weights is that these finite length effects can be included, although, as mentioned, with no formal proof. So, we will use here the quenching weights [47] as a theoretically motivated and phenomenological sound approach. In addition, this procedure is easy to implement, allows to fit the quenching parameterq once the geometry of the medium is known and, moreover, it has been extensively tested, mainly for RHIC phenomenology [48,54,59,60] but also for the LHC, e.g. [49].
The formalism relies on two basic assumptions: i) the subsequent medium-induced gluon emissions are independent and ii) the fragmentation functions are not modified, i.e., fragmentation takes place in vacuum. These two assumptions find strong theoretical support in the recent analyses of coherence effects in the medium. Starting with the simplified setup of a QCD antenna [35,38,39,42], a pair of color-correlated partons with opening angle Θ emitting a soft gluon, a simple picture of jet quenching arises [40]: a medium of length L and jet quenching parameterq has a typical transverse momentum scale for color correlations Λ ⊥ ∼ 1/ √q L; when the typical transverse size of the jet, r ⊥ ∼ ΘL, is smaller than this scale, r ⊥ < Λ ⊥ , the medium cannot resolve the inner structure of the jet, which remains unchanged, but the whole jet radiates medium-induced gluons with the total charge of the jet. This is the totally coherent case. Clearly, while this picture indicates that the fragmentation function remains basically unmodified if color coherence is maintained, it still depends on the fraction of momentum z and only a global energy loss affects the production of the fragmenting particles. This picture of jet quenching dictated by color coherence is in qualitative agreement with the experimental findings at the LHC [61][62][63][64][65][66] -see [67] for a quantitative analysis of some data.
Regarding the first assumption, as mentioned before, similar color coherence arguments can be used to prove that interference effects, which would break independent emission, are absent when the formation time of the medium-induced gluons, τ form ∼ ω/q, is much smaller than the total length of the medium [25]. So, for soft radiation, ω ω c ≡ 1 2q L 2 , τ form ∼ ω/q L, and independent gluon emission is a good approximation. Notice that the quenching weights and the rate equations are equivalent for the case of soft radiation and when finite energy effects can be neglected, i.e., when the kernel depends neither on the energy of the parent parton nor on the medium length [58].
We use the quenching weights P ( ) tabulated in [68] to model the amount of energy loss of highly energetic partons (or better, coherent jets) which will eventually fragment in the vacuum to give a given hadron h. The corresponding cross section reads where the cross section for producing parton k is (we are neglecting here any difference from parton to hadron rapidities and we take all renormalization, factorization and fragmentation scales to be equal, We make all calculations at NLO using the code in [69] with the proton PDF set CTEQ6.6M [70], the nuclear corrections to PDFs given by EPS09 [71] and vacuum fragmentation functions DSS [72,73]. The quenching weights, P ( ), are defined as where the medium-induced radiation spectrum is given by 1 with the path integral (Semi)analytical solutions of equations (4) and (5) are only known in a limited number of cases. Two main approximations are normally employed. One of them consists of expanding the exponents in series of the opacity parameter n(ξ)σ(r) with a Fourier transformed cross section given by thermal QCD in some approximation, for example with µ the thermal mass. The first term of this opacity expansion is normally employed, as used, for example in Ref. [11].
The leading term of the Fourier transform of Eq. (6) is proportional to r 2 with a logarithmic correction. The second approximation consists of neglecting this logarithmic correction and approximate With this approximation, the path integral is that of a harmonic oscillator with imaginary (eventually time-dependent) frequency. The solutions of this path integral whenq(ξ) ∼ 1/ξ α can be found in Ref. [47]. This second approximation, sometimes known as the multiple soft scattering approximation, or BDMPS approximation, will be used in this paper. Eqs.
(5) and (7) can be considered as our definition of the transport coefficientq.
The main difference between the two approximations, at the analytical level, is the presence of perturbative, power-law, tails in the opacity expansion, which are absent in the multiple soft scattering one. Notice that the AMY approach would correspond to a resummation of the multiple scatterings with the correct cross section (6) but without an interference between vacuum and medium radiation which turns out to be very relevant for finite medium lengths, making the soft part of the spectrum non-divergent when real angle emission is impossed for the emitted gluons.

III. FROM A HYDRODYNAMICAL PROFILE TO THE TRANSPORT COEFFI-CIENT
The distributions of the energy loss, ∆E, of a fast quark or gluon, i, traversing a medium is computed through the quenching weights, P i (∆E/ω c , R) tabulated in [47,68] for the case of a static medium of finite length L and transport coefficientq, where For the expanding medium case, with proper time dependence of the transport coefficient, which, in particular, reproduce Eqs. (8) for the static case. Similar implementations of the hydrodynamical model have been used before [15,[48][49][50], so it is a rather standard procedure 2 . The production point of the parton at time τ prod is distributed according to an N coll -scaling in the transverse plane and the azimuthal angle φ is taken as a random number in [0, 2π]. As usual, each parton traverses the medium in a straight-line trajectory parametrized by the proper time ξ at each point in the transverse plane. We only need to specify the relation between the local value of the hydrodynamical variables at (x ⊥ (ξ), y ⊥ (ξ)) and the local value of the transport coefficientq(ξ). Following our previous work [48] we where K 1 would correspond to the ideal QGP (see the estimate in e.g. Ref. [51]). Other relations between the transport coefficient and the local thermodynamical quantities have been explored e.g. in Ref. [50]. The local energy density (ξ) is taken from a hydrodynamical model of the medium, for which we will consider several different options in the next sections.
The rest of the formalism follows that in Ref. [48]; we include here some details, while the complete formulation can be found in that reference.
In a dynamical medium like the one considered here, there is an ambiguity on the value of the transport coefficient, defined by Eq. (11), for values smaller than the proper time τ 0 when relativistic hydrodynamics is started. One extreme case is to takeq(ξ) = 0 for ξ < τ 0 . The absence of any energy-loss effect for these early times is a strong assumption since neither thermalization nor isotropization is necessary in the approach in which the quenching weights have been computed. To quantify this uncertainty, we consider three different extrapolations for the time from the hard parton production to the thermalization: These extrapolations range from the most extreme assumption of no effect at all before the thermalization time (case (i)) to a continuous interaction from the production time (taken to be τ prod 0.04 fm/c) and a free-streaming medium with energy density dropping as (ξ) ∼ 1/ξ (case (iii)).
2 Notice that we have slightly changed the prescription to compute R ef f which is now the second moment ofq(ξ). The results are similar with the old prescription (see e.g. Eqs. (4.2)-(4.4) in Ref. [48]) but with improved stability for functional dependences ofq(ξ) that are divergent in 1/ξ.
The production point of the hard scattering is characterized by a production weight w(x 0 , y 0 ) computed as where T A are the profile functions computed from a 3-parameter Fermi distribution at a given impact parameter b taken from [74]. This weight allows to compute the average fragmentation functions for a parton k which has traversed the medium and hadronizes in the vacuum to a hadron h as where P k (x 0 , y 0 , φ, ζ) is the quenching weight for parton k and the normalization is N = 2π dx 0 dy 0 w(x 0 , y 0 ). With this definition, the cross section, Eq. (1), can be simply computed Let us emphasize again that the formalism assumes no medium modification of the fragmentation function, an early approximation which we have motivated here by color coherence arguments. For this reason, Eq. (13) should be understood as a suitable computational tool to include the energy loss of the coherent colored object traversing the QCD medium which eventually fragments in the vacuum.

IV. HYDRODYNAMICAL MODEL OF THE MEDIUM
We obtain the space-time distribution of the local energy density by solving the relativistic hydrodynamic equations. Such simulations require the specification of initial values for the energy momentum tensor, as well as parameters that describe medium properties, neither of which are accurately known. In order to test the robustness of our results and conclusions with respect to these uncertainties, we repeat all calculations using space-time profiles from several different hydrodynamic simulations.
The first, which we refer to as "Hirano", corresponds to the calculation described in [48,[75][76][77], to which we refer the reader for details. In short, this calculation uses an optical Glauber model where the initial entropy density at initial proper time τ 0 = 0.6 fm is given by a linear combination of the number density of participant nucleons, ρ part , and binary collisions, ρ bin : to which we again refer the reader for all relevant details.
One calculation, which we refer to as "Glauber", uses for an initial condition an energy density proportional to the density of binary collisions, ρ bin , while the ratio of shear viscosity to entropy density is fixed to a constant value of η/s = 0.08.
The final calculation is referred to as "fKLN". This simulation takes its initial condition from a factorised Kharzeev-Levin-Nardi model [80], with the shear viscosity set to η/s = 0.16.
Both of the latter simulations begin at an initial proper time of τ 0 = 1 fm and use an equation of state inspired by lattice QCD calculations. Each system is assumed to be in chemical equilibrium until it reaches a freeze out temperature of T f = 140 MeV.
All of these calculations have been successfully tested against various experimental data, but use different choices for initial conditions, thermalization time, viscosity, equation of state, etc. Thus, we expect that the variation in our results from using these different models should give a reasonable indication of the uncertainty coming from the hydrodynamic background. We will see that such uncertainty is negligible with respect to our main conclusions.

V. RESULTS OF THE FITS
We restrict our study here to the case of one-particle inclusive suppression at RHIC [81] and the LHC [82], i.e., we stick to the simplest observable R AA . We do not consider here other observables previously considered in fits like this, see e.g. [48], as they may involve other effects related with fragmentation, mass effects on the energy loss mechanism, etc. Something new is, however, the centrality dependence of both RHIC and LHC data. Most of previous analyses have only studied the most central class or analyzed the centrality dependence only for one energy [15,49,50].
We have performed our study using three different hydrodynamical profiles and three different assumptions for the time prior to the equilibration time, see the previous Sections.
For economy of space we present here only a subset of the total results which, anyway, lead to the same qualitative conclusions (only the actual value of the parameter K depends on the given prescription). In Fig. 1   The uncertainty band is determined by ∆χ 2 = 1. In order to make the comparison between RHIC and the LHC, these issues need to be taken into account, although the conclusions do not change at the qualitative level. In the left panels of Fig. 3, Fig. 4 and Fig. 5 we plot the different values of the K-parameter fitted to the PHENIX data [81] for different combinations of hydrodynamical profiles and behavior before the thermalization time. The corresponding values for the LHC [82] are plotted in the right panels of the same figures.   Several comments are in order. First, the extracted values of K are compatible for the cases of either frozen energy density or free streaming before τ 0 , and the results for the three different hydrodynamical implementations are similar. This is not the case when no quenching is assumed before τ 0 ; for this assumption, the two viscous hydrodynamic models that use a common (larger) τ 0 require a larger K than the ideal hydrodynamic model that considers a smaller τ 0 , with actual values which become unrealistically large. Therefore, we do not consider the results obtained for this assumption for the discussion of the values of K, but the qualitative behavior that we find is in agreement with the two other assumptions.
In any case they clearly illustrate the importance of the treatment of early times in jet quenching computations. Second, for the most peripheral collisions at the LHC, model "Glauber" demands a much larger K than the others, while model "Hirano" returns a rather flat value of K for all centralities. Third, the trend of the results at RHIC is a slight decrease with decreasing centrality, although compatible with constant, while at the LHC the behavior is constant except for the smaller centralities, where the behavior, as it was mentioned above, depends very much on the hydrodynamical profile employed.
In the end, we would like to understand the systematics and relation of LHC and RHIC results for the K-factor that we obtain. First, we notice that, in principle, Eq. (11) de-termines how far or close the perturbative estimateq 2 3/4 is from our value fitted to experimental data. In this sense, we note that there is a clear departure from unity of this value for the case of RHIC. This fact was found several times [48,50] 3 . We also find that the corresponding value of K is smaller at the LHC, a fact which has been already found before by other groups [15] but with a smaller decrease (a factor ∼ 25 % compared to our factor 2-3) 4 . The study of the centrality dependence is, nonetheless, more surprising. The extracted value of K seems to depend mainly on the energy of the collision and much less (if any) on the centrality. This is not the behavior one would expect from a naive interpretation in which the K factor only indicates the departure from the leading order perturbative estimate determined by temperature. In this naive interpretation, a medium with a smaller temperature (RHIC) would need higher orders of the perturbative series to be included, while a medium at higher temperature would be closer to the ideal limit. This simple interpretation does not correspond, notwithstanding, to the present findings as there is an overlap on typical energy densities between central AuAu at RHIC and semi-peripheral PbPb at the LHC, so their values of K should coincide in this naive interpretation. In order to provide an estimate of this overlap, we plot, in Figure 6, the K-factors obtained for different centralities and energies versus an energy density times formation time τ 0 extracted from the experimental data using Bjorken estimates [83,84] -we have checked that the overlap is similar if we plot as a function of the maximum energy density of the hydrodynamical profiles that we have used to perform the fits.
Finally, we would like to include here the predictions of our formalism for the forthcoming data on PbPb collisions at the LHC at √ s NN =5.02 TeV. In our case, the K-factor is not fixed but fitted from experimental data. However, assuming that the perturbative estimatê q 2 3/4 is approximately correct, it would be reasonable to think that our K-factor cannot be much smaller than unity. So, assuming the same values of K as the ones we obtain from the fit to √ s NN =2.76 TeV data, K=1.133 ± 0.028 for "Glauber" hydrodynamic profile and K=1.088 ± 0.028 for "fKLN" hydrodynamic profile and considering a frozenq between production time and τ 0 , we obtain a slightly stronger suppression for the 5.02 TeV case, 3 Note that the difference of the present extraction K ∼ 2-3 and K ∼ 4 from [48] comes mainly from the new definition of R in (10), indicating, again, the important role of the geometry in the extraction ofq. 4 Nevertheless this comparison needs to be taken with caution as the values ofq/T 3 quoted in [15] are performed at a given temperature and no systematics with temperature is presented. Moreover,q is not the natural fitting parameter in the models studied in that reference but a derived quantity once the parameters of the different models are extracted from the data. In our case,q is the natural parameter, given by Eqs. (5) and (7), and the K factor has a well defined meaning. This behavior is easy to understand considering that the transverse momentum slope is not that different at the two LHC energies, but the larger energy densities at the higher collision energy imply larger absolute values ofq for the same K-factor. Thus, the energy loss is larger and the suppression becomes stronger. The effect is, however, not very large.

VI. CONCLUSIONS
We have studied the one-particle inclusive suppression of particles produced at high transverse momenta at RHIC and the LHC as a function of centrality. By defining a constant K-factor with respect to the perturbative estimateq 2 3/4 we fit the corresponding experimental data at RHIC and LHC for different centralities. The fitted value at RHIC confirms previous estimates [48,50]  the static scenario equivalent to the evolving medium, Eq. (10). For the case of the LHC, instead, the extracted value of K is close to unity. One would be tempted to make the naive interpretation that the medium created at the LHC, having a larger temperature, is closer to the ideal case than the one at RHIC, for which larger corrections or even a strongly coupling treatment, could be needed.
This naive interpretation finds difficulties to be accommodated, however, with the fact that the centrality dependences at RHIC and the LHC separately are rather flat, that is, the change in the value of K is not simply due to the different temperature (or energy density), as there is a large region of overlap between RHIC and the LHC for different centralities.
At this moment we do not have an interpretation for this finding which, in any case, should be checked by other model implementations of jet quenching. It is also worth noticing that the extraction of the value of K in the case of RHIC depends on a single set of experimental data, namely inclusive π 0 suppression measured by PHENIX. The corresponding results from STAR on π + + π − suppression [85] show a smaller suppression but the smaller range of transverse momentum studied makes our analysis to be not very reliable. For this reason we have chosen not to include this set of data in the fit. For the LHC, on the other hand, CMS [86] and ATLAS [87] have measured the suppression of inclusive charged particles with results almost identical to the ones from the ALICE collaboration 5 .
From a theoretical perspective, the formalism of the quenching weights presents several limitations which could influence the result. We quote some of them here: (i) the very definition ofq in the path integral, Eq. (5), neglects the perturbative tails of the distributions which may enhance the energy loss and even change its angular dependence; (ii) as we have mentioned, the quenching weights rely on two assumptions which could fail if color coherence of the parton shower is broken during the path of the jet though the medium; (iii) the geometrical implementation of the hydrodynamical profiles relies on the relations (10) which have been proven only for a class of profilesq(τ ) ∝ 1/τ α ; (iv) finite length corrections to the independent gluon emission are not known in any of the implementations used at present; (v) finite energy corrections to the medium-induced gluon radiation could also affect the result; (vi) the jet quenching parameterq is taken to be energy or length independent, while evolution equations have been proposed [88][89][90]; (vi) finite energy corrections could also contain collisional energy loss which is neglected in our formalism and which may have a different parametric dependence with the medium properties. In spite of these limitations, it is difficult to imagine how a more refined implementation of the in-medium parton shower could qualitatively modify our finding of a mostly flat, in centrality, value of K and different for different collision energies. At the level of the partonic spectra, the main quantity affecting the suppression is a decreasing value of the slope with increasing energy: for a simple parametrization of 1/p δ T , δ varies in the range 5 − 7 from the LHC to RHIC. This steeply falling spectrum introduces a bias in the proved energy loss distributions, so the typical energy of the partons at the LHC is larger than at RHIC for the same measured p t . In this way, a softer part of the energy loss distribution is probed with increasing δ, so that the perturbative tails neglected in the multiple soft scattering approximation used may become relevant. Moreover, the typical jet could also be more collimated at lower energies (as the typical fraction of momentum probed in the fragmentation function is also larger).
These details could modify the value of K extracted for a more realistic treatment of jet coherence and thermal cross sections. A better control on the initial times and the study of different experimental observables with refined methods will allow to clarify this issue.
As an outlook, we plan to study the effect of event-by-event fluctuations as presently done e.g. in [91][92][93].
Note added: When finalizing the preparation of this manuscript, new, preliminary, experimental data on inclusive charged particle suppression R AA of PbPb collisions at √ s = 5.02 TeV measured by the CMS collaboration were presented at the 3rd International Conference on the Initial Stages in High-Energy Nuclear Collisions (InitialStages2016) in Lisbon (Portugal) [94]. Although we have not included CMS data in our analysis, for the reasons explained above, taking at face value, the new data would indicate a suppression which is almost energy independent at the LHC. In Fig. 7 we predict a ∼ 15% larger suppression at √ s = 5.02 TeV compared to √ s = 2.76 TeV. We have checked that, if confirmed, this same suppression would indicate that the K value needed to reproduce the higher energy LHC data is ∼ 10% smaller than the ones quoted here. Obviously, the exact value would require a new fit once the data for central rapidities are available.