Medium response in JEWEL and its impact on jet shape observables in heavy ion collisions

Realistic modeling of medium-jet interactions in heavy ion collisions is becoming increasingly important to successfully predict jet structure and shape observables. In JEWEL, all partons belonging to the parton showers initiated by hard scattered partons undergo collisions with thermal partons from the medium, leading to both elastic and radiative energy loss. The recoiling medium partons carry away energy and momentum from the jet. Since the thermal component of these recoils' momenta is part of the soft background activity, comparison with data requires the implementation of a subtraction procedure. We present two independent procedures through which background subtraction can be performed and discuss the impact of the medium recoil on jet shape observables. Keeping track of the medium response significantly improves the JEWEL description of jet shape measurements.


Introduction
Hard probes in heavy ion collisions that traverse the hot and dense medium, often referred to as the quark gluon plasma (QGP), undergo significant energy loss. One of the avenues to study this energy loss is by the use of reconstructed jets. Results from experiments at the LHC, have shown consistent findings regarding the degree to which the jet spectra is quenched [1][2][3] and the recovery of the "lost" energy at large angles from the jet [4]. Simultaneously, medium induced modification to the jet structure have also been observed in these experiments. These results show an increase in the multiplicity of low p ⊥ jet constituents as seen in the jet fragmentation function [5][6][7], and that these constituents are found in the periphery of the jet [8]. An increase of relatively soft particle multiplicity is also observed at large angles from the jet axis [9]. At the same time the hard component of jets gets narrower, as observed in a slight reduction of the jet girth [10]. Recent results have also probed the so called groomed sub-jet shared momentum fraction [11], where Pb+Pb jets are split more asymmetrically than p+p jets in the lower p ⊥ bins, with the effect disappearing at higher p ⊥ . The invariant charged jet mass [12] was also recently measured in different jet p ⊥ bins and found to be largely unmodified as compared to p+p. Together these results point to a picture of modifications to the jet structure in heavy ion collisions that involve a narrowing of the jet core along with a broadening of the jet at its periphery.
It has been observed [13][14][15] that the narrowing of the jet core, accompanied by a hardening of the fragmentation function at large z, is due to the fact that broader jets with a softer fragmentation pattern are more susceptible to energy loss and are thus more likely to fail analysis cuts and disappear from the sample. The increase of soft particle production at relatively large angles from the jet axis is commonly attributed to the medium's response to the energy and momentum deposited by the jet [15][16][17], although other interpretations exist [18].
Since the latest high luminosity runs at the LHC, the focus in heavy ion jet studies has moved to detailed characterisation of jet shapes and intra-jet observables to highlight the modification to the internal structure of the jet. The impact that medium response has on such jet shape or jet sub-structure observables thus offers a possibility to observe the thermalisation of energy and momentum deposited by a jet. In this context, the advantage of jet observables is that the additional soft activity can be clearly identified inside a jet, while globally it is much more difficult to separate it from the much larger, soft, fluctuating background.
Theoretically, medium response is described using different frameworks. At strong coupling it has been known for some time that energy and momentum lost by a quark is transferred into hydrodynamic modes and leads to a Mach cone and wake [19,20] (for reviews see e.g. [21,22]). Medium response at weak coupling has also been studied in a hydrodynamic framework [23,24], supported by the observation of fast thermalisation of soft fragments [25]. While these approaches rely on a clear separation between jet and medium degrees of freedom, transport codes treat all partons, whether they are hard or soft, on equal footing. In these calculations, soft partons that interacted with a jet, undergo further re-scattering and thereby thermalise [26][27][28]. Recently, hybrid approaches have been developed, that describe the propagation of jets in transport theory, but treat the thermal medium and its response to energy and momentum depositions by jets in hydrodynamics [15,17,29]. A different kind of hybrid approach is followed by [15], where the energy lost by a weak coupling jet is calculated at strong coupling and the medium response is again treated in hydrodynamics.
In the jet event generator Jewel [30] three options for dealing with medium response are available: (i) one can ignore it, (ii) extract a source term for a hydrodynamic treatment [31], or (iii) keep thermal partons recoiling against interactions with the jet in the event and let them hadronise together with the jet. In this paper, we explore the latter option. It is organised as follows: In section 2 we discuss the handling of medium response in the Jewel generator, followed by a description of two subtraction techniques in section 3. After a brief specification of the used Monte Carlo samples in section 4 we then characterise the subtraction methods systematically in section 5. We then discuss the effects of medium response on traditional jet quenching observables (section 6) and jet shape observables (section 7) and close with a discussion of our results in section 8.

