Quenching effects in the cumulative jet spectrum

The steeply falling jet spectrum induces bias on the medium modifications of jet observables in heavy-ion collisions. To explore this bias, we develop a novel analytic framework to study the quenched jet spectrum, and its cumulative. We include many energy-loss-related effects, such as soft and hard medium induced emissions, broadening, elastic scattering, jet fragmentation, cone size, coherence effects, etc. We show that, different jet spectrum-based observables are connected, e.g., the nuclear modification, spectrum shift, and the quantile procedure. We present the first predictions for the nuclear modification factor and the quantile procedure with cone size dependence. As an example, we compare dijet and boson+jet events to unfold the spectrum bias effects. We improve quark-, and gluon-jet classification using arguments based on the cumulative. Besides pointing out its flexibility, we apply our framework with other energy loss models such as the hybrid weak-, strong-coupling.


Introduction
At high-energy particle colliders, such as LHC at CERN, the exchanged momentum is large enough to resolve and scatter partonic constituents of the matter. These violent scatterings deviate the partonic constituents of nucleons and allow for intense bremsstrahlung radiation that ultimately result in collimated bunches of hadronic particles and energy. These so-called jets open a new perspective on the understanding of strong interactions at intermediate scales between the scale of the hard partonic scattering and the hadronization scale. Jets can be calculated with high-precision within perturbative QCD both gin electron-positron and proton-proton collisions [1][2][3].
For accelerated heavy nuclei, processes involving large momentum exchanges between the incoming partons have additionally the potential to be well-calibrated probes of the hot and dense quark-gluon plasma (QGP) that is formed in the aftermath of such violent collisions. In this context, jets are particularly interesting since their typical formation time-scales overlap with the time-scales governing the creation and evolution of the QGP, suggesting potentially substantial jet-medium interactions. While the strong separation of the medium scale and the jet scale motivates a perturbative description of hard jet-medium interactions, many aspects of this processes are in the realm of non-perturbative physics and have to be modeled on the phenomenological level. A satisfactory description of jet production in heavy-ion collisions would therefore allow us to separate perturbative from non-perturbative phenomena. Besides, a well-controlled scale separation is an indispensable insight when studying, e.g., jet-medium coupling, thermalization, or medium modification of hadronization.
Jet studies have a rich history and wide selection of observables have been discussed involving fully reconstructed jets or their substructure, for recent reviews see [1,2,4]. A fundamental observable is the momentum spectrum of jets for different reconstruction parameters R [5,6]. The nuclear modification factor R AA , compares the spectrum in heavy-ion collisions (AA) to proton-proton (pp) at the same reconstructed jet p T . However, jets that interact with a surrounding medium lose energy and end up with smaller p T . Thereforeand it is not emphasized enough-the jet selection for R AA compares two jet populations originated at different p T . The equal p T selection induces a bias in the observables because the probability of creating a jet is steeply falling with p T [7,8]. We explore the bias in this paper.
Recently, there has been efforts toward mitigating such bias effects by investigating novel observables or by using machine learning techniques [9]. One alternative, that we will investigate in detail in the current work, is to introduce a quantile procedure [10] to reconstruct a p T that is closer to the initial jet p T before quenching sets in. In contrast to R AA , the quantile procedure uses the tail-cumulative of the jet spectrum and momentum ratio to reduce the bias coming from the steepness of the spectrum. We demonstrate the properties of the quantile procedure within a versatile framework to incorporate quenching effects and explain its robustness for the first time.
New measurements have been suggested, e.g., boson+jet, to have a better handle on the momentum mitigation in heavy-ions [11][12][13][14]. Bosons suffer less medium effect, and their momenta are related to the initiator of the recoiling jet. The jet spectrum in boson+jet events is slightly different from inclusive QCD jets, and we will use it to illustrate the bias on the quenching. Moreover, quark-and gluon-jet contribution in the inclusive and boson+jet samples differ. This can be used for quark-gluon jet discrimination in a modelindependent fashion [15]. Using arguments on the cumulative, we improve the statistics of the classification task. Parallel with the works mentioned above, some numerical studies also appeared using Bayesian and machine learning techniques to extract the energy loss properties from data [16].
During the evolution of jets inside a QGP, their constituents scatter elastically and inelastically on the medium. The scatterings redistribute energy to larger angles out of the jet cone, resulting in energy loss. The inelastic (or radiative) part describes the emissions induced by the medium (medium-induced emissions, or MIE for short). The MIE has wellknown limits: (i) in the multiple soft scattering limit, the scattering centers act coherently, resulting in suppression of emissions (QCD analog of the Landau-Pomeranchuk-Migdal effect), captured by the BDMPSZ formula, which describes the induced emission of soft gluons [17][18][19]. (ii) the single hard scattering limit is captured by the GLV formula, which describes emissions of harder gluons [20]. Beyond the analytic limits, MIE is also amenable to direct numerical methods, such as in Refs. [21][22][23]. Recently, there has also been a progression in the better understanding of the two regimes and the scales involved [24][25][26]. In our work, we adopt the latter strategy to explore the impact of the MIE spectrum. The MIE can be resummed accounting for multiple induced emissions. The resulting formalism is the quenching weight [27,28]. It is easy to generalize for all jet constituents, including coherence effects [29], spectrum shapes, and elastic energy loss. Within this framework, we show the appearance of the spectrum bias.
The paper is organized in the following; in Sec. 2 we define our novel framework to calculate the medium jet spectrum. In Sec. 3, for the first time, we show the properties of the cumulative spectrum, and we apply the quantile procedure. In Sec. 4 we give predictions for measurements, considering cone size dependence for the single-inclusive jet sample produced in dijet and boson+jet events (we focus concretely on Z+jet). We show how to use cumulative arguments to improve quark-gluon discrimination. Besides in App. D, we use the quenching weight formalism for elastic scatterings and in the hybrid weak-, strong-coupling model.

Quenching effects in the spectrum
The main observable considered in this work is the single-inclusive spectrum of reconstructed jets in heavy-ion collisions. In this work, we will both consider jets produced in conventional QCD processes, that is dijet events, and jets produced in conjunction with a photon or weak boson, so-called boson-jet events. In the context of high-energy collisions, it is natural to assume a factorization of the partonic hard cross-section from the subsequent medium processes. This can be justified by invoking the large separation of momentum scales involved in jet production; typically the hard scattering Q hard ∼ 10 3 GeV, is much bigger than the jet scale Q jet ∼ p T R ∼ 10 2 GeV, where p T is the reconstructed transverse momentum of the jet and R is jet cone parameter. These scales are much bigger than the typical medium scale, for instance, the temperature of the medium T ∼ 0.5 GeV. Hence, one can write the medium modification of the vacuum jet spectrum due to energy loss [27], The P > (ε) describes the probability of a vacuum jet to distribute (or lose) energy out of the jet cone. Above, dσ vac R /dp T refers to the partonic cross-section to produce a jet with R in the collinear factorization and can be calculated up to high precision [5,6]. In this work, instead, we extract the spectrum using a Monte Carlo event generator, see Sec. 4 for further details.
The vacuum spectrum is well approximated by a power-law, dσ vac R /dp T ∝ 1/p n T , and is steeply falling, i.e. n 1. Therefore, we approximate dσ vac R (p T + ε)/dp T = A (p T + ε) −n ≈ e −nε/p T dσ vac (p T )/dp T , where A is a constant and n ≡ n(p T , R) is the power index of the spectrum. In the last step, we additionally assumed that ε p T . The cone-size dependent nuclear modification factor defined as is consequently related to the quenching factor Q(ν ≡ n/p T ), which is the Laplace transform of P > (ε), i.e.
This approximation is precise within at most a few percent for realistic parameters in a wide kinematic range 10 < p T < 1000 GeV that we consider here. Therefore, in what follows, we simply identify R med = Q > . For further details on such corrections, see Appendix A.
In order to go from R med to the experimentally measured R AA , one needs to include a p Tdependent quark and gluon jet production. Other differences are mostly due to geometry (for a review see Ref. [30]), and cold nuclear effects (nPDF), that we include in Sec. 4. The energy loss distribution P > (ε)-and therefore the quenching factor Q > (ν)-depends on vacuum jet properties (like the jet p T and the cone size R) and the properties of the medium (e.g., medium length L, and the jet transport coefficientq). Its normalization condition, ∞ 0 dε P > (ε) = 1, translates to Q > (0) = 1. The assumptions underlying Eq. (2.2) are quite robust for a wide range of applications. Therefore, one is flexible in defining the precise nature of the energy loss distribution P > (ε). The introduction of a probability distribution to describe effects of quenching in heavy-ion collisions has a long history in the analysis of single-inclusive hadron [27,28,31] and jet spectra [29,[32][33][34].
We will derive the quenching factor Q > of the jet in several steps. First, we consider the induced radiation spectrum of a single color-charge propagating through the medium and how to account for their multiple emissions. We thus arrive at the quenching factor for a single parton, Q (0) > (ν). Next, we consider the effect of jet fragmentation which leads to multiple vacuum-like emissions on short time-scales inside the jets. Partons from these emissions contribute to the quenching of the full jet. This is accounted for by the so-called collimator function that provides a fully resummed quenching factor Q > (ν). The corresponding jet quenching probability distribution can then be found via an inverse Laplace transform, but we will not pursue this further in this work.
In our numerical results in Sec. 4, we will also include elastic energy loss, since it potentially can contribute to the ∼ 10 − 30% level to the final jet suppression factor. As we said, the formulation above is quite general and allows to separately formulate a) the mechanism of quenching, and b) the phase space where the jet is affected. In App. D, therefore, we show how to formulate other energy loss models in terms of quenching weights.

