On heavy-flavour jets with Soft Drop

We study hadronic jets that are tagged as heavy-flavoured, i.e. they contain either beauty or charm. In particular, we consider heavy-flavour jets that have been groomed with the Soft Drop algorithm. In order to achieve a deeper understanding of these objects, we apply resummed perturbation theory to jets initiated by a massive quark and we perform analytic calculations for two variables that characterise Soft Drop jets, namely the opening angle and the momentum fraction of the splitting that passes Soft Drop. We compare our findings to Monte Carlo simulations. Furthermore, we investigate the correlation between the Soft Drop energy fraction and alternative observables that aim to probe heavy-quark fragmentation functions.


Introduction and motivation
Beauty (b) and charm (c) quarks are often considered "heavy flavours" because their masses are above the proton mass, m b = 4.2 GeV and m c = 1.3 GeV, respectively.However, their mass is not so large, compared to the typical scale of hadron formation, Λ ≃ 1 GeV, so hadronisation occurs.This does not happen for the top quark, which has a mass 170 times bigger than the proton mass.This value is so large that its lifetime is smaller than the hadronisation scale.
Heavy flavours constitute a window on two mechanisms that provide ordinary matter with mass: electroweak symmetry breaking and the binding energy of strong interactions.On the one hand, they play a crucial role in studies of the Higgs boson and, on the other hand, they constitute a bridge between perturbative and non-perturbative Quantum Chromo Dynamics (QCD).For these reasons, they have been (and still are) the subject of detailed theoretical and experimental studies.a e-mail: scaletti@phys.ethz.chb e-mail: andrea.ghira@ge.infn.itc e-mail: simone.marzani@ge.infn.it In this work, we focus on processes in which heavy quarks are produced.From a theoretical point of view, calculations for identified heavy flavours can be performed essentially because the quark mass sets a perturbative scale for the running coupling and, at the same time, removes collinear singularities.From an experimental viewpoint, the lifetime of B (or D) hadrons is long enough so that their decay happens away from the interaction point.Dedicated band ctagging techniques that exploit this property to identify B and D hadrons, or b and c jets, are widely used in collider experiments.
Flavour physics has been studied for decades, both in the quark and lepton sectors.However, the recent development of Infra-Red and Collinear (IRC) safe flavour-jet algorithms [1][2][3][4][5] opens up the possibility of setting up a yetunexplored flavour physics program that exploits jets and their substructure at the Large Hadron Collider (LHC).
In contrast, there exists a rather extensive literature dedicated to studying the properties of a reconstructed B (or D) hadron, such as its energy or its transverse momentum.The theoretical description of these observables is usually based on heavy-quark fragmentation functions, which can be computed in perturbative QCD.State-of-the-art predictions include the resummation of different classes of logarithmic corrections, e.g.mass logarithms and end-point log-arXiv:2312.11623v2[hep-ph] 29 Feb 2024 arithms, see e.g.[27][28][29][30][31][32][33][34].In this context, two of us have recently developed a theoretical framework to consistently resum both mass and soft logarithms [35,36].In this study, we will heavily rely on these results and we will apply them to jet-based observables.This way, we will be able to obtain theoretical predictions that resum both the logarithms of the observable we want to study and the logarithms of the ratio of the heavy-quark mass to the jet transverse momentum.
In this study, we investigate the possibility of exploiting a widely used jet substructure technique, namely Soft Drop [37], to study heavy-flavour jets.The Soft Drop algorithm has been extensively studied and Soft Drop jets are routinely used in experimental analyses, both in the context of measurements and searches for new physics.In particular, we consider heavy-flavour-tagged jets and focus on the angular separation (θ g ) and momentum fraction (z g ) of the first splitting that passes Soft Drop.The former directly measures the angular resolution of the groomed jet and, therefore, it gives us direct access to the dead cone.The latter, instead, allows us to probe the heavy-quark splitting function [38,39].
These observables have been recently measured on c-jets by the ALICE collaboration [40] and our study constitutes the first step towards a first-principle description of the data.A direct comparison to the ALICE data, however, goes beyond the scope of this paper, as it would require matching the resummed expressions to fixed-order matrix elements, as well as accounting for the fact that ALICE reconstructs c-jets from exclusive D-meson decays.
The paper is organised as follows.We start by reviewing the Soft Drop algorithm in Sec. 2. In Sec. 3 we perform the next-to-leading logarithmic (NLL) resummation of the θ g distribution, including heavy-quark mass effects and compare it Monte Carlo parton-shower simulations.In Sec. 4 we study the z g distribution, both from first-principle and in Monte Carlo simulations and study its correlation to a standard momentum fraction variable used in the context of fragmentation functions.Finally, we draw our conclusions in Sec. 5. Details of the calculations are collected in the Appendices.