Treatment of medium response in JEWEL
In Jewel the background medium is assumed to consist of an ensemble of partons, the phase space distribution and flavour composition of which have to be provided by an external medium model. Partons belonging to a jet may interact with these background partons through 2 → 2 scattering processes described by perturbative matrix elements, with associated gluon emission generated by the parton shower. Further details of the inner workings and Monte Carlo implementation of Jewel are available in [32,34].
As mentioned in the introduction, there are two operational modes for event generation with Jewel concerning the treatment of background partons recoiling against a scattering with the jet (so called "recoils" or "recoiling partons"). Events can be generated with or without storing the recoil information. When run without recoils, the recoiling partons do not show up in the event. In this case no medium response is considered and inclusive and inter-jet observables can be compared to (background subtracted) experimental data. So far this was the recommended mode for jet observables.
However, jet structure observables are sensitive to medium response and hence it is desirable to include these effects in Jewel by keeping the recoiling partons in the event. After the scattering these recoiling partons do not interact further with the medium and free-stream towards hadronisation. This represents a limiting case for the recoil behavior, that can be regarded as being the limit opposite the assumption of immediate thermalisation of recoil energy and momentum made by hydrodynamic frameworks. The truth is expected to be between these two extreme cases, since one would expect that these partons interact further with the medium, but do not necessarily fully thermalise.
So far the background partons could be either (anti-)quarks or gluons. For hadronisation, however, all recoiling partons are converted to gluons. It is assumed that the recoiling parton is a colour neighbour of the hard parton it interacted with. The recoil gluons are thus inserted in the strings connecting the partons forming the jet. Therefore, the hadronic final state including recoiling partons is not an incoherent superposition of jets and activity arising from recoils. At hadron level, it is impossible to assign a certain hadron to the jet or medium response.
The four-momentum of the recoiling partons has two components: the thermal momentum it had before the interaction with the jet, and the four-momentum transferred from the jet in the scattering process. Only the latter is interesting for investigating medium response, the former is part of the uncorrelated thermal background that is subtracted from the jet. As Jewel generates only the jets and not full heavy ions events, it is not possible to use the experimental background subtraction techniques for the Monte Carlo events. Instead, a dedicated procedure for removing the thermal four-momentum components from the jets when running with recoils has to be devised. Therefore, along with the recoiling partons, we are also storing the thermal four-momenta, which constitute our background 1 . These will be systematically removed from the jets during the analysis step, as detailed in the following section.

Subtraction of the thermal component
As discussed in the previous section, in order to compare predictions for jet observables with data, it is imperative to perform a background subtraction on Jewel events when running with the recoils. This is to avoid a mismatch between the prediction and data, since the jets in data have the fluctuating underlying event subtracted. In this section, we present two independent subtraction methods for Jewel, that can be employed at the analysis level 2 .

4MomSub
This method removes the thermal momenta exactly from the jet's four-momentum. In order to determine which thermal momenta should be subtracted, an additional set of neutral particles with very small energy and momenta and pointing in the direction of the thermal momenta are added to the final state particles list. These "dummy" particles are effectively the same as ghosts that FastJet [36] uses during its clustering to determine the jet area. They can get clustered into jets without affecting the jet's momentum or structure. Thermal momenta, that are matched to a dummy (in the azimuthal angle -pseudorapidity plane) inside a jet, are subtracted from the jet's momentum. The resulting four vector constitutes the subtracted jet momentum. An algorithmic implementation of the procedure is detailed below: 1. Cluster the initial jet collection from the final state particles (including dummies).
2. Compile a list of the thermal momenta (particles in the HepMC event record with status code 3).
3. For each jet, get the list of thermal momenta that have ∆R < 1 · 10 −5 with one of the jet constituents, i.e a dummy particle.
4. Sum up the four-momenta of the matched thermal momenta. This constitutes the background.
5. For each jet subtract the background four-momentum from the jet's four momentum, this provides the corrected jet collection.
This method is easily generalised to subtraction of sub-systems of jets, such as sub-jets or annuli used for the jet profile.