Constructing the quenching weight
We construct the radiative energy loss distribution P > (ε) (or the quenching factor Q > (p T )) of a jet starting from a single parton, that we denote P (0) > (ε). For the radiation of a single medium-induced gluon, this energy loss probability is simply given by where dI > /dω is the spectrum of medium-induced gluon radiation spectrum that emerges at angles larger than the jet cone, or This is an important difference with respect to the more common use of the quenching weights applied to single-hadron spectra in Refs. [27,28]. The first term in Eq. (2.4) describes an emission, while the second term is a virtual correction and provides the normalization. Our starting point in Sec. 2.2 is to discuss the specific details of this induced-emission spectrum. Accounting for multiple such emissions in course of the medium propagation allows us to derive the single parton quenching weight in Sec. 2.3. Then, having derived how one parton contributes to the energy loss of the whole jet, we next turn to the calculation of how multiple partons in the jet, resolved by the medium during their fragmentation process, add up to the total quenching effect. This will be described in Sec. 2.4.

The medium induced gluon spectrum at finite cone
The medium-induced spectrum from multiple scattering in a QCD medium was derived independently by Zakharov [18,35] and Baier-Dokshitzer-Mueller-Peigne-Schiff (BDMPS) [17,19,36], see also Refs. [37,38]. In the limit of soft gluon emission with ω energy and k transverse momentum, this spectrum can be written as where C R = C F if the radiator is a quark (C A for gluon), the path integral is and v( is related to the elastic scattering potential in medium. The leading logarithmic behavior of the scattering potential reads for any hard Coulomb tailed elastic potential, whereq 0 (t) is a local transport coefficient and µ * is related to the infrared scale that regularizes the in-medium scattering potential d 2 σ el /d 2 q. For example, the Gyulassy-Wang scattering potential [39], describes scattering off a plasma made up of Yukawa-screened quasi-particles, reads where µ is a screening mass, and n(t) is the density of scattering centers in the medium. In this caseq 0 (t) = 4πα 2 s N c n(t), and µ 2 * = µ 2 e −1+2γ E /4. The path integral in Eq. (2.7) can be solved numerically, see e.g. in Ref. [21][22][23], but analytic solutions are available in limiting cases. Here, two of the most frequent schemes are discussed. In the "harmonic oscillator" (HO) approximation, valid when the logarithm is slowly varying around a constant scale, i.e. µ 2 * 1/x 2 ∼ Q 2 , one can absorb the logarithm in the definition of the transport coefficient so that v(x) ≈qx 2 /4. This corresponds to a purely Gaussian momentum broadening in transverse momentum given by k 2 =qt during the propagation in the medium. In this case, Eq. (2.7) describes a harmonic oscillator with imaginary frequency Ω 2 =q/(2iω) in transverse to the propagation, and whose solution is well known. This approximation, however, fails to describe the hard tail of medium emissions. On the other hand, one can also truncate the resummation of medium scatterings at a fixed order N , giving rise to the so-called "opacity expansion" [20,40]. The N = 1 approximation amounts to considering a single, incoherent scattering with the medium during the propagation.
A discussion of the regions of validity of these approaches was recently addressed in Refs. [24][25][26]41]. In the soft scattering regime, the formation time of emissions t f = 2ω/k 2 becomes modified due to Gaussian broadening, i.e. t f ∼ 2ω/q. For emissions with large formation times t f ∼ L correspond to ω ∼qL 2 /2 ≡ ω c . It follows that emissions with ω > ω c cannot be produced by soft collisions and Gaussian broadening, and are dominated by a single, hard scattering with the medium constituents. Similarly, at short formation times of the order of the medium mean free path t ∼ λ L, or ω ∼qλ 2 /2 = ω BH ω c , the spectrum is again dominated by single scattering [40,41]. Since this latter regime gives a small contribution to energy loss, we will not discuss it further here.
Recently, the contribution of hard emissions has been shown to matter for precision comparisons with high-p T single-hadron spectra at RHIC and LHC [42]. A systematic procedure to calculate the spectrum for a large range of relevant emission energies ω > ω BH was developed in the so-called "improved opacity expansion" (IOE) [24][25][26]. This framework rewrites the leading-log scattering potential form Eq. (2.9) as where Q c is a separation scale of the harmonic potential. In the limit of Q c µ * , one can then expand the solution of the path integral in Eq. (2.7) around the HO solution with an effectiveq(Q 2 c ) and treat hard scatterings with the medium, given by δv(x, s), as higher-order perturbations. This approach is systematically improvable and, up to next-tonext-to-leading order in this expansion, the effectiveq parameter is given by [26] q(Q 2 c ) =q 0 ln where the scale Q 2 c is itself found through an implicit equation, is the Lambert function on the i th branch. Since jet quenching is not sensitive to the details of very soft gluon emissions, in our numerical results we freeze the logarithms at 1, i.e. ln Q 2 c /µ 2 * ≥ 1. The IOE has so far only been developed for the spectrum integrated over transverse momenta, i.e. dI/dω, and for the momentum broadening of a single particle in the medium [43]. Since the contributions to jet energy loss rely on out-of-cone emissions, cf. Eq. (2.5), we instead have to consider the matching of the partially integrated spectrum dI > /dω between the multiple-scattering HO and single-scattering N = 1 regimes. We use the quenching parameterq from Eq. (2.12) which correctly connects theq 0 parameter from the elastic potential with the multiple scattering formalism. We propose a simple interpolation scheme that relies on a single matching scale ω , that is related to the broadening of soft gluons and will be defined below. At small gluon energies, below the matching scale, the spectrum is described by the HO approach with the effectiveq. Above the matching scale, the spectrum is given by the N = 1 spectrum.
We will treat the medium as a "brick" of constantq 0 and fixed length L. In the absence of the cone constraint, a natural matching scale is ω ∼ ω c ≡qL 2 /2. For the out-ofcone spectrum, the effect of broadening after emission cannot be neglected. For Gaussian broadening that presents in HO, a particle emitted at initial time accumulates k 2 ∼qL after propagating through the medium. This corresponds to an angle θ ∼ √q L/ω in the small-angle approximation. Demanding that this angle is larger than the jet cone R for energy loss, a cut-off in energy arises ω < ω R ≡ √q L/R, above which the HO spectrum falls rapidly. For more details, see App. B. It turns out that a relatively smooth matching between the HO and N = 1 regimes is achieved by choosing where ω c and ω R are defined with the effectiveq parameter in Eq. (2.12). For our final results, see the left and right panels in Fig. 1. The postulated matching works extremely well, up to some negligible discontinuities in the spectrum. The regime of soft gluon emissions, ω < ω , is dominated by multiple scattering where we can employ the HO approximation. The spectrum in this approximation is given by whereᾱ = α s C R /π and the two factors read [44] where Ω = (1 − i) q/ω/2. The integrated spectrum in Eq. (2.5) reads then where now We find that in the R → 0 limit, Eq. (2.17) yields which is the celebrated BDMPS-Z spectrum [18,19]. As discussed above, at ω > ω , the HO spectrum has to be corrected with the single hard gluon emission spectrum (N = 1) [24][25][26], for which Using Eq. (2.5) and the Gyulassy-Wang potential from Eq. (2.10), we immediately find that the integrated spectrum reads where ζ ≡ ωR 2 L/2 and y = µ 2 L/(2ω). Again, for R → 0, we recover the familiar form of the integrated N = 1 spectrum [20,28], as expected.
The full spectrum is therefore postulated to be well approximated by the following interpolation, where all parameters in the HO spectrum and in the matching scales contain the effectivê q, given in Eq. (2.12). Note that the spectrum dI > /dω depends on the initial energy only in the combination ω = x(1 − x)E ≈ xE.
At this point, we pause to discuss the choice of medium parameters. We choose the medium coupling to be fixed at α s = α med = 0.3, corresponding to g med ≈ 1.94. The IR cutoff scale is µ ≡ 2/3g med T 0 ≈ 0.71 GeV, and ω min ≈ 1.5 GeV. The remaining parameters are chosen to reflect the conditions in 0-10% central PbPb events at √ s N N = 5.02 TeV, see Table 1. As demonstrated below, this parameter set gives a good description of the experimentally measured jet suppression factor, see Sec. 4. Since our work does not deal with the precise description of experimental data, we have not attempted to fix these parameters from a model of the underlying medium nor fitted them to experimental data, which was done in [34]. This choice finally leads to a matching scale Q 2 c = (1 − 100)µ 2 and q = (1 − 5)q 0 , depending on ω, in the IOE-matched spectrum Eq. (2.24).
The MIE spectrum is shown in Fig. 1 for quarks (left) and for gluons (right). The matching points ω * from Eq. (2.13) are shown with bullets and below (above) the spectrum is the HO (N = 1) spectrum. The matching works very well capturing the cone size dependence, however it is not perfectly smooth (see R = 0). It is good enough to study the integral of this spectrum, presented in the quenching factor. With different colors the cone size dependence points out, it is less probable to lose energy by opening the cone (i.e., recapturing emissions). The difference in the quark and gluon spectrum is the color factor C A /C F = 9/4, and thus gluons lose more energy. The dotted line in Fig. 1, represents the energy scale which below secondary branching start to dominate [37], corresponding to ω ∼ α 2 sq L 2 , see Eq. (2.29) (for more details, see the next subsection). Finally, the grey band in Fig. 1 corresponds to emissions with ω ∼ ω BH , which are given by the Bethe-Heitler spectrum [40,41]. In what follows, we will neglect such emissions since these emissions do not contribute significantly to jet energy loss at high-p T [27, 28].