Soft Drop
Soft Drop [37] is a grooming algorithm that recursively removes soft-wide angle constituents from a jet.The Soft Drop procedure starts by re-clustering a given jet (typically an anti-k t [41] jet) with radius R 0 and transverse momentum p t with the Cambridge-Aachen (C/A) algorithm [42,43].Soft Drop then parses the resulting angular-ordered branching history, grooming away the softer branch, until a branch that satisfies the condition min (p t1 , p t2 ) is found.In the expression above, 1 and 2 denote the branches at a given step in the clustering, p ti are the corresponding transverse momenta, and ∆ 12 is their rapidity-azimuth separation.If the condition above is never satisfied, we can either remove the jet from consideration ("tagging mode") or leave it as the final Soft Drop jet ("grooming mode").
The kinematics of the first branch that satisfies (1) defines the groomed jet radius θ g and the groomed momentum sharing z g : The θ g distribution is IRC safe and it was first studied, for light jets, in [37] and then resummed to next-to-leading logarithmic (NLL) accuracy in [44].The momentum sharing z g is IRC safe for β < 0 but only Sudakov safe for β ≥ 0 [38].
The NLL calculation of the z g distribution was performed in [45].
3 The θ g distribution In this section, we describe the resummation of the θ g distribution for a jet originated by a massive quark.We follow the general strategy described in Ref. [37] and we, therefore, consider the resummation of the cumulative distribution, normalised to the Born cross section, exploiting Lund diagrams [46].

Lund plane geography
Lund diagrams are a useful way to represent the available phase space for the emission of soft and collinear gluons off a hard dipole.The Lund plane is usually represented in terms of pairs of the logarithmic variables, so that, in the soft and collinear limit, equal areas correspond to equal emission probabilities.We choose to draw Lund planes as depicted in Fig. 1.On the horizontal axis, we show the logarithm of the emission angle with respect to the jet axis, as measured in the rapidity-azimuth plane, normalised to the jet radius, i.e. θ = ∆ R 0 .On the vertical axis, we show the emission's transverse momentum with respect to the jet axis, normalised to the hard scale of the process, κ = k t p t R 0 .This way, lines of Recently, two of us have achieved a generalisation of the Lund-plane resummation formalism that includes quarkmass effects [36].The Lund plane on the left-hand panel of Fig. 1 is for a b-quark line, while the one on the right is for a c-quark.When considering heavy-quark Lund planes, two important differences with respect to the massless case arise.First, the presence of the vertical short-dashed line, in black, that represents the beginning of the dead cone at θ = m b p t R 0 , i.e. ∆ = m b p t , for the left-hand figure, and at θ = m c p t R 0 , for the right-hand one.The shaded area in grey is the dead-cone region.Here, collinear emissions do not give rise to any logarithmic enhancement, leading to a suppression of QCD activity in this region.Second, the presence of heavy-quark thresholds, which are relevant when considering the running of the strong coupling.In Fig. 1, they are represented by horizontal dotted lines in red.The lines at κ = m b,c p t R 0 correspond to k t = m b,c and therefore mark the boundaries between two regions with different numbers of active flavours: n f = 5 above and n f = 4 below, and n f = 4 above and n f = 3 below, for k t = m b and k t = m c , respectively.We also show (in green) the line κ = Λ p t R 0 that marks the k t = Λ ≃ 1 GeV boundary between the perturbative and the non-perturbative regions.For the region below this line, we need a prescription to regulate the Landau pole of the strong coupling.De-tails are given in Appendix A. Finally, we show dotted lines in blue, at κ = z c and at κ = z c m b,c p t R 0 1+β , which originate from the Soft Drop condition and the dead-cone boundary.
When performing actual calculations, one needs to establish the hierarchy between these scales.Let us consider, as an example the case of b-jets, with relatively high-p t .Because Soft Drop is usually employed with z c = 0. Λ .However, this inequality is hardly satisfied if we consider commonly used values of the Soft Drop parameters, e.g.z c = 0.1 and β = 0. Thus, we do expect non-perturbative contributions to affect the deadcone region and we will assess their size using Monte Carlo simulations.In Fig. 1, we show the Lund plane for this hierarchy of scales, on the left for b-jets and on the right for c-jets.A discussion of all the other possible hierarchies and regions can be found in Appendix B.

Recap of the massless calculation
Before computing the θ g resummed distribution for a jet initiated by a massive quark, let us briefly recall the structure of the resummation in the case of massless partons.The allorder cumulative distribution can be written as the product of two contributions The function R is the Sudakov exponent, or radiator, which accounts for angular-ordered collinear emissions.It was first computed in [37] and, for the case of a quark jet, reads where θ is the emission angle, normalised to the jet radius, k t = zθ p t R 0 and the massless q → qg (timelike) splitting function is In order to achieve NLL accuracy, the running coupling must be evaluated at two loops in the so-called CMW scheme [47].
Henceforth, we will also work in the small-z c limit.Consequently, even in the case β = 0, we will ignore flavourchanging contributions [48,49].The Sudakov form factor represents the no-emission probability and corresponds to the area shaded in pink on the Lund planes of Fig. 1.This rather simple picture essentially arises from two facts.First, for the θ g distribution, independent emissions, i.e. α n S C n F contributions, exponentiate with no corrections due to multiple emissions [37].Second, all non-Abelian contributions are unresolved and therefore captured at NLL by the running coupling in the CMW scheme.However, as it was pointed out in [44], the above description is not complete at NLL and the resummed expression must be supplemented with a correction factor, S in Eq. ( 4).Indeed simple exponentiation is broken by single logarithmic contributions that arise from two competing mechanisms.Let us start with the correction to the independent emission contribution.The C/A algorithm, which is used as the first step of the Soft Drop procedure, can cluster two soft gluons together first, if they are the closest pair in angle, instead of clustering each of them with the hard quark.This happens when θ 12 < min(θ 1 , θ 2 ).This effect introduces single-logarithmic corrections that are usually referred to as clustering logarithms [50][51][52][53].On the other hand, the clustering condition also plays a role in the correlated emission contribution.The C/A algorithm can resolve a soft gluon splitting, giving rise to a non-global contribution [54,55] that spoils the CMW picture.This is also a single logarithmic contribution that can happen if the two soft gluons are not the closest pair in angle, i.e. θ 12 > min(θ 1 , θ 2 ).
Both contributions start at O(α 2 S ) and can be computed using the expression of the squared matrix element for the emission of two strongly-ordered soft gluons with momenta k 1 , k 2 off a hard fermionic dipole with momenta p a , p b [56,57] where the C 2 F term describes the independent emission contribution, while the C F C A one, the correlated one.We have introduced the eikonal factor Following [44], we consider the case in which the softer gluon (k 2 ) is emitted at a large angle and passes the Soft Drop condition, while the harder gluon (k 1 ) is emitted at an angle smaller than θ g and so it is not subject to the Soft Drop condition.If we work in the small-angle limit, considering real and virtual contributions, we find with The integration above can be simplified by noting that, in order to capture the highest power of the logarithm of θ g , we can evaluate the momentum fraction integrals with lower limit This way, the momentum fraction integrals decouple from the angular ones and we obtain The all-order resummation of these contributions has been performed numerically, in the large-N c limit [44].In this study, we have decided to limit ourselves to studying the impact of these corrections.Thus, we approximate the clustering factor S by simply considering the exponentiation of the two-loop result: where the running of the strong coupling is taken at one loop.