GridSub
This is a generic, observable independent subtraction method. A finite resolution grid (in the φ − η plane) is superimposed on the jet and its constituents. The four-momenta of particles in each cell in the grid are then vectorially summed and thermal momenta subtracted, yielding the cell four-momentum. Finally, we re-cluster the jet with the cell four-momenta as input to the jet clustering algorithm. This method does not require dummy particles. It is also possible to first discretise the entire event, subtract thermal four-momenta cell-by-cell, and then cluster jets. The algorithms for the two variants are given below.
Jet clustering before discretisation (GridSub1): 1. Cluster the initial jet collection from the final state particles.
2. Compile a list of the thermal momenta (particles in the HepMC event record with status code 3).
3. Define the grid resolution and place grid over jets.
4. Inside each cell sum the jet constituents' four-momenta and subtract the thermal fourmomenta that fall into the cell (note: no matching is required, thermal four-momenta with distance ∆R < R from the jet axis are considered 3 ), providing a single four momentum for each cell. 5. In case a cell contains more thermal momentum than jet constituents, the cells is set to have zero four-momentum. This is deemed to be the case when the (scalar) p ⊥ of the thermal component is larger than the p ⊥ of the particle component.
6. Re-cluster the jets with the cell four-momenta as input to get the final, subtracted jets.
This version is the default.
Discretisation before jet clustering (GridSub2): 1. Compile a list of the thermal momenta (particles in the HepMC event record with status code 3).
2. Define the grid resolution and place grid over the entire event.
3. Inside each cell sum the final state particles' four-momenta and subtract the thermal four-momenta that fall into the cell (note: no matching is required), providing a single four momentum for each cell. 4. In case a cell contains more thermal momentum than particle momentum, the cells is set to have zero four-momentum. This is deemed to be the case when the (scalar) p ⊥ of the thermal component is larger than the p ⊥ of the particle component.

5.
Cluster the jets with the cell four-momenta as input to get the final, subtracted jet.
Due to the finite size of the grid, it is possible to have certain cells with more thermal momentum than particle momentum, resulting in a total negative four-momentum, which in our case is set to zero before clustering. Thus, the GridSub method systematically removes less background from the jet than 4MomSub. The smearing introduced by the GridSub method will be quantified systematically in the following section.
The use of the 4MomSub method is recommended when possible, since it does not introduce finite-resolution effects and is consequently more accurate.

Limitations of the subtraction and the issue of track jets
Since the subtraction techniques introduced above subtract the thermal momenta, which are at parton level, from the hadronic final state, they only yield meaningful results for observables that are insensitive to hadronisation effects. This is the case for most infra-red safe observables based on calorimetric jets. Examples for observables that do not fall into this category are fragmentation functions and charged jet observables. In general, all cuts on the final state particles, also p ⊥ cuts, are problematic.
A few of the recent experimental results involve the use of charged or track jets [10,12,37], i.e jets reconstructed using only tracks. When the subtraction is naively applied, the techniques end up overestimating the contribution of the four-momenta to subtract. Thus, in order to compare with such experimental results, a heuristic procedure is applied. The observable of interest is calculated for full jets and re-scaled. The re-scaling between the full and the charged jet distribution is extracted from the corresponding Jewel simulation for p+p collisions. If it is larger than the resolution of the observable, it is applied to the full jet subtracted distribution in Pb+Pb. In this way an estimate of the charged jet distribution is derived. For example, a naive way of estimating the charged jet four-momentum is by re-scaling the full jet quantity with the fraction of charged particles in the jet. The charged jet mass distribution discussed in section 7 is estimated using this technique and compared with data. In other cases, for instance the jet radial moment girth (also shown in section 7), the distributions for charged and full jets are the same in p+p collisions. In this case we compute the observable for full jets in Pb+Pb collisions as well and do not apply any re-scaling.
Obviously, this method comes with an additional uncontrolled systematic uncertainty, since it is not guaranteed that the same relation between full and charged jet distribution holds in Pb+Pb and p+p.