Single-parton contribution to out-of-cone energy loss
When the number of medium-induced gluon emissions becomes large, ∞ ω dω dI dω > 1, one needs to go beyond (2.4) and account for multiple emissions to the energy loss distribution. Assuming independent emissions, we can treat it as a Poisson process. This allows to define a probability distribution of radiating energy ε off a single parton at angles

Medium parameter
Value Table 1. Choice of medium parameters corresponding to 0-10% central PbPb collisions with √ s N N = 5.02 TeV.
-9 - larger than the jet cone. This probability distribution reads where the parton radiates n soft gluons with ω j energies summing up to ε, and there is a Sudakov exponential factor to resum virtual contributions [27]. In Eq. (2.25), we only account for the primary emissions off the leading particle and neglect any secondary splittings. Performing the Laplace transform, the quenching factor of a single parton emitting multiple gluons can be calculated using where ν = n/p T . It is here worth emphasizing the role of the hard emissions described by the N = 1 spectrum. Neglecting for the moment broadening, i.e. setting R = 0 in Eq. (2.22), the spectrum at large ω > ω c µ 2 L/2 is simply The resulting single-parton quenching factor for this regime behaves as nω c , the quenching becomes at most Q (0),N =1 (ν)| p T nωc ≈ 1 − α s C Rq0 /(2q) constant factor. Hence, we conclude that the impact of hard radiation at LHC, where we expect q >q 0 , is relatively small. However, at high-p T , i.e. p T nω c , where the leading behavior is , it becomes more important. Also, since hard radiation takes place at small angles, the effect is even smaller at R > 0.
Equation (2.26) also includes contributions from small energy gluons, and an important modification should be included to improve their description in the medium. As mentioned above, gluons emitted with energy ω BH < ω < ω s , where will thermalize quasi-instantaneously in the plasma via multiple branching [37]. Their energy will basically be redistributed randomly over a cone with characteristic opening angle R rec ∼ π/2. Hence, instead of losing energy ω out of the cone, the jet loses ω where the second power comes from the area proportionality. This process describes the thermalization of soft jet particles. Moreover, if R rec (ω, η, φ), where (η, φ) describes the jet direction with respect to the reaction plane, one could use it to describe back-propagation of the thermalized energy to the cone, and thus medium response. The dashed line in Fig. 1 shows the location of ω s ≈ 8.5 GeV for our parameters. Depending on the medium and jet parameters, this scale can be below or above the matching scale ω . In Ref. [34], the importance of R rec was studied and small dependence in the result was observed and thus we used R rec = π/2. Finally, after neglecting emissions below the Bethe-Heitler energy (we assume that it is given by the plasma temperature ω BH = T 0 ), our final form for the quenching factor of a single parton inside the jet is therefore where the cone size dependence is implicit in the integration limits of the out-of-cone spectrum. The single parton quenching factor Eq. (2.30) is shown in Fig. 2 with dashed lines for quark and gluon initiators. By opening the cone, the emitted energy gets gradually recovered, and thus the quenching factor becomes closer to 1. The difference in between quark and gluon initiators is Q resulting more quenching for gluon. We expect our description to be less valid at smaller energies.