The θ g distribution for a heavy-flavour jet
We now perform the resummation of the θ g distribution for jets initiated by a heavy quark, considering both the case of a b-quark and a c-quark.Both the Sudakov contribution R and the clustering correction function S in Eq. ( 4) need to be reconsidered.We start with the calculation of the resummed Sudakov exponent.It was realised long ago that squared QCD matrix elements with massive partons factorise in the so-called quasicollinear limit [58,59].In this approximation, both the transverse momentum k t of the emission with respect to the massive emitter, and the mass m are small compared to the hard scale but they are considered of the same order.In this limit, the squared invariant amplitude for one-gluon emission takes the form where the massive splitting function for i → ig is Following the same steps as in [36], we write where θ = k t zp t R 0 .Note that the mass dependence enters in two ways.First, for a jet initiated by a quark flavour i = b, c, the quasi-collinear phase-space and the massive splitting function depend on m i .Second, regardless of the jet flavour, the running of the strong coupling may cross the b and the c thresholds, thus inducing a logarithmic dependence on the quark masses.In order to proceed further, it is convenient to change the integration variable in Eq. ( 14) from k 2 t to θ 2 : where now k 2 t = z 2 θ 2 p 2 t R 2 0 , and we have introduced the following dimensionless variables We stress again that the Sudakov exponent for either a highp t b-jet or c-jet depends on both ξ b and ξ c because the running coupling crosses both quark thresholds.However, it depends only on the dimensional ratio θ 2 i of the corresponding flavour i = b, c.
The mass-dependent shift in the denominator of Eq. ( 15) acts as an effective lower bound of a logarithmic angular integration.This is the well-known dead-cone effect, [24,25] i.e. the fact that radiation off massive partons at angles below m/p t is not logarithmically enhanced.Therefore, we can further simplify our expression by shifting the angular integration variable to θ 2 = θ 2 + θ 2 i .Neglecting power corrections in the mass, we have Let us first note that Furthermore, as discussed in detail in [36], the mass-dependent shift in the argument of the running coupling contributes at most to NNLL corrections and, therefore, it can be dropped at the accuracy we are working at.A similar argument applies to the shift in the Soft Drop condition.This can be easily checked at fixed coupling.For instance, in the soft and collinear limits, we have where in the last step we have dropped constant terms in the θ i → 0 or θ g → 0 limits, which are NNLL.Thus, the massive Sudakov exponent at NLL becomes Remarkably, with the sole exception of the mass-dependent term in the splitting function, this result, when seen as a function of ϑ 2 g,i , has the same form of the massless case, Eq. ( 5).Consequently, we can represent the massive Sudakov form factor in Fig. 1 using the same Lund planes as in the massless case, provided that we interpret vertical lines as lines of constant ϑ 2 g,i .This way, the dead-cone line acts as a phase-space boundary because ϑ 2 g,i → θ 2 i as θ g → 0. The integrations in Eq. ( 20) are all straightforward but the presence of the Soft Drop condition, of the dead cone, and of mass thresholds force us to consider many cases.Details are given in the Appendices.
Next, we have to consider the clustering correction factor S. We adopt the same, approximate, strategy of the massless case.Namely, we consider the running-coupling exponentiation of the two-loop result.The square amplitude for the emission of two soft gluons off a massive quark-antiquark dipole was computed in [60,61].For our purposes, we are actually interested in a less general result, namely the case of the emission of two soft gluons that are strongly ordered in energies, e.g.z 2 ≪ z 1 , off a single dipole.In this limit, the massive square matrix element takes the same form as Eq. ( 7) but with the massless eikonal factor w i j,k substituted by the massive one: We choose p a , see Eq. ( 7), to be the momentum of the heavy quark and we work in the quasi-collinear limit with respect to its momentum direction.Therefore, we can set p 2 b = 0. We then integrate the C 2 F and C F C A contributions with the same constraints as in Eq. ( 9): where, as before, θ 2 As in the massless case, in order to obtain the leading contribution we can substitute θ 2 ≃ θ g in the Soft Drop condition.Therefore, the momentum fraction integrals decouple from the angular ones and we obtain The Abelian contribution reads with For the non-Abelian case, instead, we find We evaluate F 1 and F 2 numerically.However, in the limit θ 2 g ≫ θ 2 i we are able to perform the integrals analytically and we obtain lim so that the massless result Eq. ( 10) is recovered.Remarkably, in the opposite limit, we find lim so that these contributions disappear asymptotically.We note that this result is not related to the clustering conditions but rather to the requirement that one of the two gluons (k 1 ) be within the groomed jet radius, i.e. 0 < θ 2 1 < θ 2 g .Changing the integration variable to θ 2 1 , we have θ 2 i < θ 2 1 < θ 2 g + θ 2 i .Thus, we have to evaluate the integral of a (regular or integrable) function over a domain that has a vanishing measure as θ 2 g → 0.
As for the massless case, we do not perform the full NLL resummation of these effects, but we limit ourselves to exponentiate the O(α 2 S ) with running coupling corrections, as done in Eq. ( 11): In summary, the resummed θ g distribution for a heavyflavour jet reads 1 for i = b, c.The Sudakov exponent R i is given in Eq. ( 20) and the clustering contribution S i in Eq. ( 29).Note that our NLL distribution is no longer normalised, as a consequence of the fact that, in the massive case, the resummed cumulative behaves as a constant for θ g → 0, while it vanishes in the massless case, see Eq. ( 4).Indeed, while in the massless case, Soft Drop with β ≥ 0 behaves perturbatively as a groomer, i.e. within resummed perturbation theory, it always returns a jet with θ g > 0, the quark mass provides an effective cutoff so that, there is a non-vanishing probability, given by to find Soft Drop jets with θ g = 0, even if we expect fixedorder corrections, as well as non-perturbative effects, to smear this effect out.Thus, we could either use Soft Drop in grooming mode and supplement Eq. ( 30) with an endpoint, δ (θ g ), contribution that ensures normalisation or consider the algorithm in tagging mode and discard jets that return θ g = 0.
Henceforth, we choose this second option.Finally, we note that, in the calculations of the resummed exponent presented so far, we have always included virtual corrections, as well as the contributions from jets that fail Soft Drop.In tagging mode, one should discard those and, in principle, repeat the calculations.However, we observe that in resummed perturbation theory Thus, in practice, we can just use the results obtained so far, provided that we remove the θ g = 0 contribution.Numerical results and their comparison to Monte Carlo simulations will be presented in the next section.