The event sample
We generate di-jet events in the standard setup [30] at √ s NN = 2.76 TeV and √ s NN = 5.02 TeV with the simple parametrization of the background discussed in detail in [33]. The values used for the formation time and initial temperature, τ i = 0.6 fm and T i = 485 MeV for √ s NN = 2.76 TeV and τ i = 0.4 fm and T i = 590 MeV for √ s NN = 5.02 TeV, are taken from a hydrodynamic calculation [38,39]. The proton PDF set is Cteq6LL [40] and for the Pb+Pb sample the Eps09 [41] nuclear PDF set is used in addition, both are provided by Lhapdf [42]. We use the Rivet analysis framework [35] for all our studies. Jets are reconstructed using the same jet algorithm as the experiments (anti-k ⊥ [43]) from the FastJet package [36].
For our studies we use Pb+Pb samples with and without medium response and a corresponding p+p sample. In addition, we also need samples at parton level for some of the systematic studies, as they allow for a parton-by-parton separation of jet and recoil contributions.

Systematic studies
The background subtraction techniques introduced in the previous section and their effects on jets are studied henceforth in a systematic fashion.

Smearing due to finite resolution of the grid
An immediate consequence of the grid, before any subtraction is introduced, is that both the jet p ⊥ and the position of the jet in the η − φ plane are smeared. This effect is studied in Jewel with p+p events generated at hadron level to highlight the inherent behavior. All our studies of the grid are shown for a nominal grid size of 0.05 in η − φ plane, which we find to be a good compromise between resolution and under-subtraction (which is more severe for smaller cell sizes). The systematic uncertainties are estimated by varying the grid size by a factor of two and most final observables are shown to be quite insensitive to the grid size within these limits.
In each event, jets are first reconstructed from the final state hadrons. Then the event is discretised using a grid and jets are reconstructed based on the grid cells. Finally, each jet of the smaller of the two collections is matched to the one from the other set that is closest in ∆R, with the constraint that ∆R is smaller than the reconstruction radius (this is the standard CMS procedure for comparing generator level jets to jets after detector simulation). The smearing is quantified in Fig. 1 with the top panels showing the smearing in jet p ⊥ (on the left) and jet axis (on the right) as a function of the particle jet p ⊥ . The latter is broken down into the respective shifts in η and φ, which are shown in the bottom panels. The deviations are observed to be small in the p ⊥ range studied here. There is a clear trend for the grid jet p ⊥ to be larger than the corresponding particle jet p ⊥ , which is due to the fact that the effective area of the grid jets can be larger due to the discretisation. As one would expect, increasing the jet p ⊥ reduces the smearing introduced by the grid.  shown as a function of the particle jet p ⊥ . Bottom plots show the relative shifts in jet η (left) and φ (right), shown as a function of the respective particle jet η and φ. (Note that the log scales on the z axis span six orders of magnitude.) −2.4 < particle jet η < −1.44 −1.44 < particle jet η < 0.48 −0.48 < particle jet η < 0.48 0.48 < particle jet η < grid jet p ⊥ /particle jet p ⊥ 30 < particle jet p ⊥ < 50 GeV 50 < particle jet p ⊥ < 90 GeV 90 < particle jet p ⊥ < 150 GeV 150 < particle jet p ⊥ < 230 GeV 230 < particle jet p ⊥ < 500 GeV  Figure 2. Average shift in the jet p ⊥ introduced by the discretisation as a function of the particle jet p ⊥ and η, respectively.
In Fig. 2 the p ⊥ shift seen in Fig. 1 is quantified. The ratio between grid jet p ⊥ and the particle jet p ⊥ is seen to be around 1.04 and thus reasonably close to unity, and largely independent of jet p ⊥ and η for p jet ⊥ > 50 GeV. Such shifts usually are corrected in experiments [44,45] by introducing detector level correction factors as a function of the jet p ⊥ and η. In this paper, GridSub jets are not corrected for this shift in their p ⊥ , since it is reasonably small. Also, it partially cancels when looking at ratios of Pb+Pb with p+p due to its independence on jet kinematics. Furthermore, since the mismatch is related to nearby jets, increasing the jet p ⊥ cut leads to a reduction of the effect. As previously mentioned, the GridSub technique sets the cell's four momentum to zero if it contains more thermal than particle momentum. This leads to a systematic undersubtraction, that increases with decreasing cell size. We quantify this effect using the event sample with medium response included. Jets are reconstructed and subtracted using the default grid subtraction, but here we keep track of the energy of cells whose four-momentum is set to zero. For each jet we then check if it contains such cells and sum the (negative) energy that these cells originally had. The sum of the negative energy per grid jet is shown in Fig. 3 for different jet p ⊥ ranges. The contribution of negative energy, i.e the amount of thermal energy that remains un-subtracted from the jet, is largely independent of the jet p ⊥ (except for the lowest p ⊥ bin) and small compared to the jet p ⊥ over most of the covered p ⊥ range.