Quenching the whole jet
Having derived how one parton contributes to the energy loss of the whole jet, we next turn to the calculation of how multiple partons in the jet, resolved by the medium during their fragmentation process, add up to the total quenching effect.
Due to the large phase space for radiation between the jet scale ∼ p T R and the hadronization scale ∼ Λ QCD , the jet forms through multiple emissions. It can be estimated, from formation time arguments, that many of these emissions occur while the parton(s) are interacting with the surrounding medium [29]. In the limit of complete decoherence, the constituents are expected to lose energy independently, following an incoherent superposition of single-particle quenching factors in Eq. (2.30). However, one has account for coherence effects leading to a finite resolution power of the medium. As long as two partons are closer to each other than the medium resolution length, the medium cannot resolve them individually. The two partons are affected coherently (as a whole color charge) by medium interactions and, in particular, by induced energy loss [45]. The relevant time-scale can be estimated by comparing the size of a dipole, that in the small-angle approximation roughly scales as r ⊥ ∼ θt (where θ is the angle of the dipole and t is the propagation time), to the resolution length of the medium, that scales as where we assume Gaussian transverse-momentum broadening. The two transverse sizes are equal at the decoherence time t d ∼ ( This condition can be translated to emission times: emissions with formation times smaller than the medium decoherence time t f < t d are vacuum-like. In other words, they are generated according to the probability distribution to split in the vacuum. The core constituents should ultimately be resolved while they are still in the medium, i.e. t d < L. Therefore, jet constituents produced in the phase space delimited by t f t d L, will be resolved by the medium and are affected by quenching. The rest of the phase space stays unaffected. These two effects (vacuum fragmentation and medium resolution) are captured by the collimator function C(p T , R) [29], which is a function of the jet and medium scales. It takes into account the additional energy loss of resolved vacuum-like emissions in the medium (see also Ref. [46] for an application to heavy-quark jets). The total quenching of the jet is therefore given as a product of the quenching of the total charge of the jet and the collimator, that is where i = q, g indicates the dependence on the color charge. This is what we refer to as the fully resummed quenching factor of a jet. The functions C i obey a set of coupled, non-linear evolution equations, see in Ref. [29]. Here, we use its linear approximation, where the quark and gluon solutions decouple. This allows to write the solution explicitly as where α s is the 1-loop running coupling, the relative transverse momentum is k ⊥ = z(1 − z)p T θ, and P gi (z) is the Altarelli-Parisi LO splitting function. The finiteness of the integrals is ensured by the phase space measure In the opposite case, p T q 0 L 2 , the angular integral is regulated by θ (q/p 3 T ) 1/4 . Finally, if θ c > R the jet is completely coherent and C(p T , R) = 1. This linearized version of the collimator function is analytically calculable, which is a big advantage in comparison to the full, non-linear version. We also tested against the full non-linear solution, which resulted in small, <10%, deviations even for big R ∼ 1 cones.
The resummed quenching weight Q >,i (p T , R) is shown in Fig. 2 with solid lines for quarks (left) and gluons (right) for different cone sizes. It results in more quenching (with the same medium parameters) compared to the single-parton quenching factor because there are more jet constituents that contribute to the total energy loss. The R-dependence is a result of two competing effects: a) the recapture of medium-induced gluons by opening the cone, and b) the opening of phase space for vacuum-like emissions that source additional energy loss. The combination of these effects balances out, leading to a very mild cone size dependence. For a full discussion of the uncertainties related to the choice of medium scales and parameters involved in the quenching, see also Ref. [34].
We would like to emphasize the flexibility of the collimator function C i , which is independent of the particulars of the model of energy loss employed in the previous sections. One could start with other models for single-parton quenching Q (0) and the resolved phase space Θ res , and then use the collimator to describe the quenching of the multiple resolved sources inside the jet. As a concrete example, we provide an alternative calculation in the context of the hybrid weak-and strong-coupling model [47] in App. D.

Quenching effects in the cumulative spectrum
The jet suppression factor, defined in Eq. (2.2), compares the jet spectra in heavy-ion collisions (medium), to that in proton-proton collisions (vacuum) at the same final p T . In this section, we turn to the discussion of other observables that could be constructed from the inclusive jet spectra measured in these two systems. One alternative is to compare the cumulative of the jet spectra (i.e., integrated above a fixed p T cut). The cumulative is less affected by the initial shape of the hard spectrum and has better statistical uncertainties. Another approach is to estimate the p T -shift necessary to match the pp and AA spectra. This is closely related to the typical amount of energy lost by a jet. Both of these procedures are straightforwardly related to the quenching factors entering the jet spectrum and will be discussed in further detail in Sec. 3.3. Now, we turn to an observable that combines the strength of the two examples given above, namely the quantile procedure.

The quantile ratio
The quantile procedure was introduced in Ref. [10] and aims to unfold the average momentum shift ε between vacuum and quenched jets. First, let us introduce the tail cumulative of the spectrum, which is a probability after dividing with the full integral. The quantile procedure compares the medium and the vacuum spectrum at equal probabilities, Σ med (p q,med T , R) = Σ vac (p q,vac T , R). In heavy-ion collisions, for a fixed p q,med T , this condition allows identifying the corresponding p q,vac T . Finally, the quantile momentum ratio is defined as Therefore, p q,vac T is the momentum of vacuum jets above which vacuum and medium jets have equal probability to be produced.
For a quick estimate, let us assume a steeply falling spectrum with a fixed power n = const, and neglect the R-dependence of the quenching. The tail cumulative crosssections in vacuum and in medium (see Eq. (2.1)-(2.2)) are simply This results in the quantile momentum ratio (3.4) To get the feeling for this quantity, it is instructive to consider a few simplified scenarios for the quenching factor Q(n/p T ). First, for a constant quenching factor Q(n/p T ) = Q 0 , the quantile ratio is a trivial function of the quenching factor Q med = Q 1/(n−1) 0 . Next, we will consider the single-parton quenching factor obtained by using the soft limit (ω ω c in Eq. (2.20)) of the BDMPS-Z spectrum, which is derived in App. B. The interplay between the jet cone and the broadening introduces a characteristic energy scale ω R = √q L/R which defines two regimes that we discuss below: • For p T nω R , it scales parametrically as Q (0) > (p T ) = exp[− ω 1 n/p T ] where ω 1 = 8ᾱ 2 πω c is a characteristic energy scale of the medium. In this case the quantile ratio becomes √ 2ω c ω R . This p T dependence is similar to that of medium-induced single hard scattering, cf. Eq. (2.28), and elastic drag, cf. Eq. (4.6). The quantile reads where y = nω 2 /p T is the scaling variable and γ(s, x) = Γ(s) − Γ(s, x) is the lower incomplete gamma function.
Both the quantile ratio Q med , and the jet suppression factor R med depends identically on a dimensionless ratio of a medium scale over the jet transverse momentum. Strikingly, the main difference resides in the n dependence. It turns out that the relation between quantile and quenching factor Q med Q 1/(n−1) > holds approximately also for p T -dependent quenching factors-at least for the case of fixed n. In particular, given that ln Q nω R , we should expect a reduced sensitivity of the quantile to the details of the initial spectrum at high-p T , i.e. ln Q med ∝ 1 − ω 2 /p T .
Here, we have mostly focused on the contribution from the out-of-cone, soft radiation spectrum to quenching. However, both quenching by hard emissions, see Eq. (2.28) and discussion below, and elastic energy loss, see App. D.1, behave in a similar fashion. For a single parton species, we should therefore expect to see a universal behavior, independent of the hard spectrum of the quenching at high-p T .
We study the onset of the independence of the power-index of the hard spectrum n in Fig. 3. In these plots, we have computed the single-parton quenching factor Q (0) > for a single parton species, in this case quarks, using as input the full radiative spectrum from Eq. (2.24), which encompass both soft and hard medium-induced emissions (dashed lines). We also plot the resummed quenching factor Q > , given in Eq. (2.31), (solid lines). On the left in Fig. 3, we plot quenching factor R med (p T ) for a wide range of constant n values, which only makes sense for pure samples of quark and gluon jets. 1 For realistic situations, e.g. dijet or boson-jet events, one has first to add up these contributions to the total spectrum (cf. Eq. (4.2)), before computing the cumulative. In this case, the scaling features of the quantile ratio are not necessarily as transparent. We will discuss these issues in more detail in Sec. 4. In order to make contact with the main objective of this paper, namely to provide predictions for the quantile in heavy-ion collisions at the LHC, we round off this section by studying the R-dependence of the quantile for pure quark and gluon jets in Fig. 4 for fixed n = 5 and medium parameters given in Tab. 1. As before, the dashed lines correspond to using Eq. (2.30) for the quenching factor, which assumes that the whole jet is quenched as a single parton, while the solid lines employ the resummed quenching factor Eq. (2.31). The single parton quenching trivially results in less modification and thus smaller momentum shift and quantile ratio for the same medium parameters. Generally, the qualitative features follows the naive expectation Q med Q 1/(n−1) > , cf. Fig. 2. The R-dependence is analogous to our previous discussion; less quenching results in a quenching factor closer to 1, and thus a smaller difference between the quantile momenta.