Numerical results and comparison to Monte Carlo simulations
We now provide numerical results for the θ g distribution according to Eq. ( 30), and compare them to Monte Carlo (MC) simulations.To this purpose, we generate events using HERWIG, version 7.2.2 [62].To test the all-order behaviour of the observables we are interested in, we consider LO matrix elements for the hard process, dressed with the parton shower. 1 Hadronisation effects are included using the HERWIG cluster hadronisation model, when explicitly declared, and the CT14 [66] set of PDFs has been used throughout the paper. 1 In this paper we concentrate on an angular-ordered parton shower, see [63] for detailed comparisons with analytic resummation.We have also tested our analytic predictions against the HERWIG dipole shower [64] and the PYTHIA8.3 [65] one and we did not find significant differences.It would be interesting, in the future, to also compare to stateof-the-art simulations that also include matching to fixed order.
We simulate the inclusive production of a pair of oppositely charged muons in association with a jet, in protonproton collisions at 13 TeV centre-of-mass energy.The muon pair is required to have an invariant mass between 70 and 110 GeV.Jets are clustered using the anti-k t algorithm [41] with R 0 = 0.4 and then ordered in transverse momentum.The hardest jet containing a b (c) quark, or a B (D) hadron, is considered.To match the heavy quark/hadron with a jet we look at the closest (with respect to the jet axis) flavoured particle, starting from the highest in transverse momentum, with p t > 5 GeV.Properties of the flavoured particle and the jet are extracted and analysed using RIVET [67] and FASTJET [68].The fiducial phase-space for the muons is defined by the following cuts: p t,µ > 26 GeV, |η µ | < 2.4, while the jets are selected in the region defined by |η J | < 2.4.We consider three different transverse-momentum regions, i.e. p t ≥ 50, 150 and 300 GeV.However, only p t ≥ 150 GeV is shown in the main text, and the other cases can be found in Appendix C.
The θ g distribution is shown in Fig. 2. Each plot includes parton level, i.e. with parton shower effect only, and hadron level, i.e. with hadronisation and the Underlying Event (UE) included.Together with the MC prediction, we show the corresponding NLL result of Eq. (30).The latter also exhibits a theoretical uncertainty given by the variation of the renormalisation scale, for which the central value is set at the hard scale p t R 0 , by a factor of two, as customary.Plots in different rows correspond to results for the identified leading jet in Z + b (B), Z + c (D) and Z + light quark/hadron production, respectively.We always consider Soft Drop with z c = 0.1, and β = 0 and 1, for left and right plots, respectively.For b and c jets, we also indicate, in green, the expected dead-cone region θ g < θ i = m i p t R 0 .All curves are normalised to have unit area.For β = 0, i.e. when Soft Drop is more aggressive, we find good agreement between the NLL prediction and the MC.They are fairly different, instead, for β = 1 at high θ g .However, in this region, we expect fixed-order corrections, not considered here, to be important.
In order to highlight possible dead-cone effects, it is useful to consider the ratio between the heavy-flavour jet θ g distribution and the corresponding one for light-quark jets.We do this in Fig. 3.We concentrate on b-jets, but we show the results for three different transverse momentum cuts, namely 50, 150 and 300 GeV (from top to bottom).We show results for both β = 0 (left) and β = 1 (right).In every plot, we show the ratio obtained with our NLL resummed prediction and with the simulation performed with HERWIG, both at parton-level and hadron-level.
As expected, because θ i = m i p t R 0 , the dead cone is more visible at lower values of p t .Interestingly, deviations between the b and the light quark distribution start at angular scales bigger than θ i .However, as p t is increased, this transition is pushed to a region that is likely beyond exper-imental resolution.We also note that mass effects are more pronounced for the β = 0 case than β = 1.Considering that we have already noted that the former is under better theoretical control than the latter, the Soft Drop jets with β = 0 appear to be an interesting choice to study the dead cone.