Comparison of two GridSub versions
Jet clustering before discretisation vs. Discretisation before jet clustering As discussed in section 3.2 we have implemented two versions of the grid based subtraction, that differ in the order of jet clustering and discretisation. It is to be expected that the two versions yield different results, as there is no reason why the two operations should commute. Using again the hadron level event sample with medium response included we quantify the differences between the versions. To this end, we find and subtract jets with both versions and event-by-event match the jets following the procedure detailed above. The relative difference in jet p ⊥ and shift of the axis due to the different ordering of operations is shown in Fig. 4. Both these effects are determined to be quite small, but the jet p ⊥ is consistently larger, when the initial jet clustering is performed before the discretisation of the event.  A final check of the two subtraction methods (4MomSub and GridSub1) is done at parton level, where the same jets can be reconstructed with and without recoiling partons. The subtraction is performed with either of the two methods and the matching procedure is again the same as before. Fig. 5 shows the relative p ⊥ difference between jets reconstructed with and without recoiling partons. As expected, the 4MomSub distribution is narrower compared to Grid-Sub1, due to additional jet smearing introduced by the discretization of the event into cells of a finite size. Additionally, the 4MomSub distribution has a tail on the positive side. This is a momentum conservation effect: the thermal distribution is isotropic (except for the longitudinal boost), while the recoiling partons have a net momentum in direction of the jet due to momentum conservation. Therefore, when including medium response more momentum is added to the jet than is subtracted. This is a physical effect that is independent of the subtraction method, but for the GridSub method the shoulder is towards the negative side. This is due to the aforementioned nature of the GridSub to under-subtract the jets, which overcompensates the momentum conservation effect.

Application to traditional jet quenching observables
Observables built from the jet p ⊥ and axis, such as jet R AA or the di-jet asymmetry A J , for smaller radii jets typically show a rather mild sensitivity to medium response. The jet axis is dominated by the hard jet components and for the jet p ⊥ the only effect of medium response is a partial recovery of lost energy. For small reconstruction radii, this is at best a moderate effect, while for very large radii, such as R ≈ 1.0, the effect becomes sizable. For such large radii also the systematic uncertainties related to the subtraction become large. Experimentally, the study of such large jets in a heavy ion environment constitutes an almost impossible task of discriminating between underlying event and the jets. For small radii jets at small momenta the same problem persists, which is why different experiments utilize different procedures to remove the effect of the underlying event in the jet collection of interest [2,46,47]. As our primary validation, Fig. 6 shows the nuclear modification factor R AA of jets, i.e the ratio of jet yield in Pb+Pb over binary collisions scaled p+p, for a moderate radius of R = 0.4. As expected, including medium response leads to a small increase of R AA over the entire jet p ⊥ range. The grid based subtraction leads to a significantly larger increase. This reflects the under-subtraction of the GridSub method discussed in section 5. Increasing the cell size leads to a reduction of R AA . There is good agreement between the two versions of the grid subtraction.
The jet radius dependence of R AA is shown in Fig. 7 with medium response and 4Mom-Sub. The expected increase of R AA with R, because with increasing jet radius more and more of the lost energy is recovered, is indeed observed 4 .

Figs. 8 and 9 show the di-jet momentum asymmetry
This is in contrast to the behaviour observed in [15], where RAA decreases with increasing jet radius because wider jets are more easily lost and medium response cannot compensate this loss.  and relative azimuthal angle ∆φ 12 , respectively. Here, the leading jet is required to have p LeadJet ⊥ > 100 GeV and the cut on the sub-leading jet is p SubLeadJet ⊥ > 30 GeV 5 . The momentum asymmetry A J is calculated without ∆φ 12 cut. The jet axis and thus ∆φ 12 are unaffected by medium response, while in the case of A J it leads to a mild reduction of the medium modification obtained without medium response.