Statistical advantage of the cumulative spectrum
Using the tail-cumulative distribution has certain advantages in statistical analysis. For a steeply falling distribution f (p T ) ∼ p −n T (like the hadronic or jet spectrum), there are fewer and fewer entries in the higher p T bins. This results in a rapid increase of the relative statistical uncertainty with p T such as ∆f /f (p T ) ∼ p n/2 T (N ∆p T ) −1/2 , where N is the total number of hits, and ∆p T is the size of a bin. By using The scaling of the "gluon quantile" follows the same trends as for the quarks in

tion, the sum of the higher bins results in less uncertainty
. This can be further slowed down using ∆p T (p T ) ∼ p T logarithmic binning, resulting in ∆f /f ∼ p It is also true that the cumulative is equivalent to the original distribution, therefore for a given set of statistical samples the tail-cumulative distribution could be advantageous. 2 We would like to note, however, that one also has to consider systematic uncertainties. Some of them cancel in ratio observables, such as the nuclear modification factor R AA . This cancellation is less trivial in the cumulative case.

Other types of observables based on the spectrum
So far we have discussed the nuclear modification factor and the quantile procedure. But other observables related to the jet spectrum could also be defined (see also in Refs. [10,27,48]). We demonstrate the relation between these observables within the quenching weight formalism and show how these observables are related to the quenching factor Q(n/p T ) (or R med ) and the quantile ratio Q med below.
• Pseudo-quantile [10]: is a version of the quantile procedure which matches directly the spectrum instead of the cumulative, dσ med /dp q,med T ≡ dσ vac /dp q,vac T . The condition relates the two momenta, i.e.p q,vac T (p q,med T ). We obtain then, Because of the momentum ratio, this observable has similar n dependence to the quantile ratio, and therefore it is more robust against the initial spectrum. The statistical uncertainty, however, is similar to the R med since bins are not summed. It is equivalent with the momentum shift parameter.
• Momentum shift [27,48]: is another interpretation of the pseudo-quantile, defined by demanding dσ med (p T ) dp T ≡ dσ vac (p T + S(p T )) dp T . (3.8) The trivial connection between the definitions is p vac Therefore, the pseudo-quantile and the momentum shift are equivalent. The spectrum shift can also be expressed with the quenching factor, expressing the connection between the pseudo-quantile Q med and the shift parameter S(p T ).
A different definition of the momentum shift parameter was used by the PHENIX collaboration S loss (p T ) [48]. Their definition, however, assumes n = const and postulates S loss (p T ) = S 0 p T . The S(p T ) and Q(p T ) are more general. We found to be necessary to consider changing power n(p T ) in the spectrum, see Sec. 4 and App. A.
• Cumulative-R med (or pseudo-ratio [10]): is similar to the R med but uses the ratio of the cumulative spectrum instead, where the cumulative of the spectrum is defined in Eq. (3.1), and where we have suppressed the R dependence for now. The integral reduces bias effects from the initial spectrum and improves the statistics as we showed in Sec. 3.2. For a p Tindependent quenching weight, as in the single-parton P(ε) in the soft BDMPS-Z limit, the cumulative R med can be written as A deviation from this expectation indicates a p T dependence of the quenching weight. Furthermore, for a spectrum with constant n, where Σ vac (p T + ε)/Σ vac (p T ) = (1 + ε/p T ) 1−n , we simply get that R med = Q((n − 1)/p T ). Finally, a trivial connection between the R med and the quantile ratio for constant n is R med (p T ) ≈ Q med (p T ) n−1 .
In conclusion, we demonstrated the relationship between the different observables one can construct from the jet spectrum. It turns out that the quantile procedure is particularly appealing due to the reduced sensitivity to the hard spectrum and because of the improvement of the statistical uncertainties.

Predictions
In the previous sections, we introduced the quenching weight formalism and showed the way one can construct the quenched jet spectrum and its cumulative. In this section, we extend our formalism to provide more reliable predictions for measurements in heavy-ion collisions.

Including the realistic spectrum and elastic energy loss
A realistic calculation of quenching effects has to include the partonic cross-section for jet production. This is evaluated at the hard scale of the collision, Q hard ∼ p T . The jet spectrum at a given cone size is computed then by including a DGLAP evolution to the scale Q jet ∼ p T R [5,49,50]. Alternatively, the partonic cross-section that results in a jet with a given cone can be parameterized using a Monte Carlo event generator in which the partonic cross-section and the parton shower are matched. In our study, we generate events with Pythia8 [51] to fit the p T -dependence of the spectrum with reconstruction parameter R, for quark-, and gluon-initiated jets. The large-angle DGLAP radiation results in additional R-dependence of the R AA , through the recapture of vacuum radiation. The spectrum parameterization, therefore, includes the vacuum radiation recaptured by the cone, resulting in an additional R dependence [52], see App. C and Fig. 10, in particular. Currently, we restrict our study to inclusive jets in dijet samples, generated in pp, and 0-10% central PbPb collision at √ s N N = 5.02 TeV, with p T = 20 − 1000 GeV and |η| < 2.8, similar to the kinematics used by ATLAS [53]. 3 The details of the event generation (excluding ISR, MPI and including the effect of nuclear PDFs in Pb, jet selection, and quark/gluon flavor assignment procedure) are described in App. C. We also describe the proposed functional form, following Ref. [32], to fit the spectral indices n q (p T , R) and n g (p T , R), that automatically parameterizes the (p T , R) dependence of the quark-gluon fraction.
The extracted spectra dσ pp/AA i (p T , R)/dp T , (i = q, g) were included to calculate the nuclear modification factor R AA , defined as R AA (p T ) = dN AA (p T , R)/dp T N coll dσ pp (p T , R)/dp T , (4.1) where N coll gives the number of collisions in the nuclear overlap at a given impact parameter and we identify dN AA /N coll = dσ AA . The main difference between this ratio and the previously defined R med , defined in Eq. (2.2), is the addition of both quark and gluon jets with their respective quenching factors and the nPDFs. As a result, in our framework, we get where we have explicitly written out the dependence of the quenching factor on the spectral index n AA i (p T , R). Similarly, the cumulative spectrum in vacuum and medium are given by The quantile is finally defined as, where the two momenta are determined from the condition Σ pp (p q,pp T , R) = Σ AA (p q,AA T , R). In Sec. 2, we focused our discussion on the quenching effects emerging from mediuminduced radiation and broadening. However, for realistic predictions we should also include quenching from elastic scattering. Elastic energy loss is described by the transport coefficientê, which is related toq through Einstein's fluctuation-dissipation relationê g =q/(4T ) for gluons andê q =ê g C F /N c for quarks [54,55]. Here, T = T 0 is the local temperature of the plasma. We model the single-particle energy loss distribution simply as P(ε) = δ(ε−êL), where the flavor index is suppressed. Assuming that the energy lost in elastic processes thermalize instantaneously, we also build in the possibility to recover part of this energy through the phenomenological parameter R rec , see Eq. (2.30). This finally results in a single-parton quenching factor from elastic energy loss, given by (4.6) The criteria for resolving the partons in the jet are based on geometry and are therefore assumed to be identical for elastic and radiative processes. Therefore, the complete singleparticle quenching factor Q