The z g distribution
We now discuss the second variable that characterises Soft Drop jets, namely the momentum fraction z g .In this discussion, we want to provide a simple theoretical prediction for z g for jets initiated by massive quarks.We also compare, both analytically and in Monte Carlo simulation, the variable z g to a widely-used fragmentation variable that measures the transverse momentum fraction of the B(D) hadrons (or b(c) quark at parton level) with respect to the jet p t .

Recap of the massless calculation
We start by briefly reviewing the calculation of the z g distribution for light jets.The value of z g is fixed by the first declustering of the jet that passes the Soft Drop condition.Because we are completely inclusive over the splitting angle, we must integrate over all possible values of θ g , including configurations where the two emissions become collinear, for which the integral diverges.If β ≥ 0, collinear splittings always pass the Soft Drop condition and these divergent configurations are not cancelled by the corresponding virtual corrections, for which z g is undefined, and heralds the fact that the observable is not IRC safe.
However, z g belongs to a wider class, i.e.Sudakov safe observables [37,38,69,70], that despite being IRC unsafe, can be computed in perturbation theory, provided that we use resummation.For this purpose, we need to introduce a safe companion observable.The Soft Drop procedure itself suggests using the groomed angle θ g , which we have discussed in the previous section.Therefore, we imagine to measure a value of z g , given a finite angular separation θ g .Using the language of conditional probabilities, we have [38]: where p(θ g ) = 1 σ 0 dσ dθ g and p(z g |θ g ) = p(z g ,θ g ) p(θ g ) is the conditional probability for measuring z g given a value of θ g .If β < 0, z g is IRC safe and the integral in Eq. ( 33) can be computed order by order in α S .This is no longer true when β ≥ 0, which is the standard configuration in which the Soft Drop algorithm is used and, therefore, our case of interest.In this situation, the integral (33) diverges order by order in the strong coupling because of the 1/θ g behaviour of the integrand.However, if we take p(θ g ) to be the resummed distribution, i.e. the derivative of Eq. ( 4), then the Sudakov form factor regulates the θ g = 0 singularity, providing us with a finite answer for Eq.(33).
Because in this work we are not concerned with logarithms of z g and z c , we evaluate the conditional probability at fixed-order, focusing on the case of a quark-initiated jet:2 p(z g |θ g ) = P gq (z g )α S (z 2 g θ 2 g p 2 t R 2 0 ) The Soft Drop condition (1) dictates z g < 1 2 and, if β = 0, z g > z c .
The integral in Eq. ( 33) with running coupling must be performed numerically.However, it is interesting to study its fixed-coupling and lowest-order approximation, see Eq. ( 19).Thus, we consider For β ≥ 0, we find [38] 1 where erf(x) is the error function and a = min[z g , z c ].If β > 0, the above result is non-analytic in α S : Even more interesting is the β = 0 case.At fixed coupling, the conditional probability, Eq. ( 34), becomes independent of θ g and factors out of the integration to give 1 σ 0 dσ (f.c.) dz g = P gq (z g ) In all plots, we show our NLL resummation, as well as the results obtained with the Monte Carlo event generator HERWIG, both at parton and hadron level.Jets are selected with the anti-k t algorithm with R 0 = 0.4, with p t ≥ 150 GeV, and groomed with Soft Drop with z c = 0.1, and two values of the angular exponent: β = 0 (on the left) and β = 1 (on the right).The uncertainty bands for the analytic predictions are obtained by varying the resummation scale by a factor of two above and below the hard scale p t R 0 , i.e. µ R ∈ p t R 0 2 , 2p t R 0 .It can be shown that the β = 0 case does have a valid perturbative expansion in α S , despite being α S -independent at lowest order.Eq. (38) shows that the β = 0 z g -distribution is essentially driven by the QCD splitting function [38].This observation has initiated numerous theoretical [73][74][75][76] and experimental [77][78][79][80][81][82][83] studies (see also [71,72]) that aim to use z g as a probe of QCD dynamics, both in pp and heavyion collisions.