Application to jet shape observables
In contrast to the observables discussed in the previous section, that aim at characterising global properties of jet events, jet shape observables are sensitive to the momentum distribution inside the jet. The latter are thus more affected by medium response. The energy in QCD jets is very much concentrated towards the jet axis, while medium response leads to a much broader distribution of relatively soft activity. Also the fluctuations of the two components are different. In this section we discuss a number of jet shape observables and how they are affected by medium response in Jewel. The reconstructed jet mass is a good probe of medium induced jet modifications and medium response, since it is sensitive to the soft sector. Fig. 10 shows the Jewel+Pythia results for the jet mass distribution. The Monte Carlo shows a shift towards larger masses when medium response is included, whilst for events generated without recoils, a smaller jet mass is observed for jets belonging to the same kinematic range. The latter is due to the known narrowing of the hard jet core. The partial cancellation between two competing effects -the narrowing due to energy loss and the broadening due to medium response -is typical for this kind of observables and also seen in other jet shapes (e.g. the jet profile and girth). We observe a large difference between 4MomSub and GridSub subtraction in this observable, but good agreement between the two versions GridSub1 and GridSub2. In fact,   Figure 12. Left: Jet mass distributions from Jewel+Pythia for p+p. The blue histogram shows the full jet distribution, the red the one for charged jets, and the green histogram is the re-scaled blue histogram. Right: Comparison of the re-scaled full jet mass distribution with recent ALICE data [12] for the charged jet mass.
the jet mass is very sensitive to the details of the grid subtraction. In Fig. 11 we compare two different cell sizes and two ways of computing the cell momentum. One is the default, which consists of summing the four-momenta of the particles in the cell (and subtracting the thermal momenta), and the other sums the particles' energies and assumes the cell fourmomentum to be massless and to point in the direction defined by the cell centre. Both variations lead to large differences in the jet mass distribution (which is not observed in any other observable we studied). We therefore strongly discourage the use of grid subtraction for the jet mass and from here on show results only for 4MomSub subtraction.
As discussed in section 3.3, in order to be able to compare the Jewel+Pythia results to the ALICE data, the charged jet mass has to be estimated from the full jet mass. We do this by re-scaling the full jet mass with a constant factor 2/3 and the jet p ⊥ with a factor 3/4 (this is needed to match the p ⊥ cuts in the charged jet sample). The scaling factors are extracted from the Jewel+Pythia p+p sample. The left panel of Fig. 12 shows the charged jet, full jet and re-scaled full jet mass distributions in p+p and gives a lower bound on the related systematic uncertainties. We would like to stress once more that this is an ad hoc procedure and that there is no guarantee that it yields meaningful results. The right panel of Fig. 12 shows the comparison of the re-scaled full jet mass distribution from Jewel+Pythia to a recent ALICE measurement [12]. The Monte Carlo predicts significantly larger jet masses, but given the uncertainties involved in obtaining the charged jet distribution it is difficult to interpret this comparison with data. Intra-jet fragmentation function [5][6][7] in p+p and Pb+Pb collisions are also an important jet sub-structure observable. However, in Jewel there is no way of doing the subtraction for individual hadrons or, as in this case, tracks. In Fig. 13, which shows the modification of the fragmentation function in Pb+Pb collisions compared to p+p, we therefore in the sample with medium response correct the jet p ⊥ , but all tracks enter the fragmentation function. It is thus expected that Jewel+Pythia overshoots the data in the low z or p ⊥ , corresponding to high ξ, region. The sample without medium response in this region shows a suppression as opposed to the enhancement seen in the data and the sample with recoiling partons, confirming the interpretation that the low p ⊥ (high ξ) enhancement seen in the data is due to medium response. The enhancement at high p ⊥ (low ξ) region is caused by the already mentioned narrowing and hardening of the hard jet core, and is more pronounced in Jewel+Pythia than in data. It is stronger without medium response, because the latter does not affect the hard fragments, but slightly increases the jet p ⊥ .  Figure 14. Ration of the differential jet shape (or jet profile) in Pb+Pb and p+p measured by CMS [8] (black points) and compared with Jewel+Pythia results with (blue line) and without medium response (green line). The data systematic uncertainties are shown in the yellow band around unity.