Numerical results for dijet events
The single-inclusive jet R AA , generated from a sample of dijet events, for a set of cone sizes, 0.2 < R < 1, is shown in the left panel of Fig. 5 (solid curves). In the current work, the medium is treated as a static brick with fixedq 0 and length L, see Tab. 1 for details that is generally a good approximation even for expanding media, see Ref. [56]. There is a notable change of the curves at high p T due to the inclusion of nPDFs (see also in Fig. 10 in App. C for only the nPDF effects). The overall R dependence is very modest and will be discussed in more detail shortly.  [57] and CMS [58] are also shown.
The parameters of the calculation are tuned to the measured inclusive jet data from ATLAS [53] at p T 100 GeV and R = 0.4 with |η| < 2.8, cf. Tab. 1, resulting in good agreement between data and theory for the whole p T range. The measured inclusive jet R AA from ALICE [57] for R = 0.2 and R = 0.4 are also shown in Fig. 5, where the rapidity range for the jet selection |η| < 0.5 is slightly different. The recent CMS [58] results are also shown in Fig. 5 for various R, where the rapidity is |η| < 2. We would like to note, there is a disagreement between the ATLAS and CMS data that was not pointed out in the CMS' latter publication. We would like to also note that the magnitude of R AA (p T ) in any BDMPS-Z type of calculation is mostly sensitive to the combination ω s ∼ α 2 medq 0 L 2 , as observed in Refs. [27,59]. The slope of the R AA (p T ) is quite robust to changes in the parameters. 4 Given our simplified modeling of the medium, we do not attempt to reproduce the centrality dependence of the jet R AA at high-p T which will be left to future work, see also in Ref. [34].
Having constrained the medium parameters with R AA , we now turn to the predictions for the quantile momentum ratio which is showed in the right panel of Fig. 5. There is no drastic change due to the nPDFs in comparison to the R AA at high-p T , demonstrating the robustness of the quantile procedure against modifications in the partonic cross-section. The shape of the curves and even the R dependence is very similar to the R AA and is well captured by the approximate relation Q AA ∼ R 1/(n−1) AA . To date, there are no experimental measurements of the quantile ratio.
Our formalism accounts for the cone size dependence of jet quenching through competing effects related to the early vacuum shower and medium-induced elastic and radiative processes. On the one hand, by opening the cone, one captures more of the particles that are affected by medium interactions and thus recover the lost energy. This is manifested as a suppression of the spectrum dI > /dω at ω > √q L/R in Fig. 1, meaning it is less probable to lose energy. For R ∼ R rec all quasi-thermalized modes, both in the radiative spectrum (ω < ω s ) and due to elastic collisions, are recaptured within the jet cone. On the other hand, opening the reconstructed jet cone results in more phase space for vacuum fragmentation at an early stage of the evolution. This leads to a higher multiplicity of vacuum-like emissions and, therefore, more sources for energy loss.
The overall effect is a relative cancellation of the R dependence (see Fig. 5). Figure 6 shows the R AA (left) and quantile ratio (right) at a given R divided by the R = 0.2 results. As an illustration, here the dashed curves are the results obtained by using the singleparton quenching factors Q (0) > , which amounts of treating the whole jet as a completely coherent single parton that is not resolved by medium interactions. Their R-dependence reflects directly the effect of recapturing energy at large angles without sourcing more energy loss though vacuum fragmentation. The full curves are for the full quenching factors which leads to a more complicated R dependence where, several effects contribute, such as vacuum fragmentation and its recapture, color coherence effects and the thermalization of the medium-induced emissions. This ratio was measured by ALICE [57] (for |η| < 0.5), and CMS [58] (for |η| < 2), showing great agreement with our model. Reference [58] compares many theory prediction and we can say, currently, our model has the best agreement. The main differences between the "bare", and "resummed" quenching factors is at lower p T ∼ 50 − 100 GeV. As we mentioned, in our model the R AA is mostly sensitive to the ω s ∼ α 2 medq 0 L 2 combination of the parameters. A precise measurement on the R-dependence would help to constraint more parameters. The right side of Fig. 6 shows the R-dependence of the quantile ratio. It is much less sensitive to the jet cone.

Differences between dijets and Z+jets
At first glance, there should be no difference in the quenching mechanism in dijet and boson+jet events. 5 However, their R AA are not expected to be the same. In this subsection we explain why, and what are the consequences for quark and gluon classification. The Z+jet process is used as an illustration. We also show how can one use the cumulative distribution for quark-gluon discrimination to perform better statistics.

Difference in the R AA
In recent years, a lot of effort has been put to measure and understand the boson+jet processes both in pp and AA collisions. Firstly, it is a favorable process for quark-, and gluon-jet discrimination, because it provides a natural definition of the initial jet flavor. Secondly, measuring the boson momentum, one gets a label on the initial momentum of the recoiling jet. This is especially advantageous in heavy-ion collisions, where the quenching of bosons is suppressed, and one, therefore, gains knowledge about the jet before final-state interactions with the medium.
Here, we focus on Z+jet process, but the arguments are valid for other boson+jet processes. We generated the Z+jet events with Pythia and identified jets using the same cuts as in the dijet study. Further details on the event generation and parameterization are summarized in App. C. We include the corresponding spectrum using Eq. (4.2) and the result is shown in the right panel of Fig. 7. Dijet and boson+jet processes have different spectra and thus their bias is different on the energy loss distribution. As discussed before, this appears through the spectrum power n(p T ) in the quenching factor. For dijets n ≈ 6, for Z+jets n ≈ 4, therefore based on the approximate R AA ∼ exp(−2 πᾱ 2 ω c n/p T ), the bigger the power, the stronger is the quenching. To stress this, the dijet-, and Z+jet-R AA ratio is also shown on the lower panel. The ratio of two exponential remains to be exponential that looks linear on the semi-log scale. Moreover, dijets are gluon dominated, while Z+jets are quark dominated, see in Fig. 10 in App. C. Therefore we expect Z+jets to have smaller quenching, which is in agreement with Fig. 7. The p T dependence of the power is also different, resulting in a different slope of R AA as high-p T . Note also that the Z+jet quark contribution is relatively constant in p T in contrast to the increasing quark contribution in dijets. Since the quenching roughly scales with the color charge (gluons have more quenching), at higher p T , Z+jets has stronger quenching, which is in agreement with the smaller slope in the plot. Although the Z+jet and dijet spectra are different, the relative R dependence is similar.
The quantile ratio is also shown in the right panel of Fig. 7. The difference is smaller between dijets and Z+jets than for the nuclear modification factor, pointing out the reduced sensitivity to the details of the initial spectrum. Using the cumulative spectrum, and the ratio of momenta, one gets much less sensitive to the initial shape of the spectrum. The degree of scaling the initial spectrum of the quantile is nevertheless not as ideal as for the results in Fig. 3. This can be traced to the fact that the single inclusive jets in neither dijet nor Z+jet events are pure samples of quark-, or gluon-initiated jets. The different admixture of parton species, as well as the different level of quenching of the two both, contribute to delaying the onset of scaling effects up to higher transverse momenta. One also has to point out that Fig. 3 was obtained by assuming n =const, while for realistic calculations n q and n g are complicated functions of p T and differ significantly for dijet and Z+jet events, see Fig. 10.
All in all, our results are also qualitatively similar to the ones observed in Ref. [10]. Similarly to dijets, the R-dependence is very similar to the Z+jet R AA and is to a great degree captured by the approximate relation Q AA ∼ R 1/(n−1) AA .