The z g distribution for a heavy-flavour jet
We now move to the case of b and c jets and we start by considering the simple O(α S ) calculation of the z g distribution, in the quasi-collinear limit: where, as before, we have dropped the mass-dependent shift in the Soft Drop condition.Even before performing the integral, it is clear that the mass of the heavy quark, as one might have expected, regulates the collinear singularity and so the z g distribution is IRC safe, for every value of β .The computation of the integral is straightforward but the presence of the θ 2 dependent contribution in the splitting function complicates the result.However, to illustrate our point, it is sufficient to work at LL. Therefore, we approximate P gi = 2C F /z g , and we find 1 Note that in the case β = 0, we have z g > z c and only the first term survives.Albeit finite, this expression contains logarithms of θ 2 i , which become large in the boosted regime.The all-order resummation of logarithms of z g partially addresses this problem.Indeed, keeping our focus on the LL fixed-coupling approximation, we find the following result for the normalised cumulative distribution: log The β = 0 case is rather simple We note that these expressions indeed resum those logarithms of θ i that are associated with logarithms of z g .However, the θ i → 0 limit is not smooth and we do not recover the massless result of Eq. ( 36).This is related to the noncommutativity of the soft and massless limits, discussed at length in [35,36].
Another way of resumming logarithms θ 2 i is to resort to the conditional probability procedure described above for the massless case.To illustrate the procedure, we repeat, for the massive case, the calculation that led to the fixedcoupling result in Eq. (36).At LL, the splitting function can be approximated by its soft contribution.Moreover, at this accuracy , with R given by Eq. ( 35) and ϑ 2 g,i ∈ [θ 2 i , 1].Therefore, for β ≥ 0, we find In particular, for β = 0, we find The first-order expansions of the above expressions agree with Eq. ( 40) and large logarithmic corrections in θ 2 i are resummed.Furthermore, if we take θ i → 0, we recover the massless distribution of Eq. (36).We further note that Eq. ( 44) is the same as its massless counterpart, but for the normalisation factor in brackets.This makes sense because in the LL fixed-coupling approximation, both distributions are driven by the most singular part of the splitting function, which is the same for P gq , P gb , and P gc .
Because it has the correct θ i → 0 limit, we decide to use the conditional-probability approach to compute the z g distribution of a heavy-flavour jet: where the resummed θ g distribution is given in Eq. ( 45).We can make some further simplifications.The derivative in Eq. ( 45) gives a factor that simplifies, within our accuracy, the denominator of the conditional probability.Thus, we obtain Numerical results and their comparison to Monte Carlo simulations will be presented in the next section.We close this discussion by noting that Eq. ( 45) does not systematically resum logarithms of z g and z c .This is acceptable for our purposes because we are mostly interested in the β = 0 case, for which z g > z c = 0.1.However, it would be interesting to extend the full NLL resummation for z g , performed for light jets in [45] to the massive case.Because the calculation in [45] is based on the resummation of the double differential (θ g , z g ), it may overcome the difficulties about the θ i → 0 limit of Eqs. ( 41) and ( 42), previously discussed.

Numerical results and comparison to Monte Carlo simulations
In this section, we compare our resummed results for the z g distribution to the ones obtained with Monte Carlo event generators, both at parton-and hadron-level.We limit ourselves to the case β = 0 and we simulate events using HERWIG, with the same settings as the ones described in Sec.3.4.
In Fig. 4, we compare resummed results for two different values of the jet transverse momentum cut, namely p t ≥ 50 GeV, on the left, and p t ≥ 150 GeV, on the right.The plots at the top are for b-jets, the ones in the middle for c-jets, and the ones at the bottom for light-quark jets.All curves are normalised to have unit area.As already pointed out, the results all look very similar, because the shape of the distribution is mostly driven by the singular part of the splitting function for the emission of a gluon off a (massive) quark.However, we also note that differences between b,c, and light flavours are larger in the Monte Carlo results than in the analytic ones.This could be due to quark masses influencing the kinematics, thus causing mass-dependent recoil effects.These are, to a certain extent, accounted for in the parton shower, but neglected in the analytic resummation.
Overall, we find good agreement between our calculation and the Monte Carlo results, although our results are strangely closer to the full simulation than to the partonlevel one.We note that our predictions undershoot the Monte Carlo at large z g .We have traced this back to the fact that our calculation does not take into account symmetrised splitting functions, as discussed in the footnote of Sec.4.1.Indeed by rescaling the z g distribution by the symmetrised splitting function, the tail of the analytic calculation moves closer to the Monte Carlo results.