Fragmentation functions
The differential jet shape or jet profile ρ(r) measures what fraction of the jet p ⊥ is found at what distance from the jet axis. It is defined as where the sum runs over all particles in the jet. The CMS measurement [8] was performed using the full jet p ⊥ , but ρ(r) was built only from tracks. Therefore, as is the case of the fragmentation function, we can do the subtraction for the jet p ⊥ , but not for the charged particles. In this case, however, this is not a problem, since the jet profile built from tracks and the one built form all particles differ only by a constant factor. Assuming this factor to be the same in p+p and Pb+Pb, it will cancel exactly in the ratio of the jet profiles. We can therefore compare Jewel+Pythia results for full jets directly to the CMS data on the jet profile ratio. A more serious problem is that in experimental analysis only tracks with p trk ⊥ > 1 GeV are included. Since we can only subtract for the inclusive final state, this leads to a small mismatch, that becomes visible only at large r and reaches up to 10 % in the highest r bin. Fig. 14 shows the Jewel+Pythia result compared with CMS data [8] for the modification of the differential jet shape ρ PbPb /ρ pp in Pb+Pb collisions compared to p+p. Including medium response and after performing the subtraction using the 4MomSub method, we are able to reproduce the general trend of the data. Jewel+Pythia with recoiling partons describes the enhancement of the jet shape at large radii mostly due to soft particles (p ⊥ < 3 GeV), while without medium response the enhancement is entirely absent.  [10] for charged jets (R = 0.2 and 40 GeV < p jet ⊥ < 60 GeV) compared with Jewel+Pythia for full jets (with adjusted p ⊥ range). The yellow shaded region around unity represents the data systematic uncertainties.
The first radial moment of the jet profile is called girth [48] and is defined as where the numerator sums the distance from the jet axis weighted with p (k) ⊥ of each constituent k of the jet. It characterises the width of the p ⊥ distribution inside the jet.
Jewel+Pythia results for girth using GridSub1 subtraction for fully reconstructed jets in central Pb+Pb collisions are shown in the left panel of Fig. 15. We find a shift to smaller values of g due to narrowing of the hard component, which is partly compensated by a broadening of the jet due to medium response. We also compare our results with preliminary ALICE data [10] for charged jets in the right panel of Fig. 15. Following the same argument as above for the jet profile, the girth of full and charged jets should be the same, provided the p ⊥ range is adjusted accordingly. We confirmed this in the Monte Carlo for p+p collisions. We therefore in Fig. 15 compare Jewel+Pythia results for fully reconstructed jets at a correspondingly higher p ⊥ with the ALICE data. We find reasonable agreement, but the Jewel+Pythia distribution peaks at slightly higher values than the data.   [11]. Note that the data is not unfolded, but the resolution is not published so no smearing is applied to the Monte Carlo events. A comparison to properly smeared Jewel+Pythia results can be found in [11]. The yellow shaded region around unity in the left panel highlights the statistical uncertainty in the p+p reference and on the right represents the the data systematic uncertainties.