Improving quark-gluon discrimination
As mentioned above, jets recoiling from a boson (γ or Z/W ) is the preferred process for quark-gluon (q/g) discrimination, because the hard scattering naturally defines the initial flavor of the jet at leading order in perturbation theory. Most q/g discriminators apply cut on jet substructure observables, e.g., jet mass or soft-drop multiplicity, to classify jets, see in Refs. [1,60]. While the best performance is achieved by machine learning algorithms, these nonetheless rely on training sets resulting in model dependence. In pp, however, the main description of jets is predominantly perturbative or it can be made perturbative with grooming. Model dependence, therefore, is not as crucial as for heavy-ions, where the description is not yet unique and where non-perturbative effects are more prominent.
The recently introduced topic modeling [15, 61, 62] is a data-driven method that is largely model independent, and its outstanding performance was demonstrated for event generator samples both in pp and AA. For this reason, it is also applicable to heavy-ions. There are some caveats of the classifier; (i) it works only with certain observables for which quarks and gluons are mutually irreducible (usually counting-type observables [63,64]), (ii) it is heavily limited by statistical uncertainty, (iii) the performance is limited by the cuts on the phase space, and (iv) it works on statistical samples.
Topic modeling aims to un-mix the sample probability distributions (e.g., dijet and boson+jet samples) to a common basis (quark and gluon distribution). Consider, for example, where p(x) is a probability density of some observable for which Eq. (4.8) is true (mutual reducibility), and f Zj/2j q and f Zj/2j g are the weight factors. To un-mix, one uses the fact that phase-space of x exists, where either the dijet or boson+jet (and thus the quark or gluon) dominates the distribution. This usually happens on the domain border of the observable (e.g., small/big jet multiplicity) [63,64]. With this, called anchor, bin one can statistically decouple the basis using p Zj (x) − κp 2j (x) ≥ 0, where the reducibility factor is where i, j = Zj, 2j. Equation (4.9) is true if inf x p q (x)/p g (x) = 0, however, for real data it has a finite minimum. The minimum is typically on the edge of the x distribution, and thus the extraction of κ is limited by the statistical uncertainty of this corner bin. Because of the linearity, one could integrate both sides of Eq. (4.8), and rewrite Eq. (4.9) using the cumulative distribution of p(x),κ . (4.12) Our cumulative method improves the statistical uncertainty by definition (see Sec. 3.2) that can be trivially tested with arbitrary combined distributions. We would like to note the cumulative in this subsection refers to the p(x) distribution, and has nothing to do with the cumulative of the jet spectrum. Unfortunately, the jet spectrum is not mutually irreducible. One can see this from Eq. (4.2), where the quenching factor depends on the quark/gluon spectrum through their n indices. However, the medium modified quark-gluon ratio of the jet spectrum is important for any quark-gluon discriminator, and thus we provide it in Fig. 8, for both dijets and Z+jets. We observe that at high-p T the ratio barely changes, however, at lower p T quarks start to dominate. This is the result of the stronger quenching of gluons, effectively suppressing them in the samples. For future quark-gluon classification, this suppression effect indicates that the Z+jet quark-gluon ratio will be less different from the dijet ratio, making the separation, unfortunately, harder, in line with what was observed in Ref. [65].

Conclusion
In heavy-ion collisions, the steeply falling jet spectrum, convolved with the probability for quenching, biases the measured jet observables. In this paper, we studied the origin of this bias and its presence in recently suggested observables, e.g., the quantile ratio and its comparison of single-inclusive jet spectra in dijet and boson+jet events in heavy-ion collisions.
We developed a novel analytic framework based on quenching factors to construct the jet spectrum in heavy-ion collisions. Starting from single parton energy loss, we showed the necessity to include the full medium-induced radiation spectrum, including both multiplesoft (described by the BDMPS-Z spectrum) and single-hard (included in the GLV spectrum) scattering regimes. By keeping track of the jet cone size, the energy is only lost if the emissions propagate out of the cone. We also account for the quenching of multiple jet partons resolved by the medium interactions through the collimator function, see Eq. (2.31). Therefore, opening the cone, vacuum-like jet fragmentation sources more partons to quench, resulting in a relative cancellation of the cone-size dependence. It is also important to use realistic jet spectrum for predictions by including the p T dependent spectrum power n(p T ), p T -dependent quark-gluon jet ratio, and nPDFs. We revealed the importance of elastic energy loss, and also included it in our quenching framework for jets.
Different observables have been introduced to study the jet spectrum in heavy-ion collisions, e.g., spectrum shift, cumulative-R AA , quantile or pseudo-quantile procedure. We showed their relation to the nuclear modification factor, to each other, and elucidated their dependence on the jet momentum p T , the spectrum power n, and the jet cone size R. The cumulative-based observables reduce statistical uncertainty, and the momentum-ratio-based ones change the spectrum power dependence. The quantile momentum ratio maximizes both of these advantages, and we provided predictions for its cone-size dependence for the first time. Our approximate formula for the quantile momentum ratio Q AA ∼ R 1/(n−1) AA captures the rough properties of the observable: it is similar to R AA , with reduced spectrum power index n.
Finally, we demonstrated the bias effect by comparing quenched dijet and Z+jet spectra. For a pure sample of quark or gluon jets, the quenching factors for single-partons scale like where ω R = √q L/R. This qualitative behavior is also numerically confirmed for the resummed quenching factor Q > . As a direct manifestation of the bias effect, i.e. the ndependence in Eq. (5.1), the jet spectrum in Z+jet events, which is less steep, results in a smaller nuclear modification factor R AA than for dijet events, see Fig. 3 (left). In other words, the steeper the spectrum the stronger the effect of quenching.
The quantile ratio Q AA is much more resilient to details of the hard spectrum, both the actual values of n and of the relative admixture of quarks and gluons. For pure samples of quark-, and gluon-jets, the behavior in Eq. (5.1) predicts an almost ideal scaling at p T nω R , which was largely confirmed in Fig. 3 for n =const. For realistic samples of jets in dijet and Z+jet events, however, the universal behavior of the quantile ratio can only be expected to be approximate, see Fig. 7 (right).
The dijet and Z+jet events are also useful for quark and gluon discrimination. Following our cumulative spectrum experiences, we improved quark-, and gluon-jet discrimination based on topics modeling in general. However, due to the n-dependence of the quenching factors, quarks and gluons in the jet spectrum are not mutually irreducible. We will return to the challenging of quark/gluon discrimination in heavy-ion collisions in future work.
Many of the assumptions underlying the concrete realization behind our numerical results in Sec. 2 are already implemented in varying degrees in various phenomenological Monte Carlo models. Our framework, therefore, provides theoretical tools to organize the effects of quark/gluon contributions, jet fragmentation, and finally, the details of medium interactions. In the BDMPS-Z framework, these are, to a first approximation, all controlled by a single transport coefficientq. also supported by the MCnetITN3 H2020 Marie Curie Initial Training Network, contract 722104, and wishes to thank the Institut de Physique Theorique (IPhT) and Gregory Soyez for the hospitality.
A Corrections to the Laplace transformation and to the running power Our objective in this paper, is to compute the ratio of medium to vacuum inclusive jet spectra, which can be written as (see Eqs. (2.1)-(2.2)) where we omit the R dependence for now, and assume that n = const. In Eq. (2.3), we took the (1 + ε/p T ) −n ≈ exp(−nε/p T ) approximation. We can, however, easily include corrections to this by noting that where Q(ν) ≡ ∞ 0 dε P(ε)e −νε is the Laplace transform of the energy loss distribution and ν = n/p T . In this appendix, we investigate the impact of these higher-order corrections for a concrete example that can be solved analytically, namely the energy loss distribution obtained in the strictly soft limit of the BDMPS-Z spectrum. It is given by which only depends on the energy scale ω s and is properly normalized. In this case, its Laplace transform is Q = exp(−2 √ πω s ν). We can, in fact find any of the terms in Eq. (A. 2) by noticing that where K m (x) is the modified Bessel function of the second kind and I 0 = Q(ν). We can therefore write where, c 0 = 1, c 1 = 0, c 2 = ν 2 /(2n), c 3 = ν 3 /(3n 2 ) and so forth, by assuming R med (p T ) to be analytic function. On the left of Fig. 9, we study the corrections by comparing to the exact value from Eq. (A.1). The parameters we use are ω s = 5 GeV and ω s = 10 GeV and n = 5. For the realistic choice of ω s 5 GeV, the leading behaviour is already of the order of O(10 −2 ), even at low p T ∼ 100 GeV. Finally, the effect of the corrections from Eq. (A.2) for the realistic quenching weight employed in the main body of the paper is showcased in Fig. 9 (right) up to the leading correction O(∂ 2 ν ). The sign of the correction depends on the concavity/convexity of  the quenching weight (note the second derivative in p T in Eq. (A.2)). For our choice of parameters, the correction is tiny < O(10 −3 ), and can safely be neglected.
The conclusions from above also hold if we consider a p T dependent power n(p T ) of the hard spectrum. One can trivially replace n → n(p T ) in the quenching weight capturing most of the effects of the running power. The correction to this simple replacement is B p T -scaling of the out-of-cone emission quenching factor Consider the BDMPS-Z spectrum in the soft limit. In the short formation-time approximation the emission and subsequent broadening of a soft gluon factorizes, and we can write where P(k ⊥ , L − t) is the probability for a particle emitted at time t to acquire transverse momentum k ⊥ = ωθ through elastic scattering up to the end of the medium L. In the Gaussian approximation, it reads For this simplified ansatz, the out-of-cone spectrum reads where ω R = √q L/R and the broadening factor B is given by where y = ω 2 /ω 2 R . Therefore, Eq. (B.3) falls rapidly for ω > ω R , and thus we used ω R in the matching definition in Eq. (2.13). The single-parton quenching factor reads then where we changed variables to x = ων and x R ≡ ω R ν. We will solve this integral in two limiting cases, namely x R 1 and x R 1. In the former case, which corresponds to the high-p T regime where p T ω R n, we can expand the terms in the bracket and find ln Q In the opposite limit, for p T ω R n, we can extend the upper integration limit to infinity, to obtain which is independent of the jet cone.
C Generating and parametrizing the jet spectrum As we mentioned in Sec. 2, the partonic cross-section to produce a jet with a given p T and R is perturbatively calculable in the collinear factorization up to high precision [5,6]. In our work, instead, we extract the spectrum using the Pythia8.235 event generator [51].
To generate dijet events we used default settings and tunes with HardQCD:All both in pp and in 0-10% PbPb collision at 5.02 TeV. This results in LO 2 → 2 matrix elements. The nPDF was EPS09LO which has a relatively important effect on the R AA , see in Fig. 10. The ISR, MPI, and hadronization were turned off to focus on final state radiation only. We reconstructed jets using anti-k t algorithm with FastJet3 [66] for R = 0−1, p T,jet = 10−1000 GeV and |η jet | < 2.8, similar to the kinematic cuts of ATLAS. 6 To label the flavors of the jets, we compared them to the outgoing partons from the hard scattering, and we kept the closest in angle if it was less than 2R. We only associated one jet (the hardest) with an initiator, and thus we only kept the two hardest associated jets. This selection detail becomes important for small R jets, where more than one jet can be reconstructed for one initiator. We kept those events in which there no jet passing the criteria, which are important in the proper R AA ratio (before quenching). We parametrized the spectrum where i is the flavor of the initiator of the jet and {p 0 , c n } are (p T , R, i, pp/AA) dependent fitting parameters. We kept terms up to n = 3, achieving < 3% relative deviation. The resulted parametrization is showed in Fig. 10. On the left, there is the quark contribution, which increases with p T . The cone size dependence shows, gluons are emitted at larger angles even in the vacuum. On the right the R AA is shown, resulted by the nPDF (no quenching on the plot). The inclusive jet spectrum would be similar to the dijets keeping not only the two hardest jets, however, the flavor assignment would be less trivial especially for smaller cone sizes, therefore we preferred to use the dijet samples. 7 For the Z+jet samples, we followed the logic from previous. We used however the WeakBosonAndParton:qg2gmZq and WeakBosonAndParton:qqbar2gmZg processes and we forced the Z-boson to decay invisibly. In this case, we kept only the hardest recoiling jet if it fulfilled the same criteria as before. In the Z+jet case, the spectrum power is smaller in comparison to dijets n Zj ≈ 4 < n 2j ≈ 6. The quark and gluon ratio is different, in comparison to dijets, but the cone size dependence is similar, see on the left in Fig. 10.
In both cases, opening the cone, one captures more gluons meaning, they are radiated at larger angles. The R AA resulted from the nPDF is also shown in Fig. 10, relatively similar for both dijet and Z+jet processes.