Comparison to fragmentation functions
The study of heavy-flavour production at high energies is a multi-scale problem and, as already pointed out, logarithms of the ratio of the quark mass to the hard scale of the process can spoil the convergence of the perturbative expansion.Because these logarithmic corrections are related to collinear dynamics, heavy-flavour production cross-sections obey a factorisation theorem and they can be written as the convolution of process-dependent partonic (massless) coefficient functions with universal heavy-flavour fragmentation functions.Fragmentation functions obey DGLAP evolution equations with timelike splitting functions, which allow one to resum these large logarithmic corrections to all perturbative orders.
The initial condition for heavy quark fragmentation functions can be computed in perturbation theory, as originally pointed out in Ref. [27,28], where the next-to-leading order (NLO) computation in QCD was presented.The NNLO corrections were computed later in Refs.[29,30].Furthermore, the initial condition of the evolution is affected by soft logarithms, that can be resummed to all orders too.[31][32][33][34]36].Fragmentation functions for b (or c) quarks are usually then supplemented with non-perturbative corrections before being compared to experimental data.In the context of heavyflavour production in e + e − collision, one of the most widely studied observables is perhaps the energy fraction x of the heavy quark (or hadron) with respect to the energy of the incoming beam, in the centre-of-mass frame, see [84] (and references therein) for a recent review.
Collinear factorisation ensures that fragmentation functions are universal objects and so they can be used, in principle, to describe heavy flavours in hadronic collisions.However, a few changes in their definition are usually required from a practical point of view.For instance, although is possible to study directly the properties of identified B or D decays, the LHC experiments typically measure properties of the B and D hadrons within jets [85][86][87][88].Theoretical studies of heavy quarks fragmenting in jets have been performed using effective field theories [89][90][91][92].
In hadronic collisions, one usually measures projections of the hadron momentum with respect to the jet one (see e.g.[85]) or the transverse momentum fraction of the heavy flavour with respect to the jet transverse momentum (see e.g.[86]): In all plots, we show our resummed result, as well as the one obtained with the Monte Carlo event generator HERWIG, both at parton and hadron level.Jets are selected with the anti-k t algorithm with R 0 = 0.4, and groomed with Soft Drop with z c = 0.1 and β = 0.The plots on the left are for p t ≥ 50 GeV, while the ones on the right are for p t ≥ 150 GeV.The uncertainty bands for the analytic predictions are obtained by varying the resummation scale by a factor of two above and below the hard scale p t R 0 , i.e. µ R ∈ p t R 0 2 , 2p t R 0 .In what follows, we would like to study the correlation, if any, between the variable ζ , which measures the departure from Born kinematics of the identified heavy flavour within a jet, to z g (with β = 0), which probes the kinematics of the splitting that passes Soft Drop, in a heavy-flavour tagged jet.
Let us start with the O(α S ) calculation, in the quasicollinear limit.The result for z g can be immediately read off Eq. ( 39) by setting β = 0: The computation of the momentum fraction ζ , in the same approximation, is rather straightforward.We consider the emission of a gluon with momentum fraction z in the quasi-collinear limit.If the gluon is emitted within the jet, then z = ζ , i.e. the b(c) quark has momentum fraction 1 − z = x.If the gluon is emitted outside the jet, it does not contribute to the observable and x remains equal to unity. 3 We have 1 where the second equation holds up to power corrections in the mass.Eq. ( 49) coincides with Eq. ( 48), for ζ > z c , i.e.
x < 1 − z c .We conclude that 1 − x and z g are fully correlated at this perturbative order.It is interesting to investigate whether this remains true beyond that.We study this problem by performing a Monte Carlo study.
Before doing so, we study the ζ = 1 − x spectrum using the same HERWIG simulation setup described for our previous analyses.The results are shown in Fig. 5, on the left, where we plot the ζ distribution at parton-and hadron-level, for b-jets with p t ≥ 150 GeV.We note that non-perturbative effects are large.This is interesting because, while fragmentation functions are under better perturbative control than Soft Drop observables, the latter ones seem to be more robust against non-perturbative corrections.
Finally, the right-hand plot of Fig. 5 shows the two dimensional distribution for z g (with z c = 0.1 and β = 0) and ζ , again for p t ≥ 150 GeV, at hadron-level.We note that the O(α S ) correlation between these two observables is not maintained when higher-order corrections and non-perturbative effects are included.Therefore, we conclude that z g and ζ = 1 − x can offer different handles to study heavy-flavour dynamics at the LHC.

Conclusions and Outlook
In this work, we have considered heavy-flavour jets, namely band c-jets.In particular, we have studied heavy-flavour jets in conjunction with the Soft Drop grooming algorithm.Exploiting resummed perturbation theory, we have computed the NLL distribution for the groomed jet radius θ g , for both band c-jets.Our calculation accounts for the dead-cone effect due to the heavy-quark mass, as well as m b and m c thresholds in the running coupling, which is frozen at a non-perturbative scale Λ = 1 GeV.We also discuss the role of clustering logarithms, first computed in [44].The presence of many scales originates many cases and regions to be considered.We have listed all of them and shown numerical results for representative cases.We have compared our findings for the θ g distribution with simulated data obtained with Monte Carlo event generators, assessing the role of nonperturbative corrections, such as hadronisation and the Underlying Event.We have found good agreement between resummed perturbation theory and Monte Carlo simulation, for the β = 0 case, even for relatively low values of the hard scale, i.e. p t = 50 GeV and R 0 = 0.4.The agreement worsens if positive values of the angular exponent β are considered.
We have also considered the momentum fraction z g of the first emission that passes Soft Drop.Our calculation generalises to the massive case the conditional-probability approach of [38] and it allows for a resummation of mass effects.While our results hold for any β ≥ 0, we have focussed our numerical investigations on the β = 0 case.This choice is motivated by a few reasons.First, this is the value for which, at least to lowest order, one has a clear factorisation of the conditional probability expression, leading to a distribution proportional to the splitting function.For β > 0, the dependence on the splitting function is smeared out.Second, it is for the β = 0 case that one expects to find more similarities to the fragmentation variable x, and third, from the study of the θ g distribution, we have found that the β = 0 is under better theoretical control and more sensitive to dead cone effects.We have compared our numerical results to Monte Carlo simulations.As for the massless case, the result is driven by the QCD splitting function and it is largely insensitive to non-perturbative effects.
Finally, we have compared the z g distribution to the fragmentation variable x, which measures the ratio of the heavy quark (hadron) transverse momentum to the jet p t .We performed the analytic calculation of both distributions at O(α S ), showing that they lead to the same results for z g = 1 − x > z c .We have investigated the all-order behaviour of these observables using Monte Carlo simulations and discovered that the parton-shower significantly dilutes this correlation.However, we have found that, in contrast to z g , the x receives sizeable non-perturbative corrections.In the future, it would be interesting to repeat these studies with the modified version of the declustering procedure [93], whereby one follows the flavoured branch.
In this work, we have performed a detailed theoretical study of the effect of Soft Drop grooming on heavy-flavour jets.Before being able to compare to experimental data, such as the ones collected by the ALICE collaboration [40], a few steps need to be taken.We are going to implement our calculation in the resummation plugin to SHERPA [94,95] in order to match our results to NLO theoretical predictions, with fiducial cuts, supplemented with non-perturbative corrections, as done for instance in [96][97][98][99][100][101].In this context, it would be also interesting to lift the small-z c limit and investigate, in the β = 0 case, flavour-changing contributions that may induce radiation into the dead-cone region.Finally, it would be important to improve the accuracy of the z g calculation by including the resummation of z g and z c logarithms, as done in [45].t > Λ 2 , the integral in Eq. (A.12) can be written as: with: .Thus, we have only 3 distinct subcases for (a) and (b), 2 for (c) and still 1 for (d), totalling 9 cases for b-jets.
In this work, we have fixed the value of z c = 0.1 and considered only two distinct values for the angular exponent β = 0, 1.We have also worked with just one value for the anti-k t jet radius, R 0 = 0.4, and fixed the non-perturbative scale Λ = 1 GeV.This results in only two different hierarchies to be considered.In fact, for p t = 150, 300 GeV we have:

Fig. 1 :
Fig.1: Lund plane representation of the soft and quasi-collinear phase space for emissions off a b quark, on the left, and off a c quark, on the right.The dead-cone region is indicated in grey.Soft Drop is applied and the groomed-away region appears in light purple.The vertical dashed line in black indicates a measurement of the groomed radius θ g and the corresponding area in pink is the vetoed region, which gives rise to the Sudakov form factor.The horizontal red lines indicate the boundaries between different flavour numbers for the running coupling, while the green one marks the boundary of the non-perturbative region.

)Fig. 2 :
Fig.2:The groomed jet radius (θ g ) distribution for b-jets (top), c-jets (middle) and light-quark jets (bottom).In all plots, we show our NLL resummation, as well as the results obtained with the Monte Carlo event generator HERWIG, both at parton and hadron level.Jets are selected with the anti-k t algorithm with R 0 = 0.4, with p t ≥ 150 GeV, and groomed with Soft Drop with z c = 0.1, and two values of the angular exponent: β = 0 (on the left) and β = 1 (on the right).The uncertainty bands for the analytic predictions are obtained by varying the resummation scale by a factor of two above and below the hard scale p t R 0 , i.e. µ R ∈ p t R 0 2 , 2p t R 0 .

Fig. 3 :
Fig.3: Ratios of groomed jet radius (θ g ) distribution of b jets to light-quark jets for three values of the jet transverse momentum: 50 GeV (top), 150 GeV (middle) and 300 GeV (bottom).The plots on the left are for β = 0, and the ones on the right are for β = 1.In each plot, we show our NLL resummation and the results obtained with HERWIG, both a parton-level and hadron-level.

Fig. 4 :
Fig.4: The groomed momentum fraction (z g ) distribution for b-jets (top), c-jets (middle) and light-quark jets (bottom).In all plots, we show our resummed result, as well as the one obtained with the Monte Carlo event generator HERWIG, both at parton and hadron level.Jets are selected with the anti-k t algorithm with R 0 = 0.4, and groomed with Soft Drop with z c = 0.1 and β = 0.The plots on the left are for p t ≥ 50 GeV, while the ones on the right are for p t ≥ 150 GeV.The uncertainty bands for the analytic predictions are obtained by varying the resummation scale by a factor of two above and below the hard scale p t R 0 , i.e. µ R ∈ p t R 0 2 , 2p t R 0 .

Fig. 5 :
Fig. 5: Monte Carlo studies of the fragmentation function variable x, obtained with HERWIG.The plot on the left shows the ζ = 1 − x distribution for b-jets, at parton-and hadron-level.The plot on the right shows the correlation between z g (with z c = 0.1, and β = 0) and ζ , at hadron level.In both cases, we have p t ≥ 150 GeV, and R 0 = 0.4.

,
defined above.In contrast, c can only assume the values: 1 1+β , − β 1+β , −1, 0 .As long as k 2 n f = 5, δ 54 for n f = 4, δ 54 + δ 43 for n f = 3. (A.15) is always the smallest one.We classify the different cases according to the value of z c .(a) The case z 2 c > ξ b has the following subcases: understood that the values in the table are ordered from big to small.(b) The case ξ b > z 2 c > ξ c has the same subcases as (a).(c) The case ξ c > z 2 c > ξ Λ has 3 subcases: Finally, for the last case, we can only have ξ Λ > z 2 have to consider 14 cases.This is what actually happens for light jets and c-jets.For b-jets, the situation is slightly simpler because the dead cone at θ 2 i = ξ b implies that we are not sensitive to the scale z 2 c ξ 1+β c