Groomed shared momentum fraction z g
The groomed shared momentum fraction z g is a measure for the momentum asymmetry in the hardest, i.e. largest angle, two-prong structure in the jet. In p+p collisions it is closely related to the Altarelli-Parisi splitting function [49]. It is defined through the Soft Drop procedure [50,51] detailed below and implemented in FastJet [36] contrib. First, jets are clustered with the anti-k ⊥ algorithm and re-clustered with Cambridge/Aachen. Then the last clustering step is undone, yielding the largest angle two-prong structure in the jet. If this configuration satisfies the Soft Drop condition where z cut and β are parameters, it is kept. Otherwise, the softer of the two prongs is discarded and the procedure of un-doing the last clustering step is repeated for the harder prong. In this way soft contaminations are systematically removed from the jet and the hardest two-prong structure is identified. Soft Drop jet grooming thus takes an inclusive jet collection and turns it into a different collection of jets with two-prong structure of a minimum momentum symmetry provided by z cut . Varying z cut up or down varies the degree of asymmetrical splitting allowed in the parton's fragmentation, while the β controls how collinear the configuration has to be. In p+p collisions, this method is has been studied in some detail [49,51], but in heavy ion collisions the exact meaning of the grooming procedure is not obvious, due to the presence of the fluctuating underlying heavy ion event and the increased soft sector, that the procedure tries to remove. Recent analytical studies [52] have shown that grooming increases the sensitivity to medium induced gluon bremsstrahlung thus experimentally opening up different avenues to directly probe the effect of the medium on a jet by jet basis. In Jewel, however, a different story unfolds.
As shown in the right panel of Fig. 16, there is an increase in asymmetrical splittings in Pb+Pb jets as opposed to p+p jets, which is observed in recent preliminary CMS results and reproduced in Jewel+Pythia. The secondary feature observed in this measurement is an apparent reduction of the effect for higher p ⊥ jets. Jewel reproduces this behavior qualitatively as well, with very high p ⊥ jets showing very little difference in the momentum fraction of the first splitting [11]. The left panel of Fig. 16 shows that in the Monte Carlo the modification of the z g distribution in Pb+Pb collisions is partly due to the narrowing of the jet, as seen in the sample without medium response. The more important contribution, however, comes from adding the recoiling partons 6 . In Jewel+Pythia we see no sign of medium induced bremsstrahlung contributing to the effect, as advertised in [52].

Discussion and conclusions
Studies of jet sub-structure modifications in heavy ions probe the intricate interactions between the medium and jets. Due to their sensitivity to medium response, they offer the power to discriminate between several models and shed light on the underlying jet quenching mechanisms as well as the thermalisation of the deposited energy and momentum.
In Jewel it is possible to study medium response in detail by keeping the partons recoiling against interaction with the jet in the event. One has to keep in mind that this is only a limiting case, since these partons do not undergo further interactions in the medium. In order to be able to compare these results to experimental data, the thermal component of the recoiling partons' momenta has to be subtracted. In this paper we introduced two methods for doing this, a four-momentum and a grid based one. With these tools we can for the first time quantitatively study jet shape observables.
We find that -at least in Jewel+Pythia -a number of qualitative feature in the data can only be explained by medium response. These are • the increase at low z of the ratio of intra-jet fragmentation functions in Pb+Pb compared to p+p, • the increase of the jet profile at large distance from the jet axis in Pb+Pb compared to p+p, • and the enhancement of asymmetric two-prong structures in Pb+Pb compared to p+p as seen in the z g distribution.
This is in line with observations by other authors [15,17]. In other observables, in particular the jet mass and girth, a non-trivial cancellation between a narrowing of the jet core due to energy loss [13][14][15] and a broadening due to medium response takes place. Also in the case of girth, including medium response leads to an improvement of the agreement between Jewel+Pythia and ALICE data. For the jet mass we find that the grid based subtraction does not yield reliable results. The 4MomSub subtration should be more robust, but without grid subtraction we do not have an independent way of cross-checking the results. We therefore recommend not to use GridSub for the jet mass and to take the comparison of Jewel+Pythia results to the ALICE data with a grain of salt.
Jet shape observables open a new perspective on jet quenching and may also help to address the question of thermalisation, and it is important to develop tools capable of quantitatively describing medium response. The present study with Jewel can only be a first step in this direction. As emphasised above, the treatment of recoiling partons is still schematic. The subtraction methods introduced in this paper are solid, but have their limitations, in particular when it comes to the description of charged jets. It is currently also impossible to perform the subtraction for particles (for instance in the fragmentation functions), due to the mix of parton and hadron level in the subtraction. The grid method also introduces systematic uncertainties related to the discretisation, that can, however, be quantified (cf. section 5). Nevertheless, the results for jet shapes obtained with Jewel+Pythia are very promising. In some cases this is the first time that they can be studied quantitatively in a consistent jet quenching model including medium response.
Upcoming measurements at the LHC will further advance the understanding of jet shapes by utilizing the jet grooming tools, amongst others. This ushers in a new era of sub-structure studies in heavy ion collisions, where correlation between different observables could point the way to the future in decoupling several of the physics features hidden in individual observables.