D Other Models
As we mentioned in Sec. 2 the quenching weight and collimator formalism are independent of the energy loss description. To illustrate this we show how to use this formalism to describe quenching of elastic scattering and within the strong coupling approximation.

D.1 Momentum shift and elastic energy loss
The simplest example for the jet energy loss probability P(ε) is to take a momentum shift of vacuum jets P(ε) = δ(ε − S(p T )), resulting equivalent definition to the old fashioned momentum shift parameter dσ med (p T )/dp T ≡ dσ vac (p T + S(p T ))/dp T [27,48]. By expanding the steeply falling spectrum in S(p T ) p T , dσ med dp T (p T ) = ∞ n ε n d n dp n T dσ vac dp T ≈ dσ vac dp T (p T + ε ) , (D.1) where we used ε i ≈ ε i . The shift parameter is roughly the mean energy loss S(p T ) ≈ ε = dε ε P(ε) ≡ ∆E . (D. 2) The quenching factor in the simple power-law case is The elastic scattering is approximated by a constant momentum shift P el,i (ε) = δ(ε − e i L), whereê = −d E /dt ≈ C iq0 /(4N c T 0 ), and T 0 ≈ 0.45 GeV [55,68]. This can be translated to an additional quenching weight compared to the form Eq. (2.30), where we included some energy recapture through the second term with R rec = π/2. Therefore in total Q >,el . This factor runs slower with p T than the BDMPS-Z, similar to the "N = 1" (or GLV) spectrum. We can also include the elastic energy loss of each jet constituent, by using the collimator function form Eq. (2.32), Q >,tot = Q >,rad Q >,el . Fig. 11 shows Eq. (D.5) with dashed lines and with the collimator with full lines, indicating the importance of elastic scattering in the overall quenching. Therefore we included this effect in Sec. 4.

D.2 Hybrid weak-, strong-coupling
In the strong coupling approximation, the energy loss is assumed to be described by a semi-classical string falling inside a black hole horizon [69]. This model does not use dI/dω because there are no emitted gluons, and thus the energy loss is directly connected with P (0) (ε). The average lost energy of a single parton traversing through an L sized, T temperature strongly coupled medium is where E is the initial energy, x s = E 1/3 /(2κ sc T 4/3 ) is the stopping length and κ sc = 1.05 g 1/3 N 1/6 c . Using the definition of momentum shift from Eq. (D.2), the single parton quenching is estimated by where the p T dependence is presented in x s . For high p T , Q (0) (p T ) ≈ 1−32nκ 3 sc L 3 T 4 /(3πp T ) the same p T dependence obtained from GLV (see Eq. (2.28) and below). We did not include here the broadening in and out of the cone due to the lack of particles, however, one could include the linearized hydro response through R rec (η, φ). We include multi parton quenching and thus the R-dependence by using the collimator from Eq. (2.32). We used Eq. (D.7) as Q (0) (p T ), and for Θ res we used the resolution condition used in the hybrid model [70], t f < t d = (θµ) −1 < L, where µ = πT /2 is the IR screening scale used in the hybrid model (the Debye mass would be gT ). Eq. (D.7) is shown in Fig. 12 with the Pythia spectrum, resulting similar quenching to the GLV assumption. The curves in Fig. 12 similar to the results in [71] without the medium response. The R-dependence is simplified in our case, because we neglected the medium response. We also estimate the quantile ratio using Eq. (D.3). The parameters are T = 0.27 GeV, L = 4 fm and κ sc = 0.4.