A consistent resummation of mass and soft logarithms in processes with heavy flavours

Perturbative calculations for processes that involve heavy flavours can be performed in two approaches: the massive scheme and the massless one. The former enables one to fully account for the heavy-quark kinematics, while the latter allows one to resum potentially-large mass logarithms. Furthermore, the two schemes can be combined to take advantage of the virtues of each of them. Both massive and massless calculations can be supplemented by soft-gluon resummation. However matching between massive and massless resummed calculations is difficult, essentially because of the non-commutativity of the soft and massless limits. In this paper, we develop a formalism to combine resummed massive and massless calculations. We obtain an all-order expression that consistently resums both mass and soft logarithms to next-to-leading logarithmic accuracy. We perform detailed calculations for the decay of the Higgs into a heavy-quark pair, and discuss the applications of this formalism to different processes.


Introduction
The physics of heavy flavours is especially important in particle physics phenomenology, for a number of reasons.The Higgs boson decays primarily into pairs of b quarks and, although this decay mode is challenging because of its large background, it plays a central role in studies of electro-weak symmetry breaking.Furthermore, many aspects of so-called flavour-physics can be scrutinised in the charm and bottom sectors.
Heavy flavours are also a valuable probe of strong interactions.Despite the fact that gluons couple to quarks irrespectively of their mass, quark masses do affect emergent phenomena, such as jet formation and their substructure.In this context, a noteworthy effect is the socalled dead-cone, i.e. the fact that colour radiation around heavy quarks is suppressed [1,2].Dedicated phenomenological strategies have been designed to expose and study this effect, e.g.[3][4][5], which has been recently measured by the ALICE collaboration at the LHC [6].Furthermore, the possibility of exploiting the imprint that quark masses leave on colour correlations has been recently investigated in the context of b-tagging [7].Moreover, while the top quark mass is so large that this particle's lifetime is shorter than the typical time scale of hadron formation, the mass of c and b, although in the perturbative regime, is not so large so that hadron formation occurs.Thus, processes involving c and b quarks can be exploited to scrutinise the mechanism that binds quarks and gluons into colour-neutral hadrons.Finally, studies of intrinsic heavy-flavour component (mostly charm quarks) in the proton requires precision calculations of perturbative cross-section involving heavy quarks.
Two main strategies to perform QCD calculations of observables with heavy flavours are usually employed.In the so-called massive scheme, heavy quarks in the final state are considered as real, on-shell, particles with a non-zero mass.The main advantage of the massive scheme is that the kinematics of the produced heavy flavour is treated exactly, because the full mass-dependence is retained.An important example of a calculation performed in this scheme is heavy flavour production in hadron-hadron collisions, which has been computed up to nextto-next-to-leading order (NNLO).[8][9][10] The fixed-order precision can be improved by various types of all-order calculations, e.g.soft-gluon-resummation, for both inclusive production cross sections and differential distributions [11,12], high-energy resummation [13,14] or even transverse momentum resummation [15].The inclusion of heavy-quark effects in generalpurpose Monte Carlo parton shower codes is also an area of active research, see e.g.[16][17][18][19].
The range of energies probed by collider experiments is typically much larger than the heavy-flavour mass, making heavy-flavour production a multi-scale problem.Theoretical predictions for these processes, even for inclusive observables, are plagued by logarithms of m 2 /q 2 , where q 2 the square of the hard scale, that can spoil the convergence of the perturbative expansion.Therefore an alternative calculational framework is often employed.This second approach exploits fragmentation functions to resum these logarithmic corrections to all orders.This is possible because these logarithmic corrections are related to collinear dynamics, which would give rise to divergencies in a massless theory.It follows that, up to corrections O m 2 /q 2 , heavy-flavour production cross-sections obey a factorisation theorem and can be written as the convolution of process-dependent partonic (massless) coefficient functions with universal heavy-flavour fragmentation functions.Fragmentation functions obey DGLAP evolutions equations (with time-like splitting functions), which allow one to resum the large logarithmic corrections we are discussing, in analogy to the initial-state collinear factorisation theorem.The initial-condition for heavy quark fragmentation functions can be computed in perturbation theory, as originally pointed out in Ref. [20,21], where a NLO computation is presented.The NNLO corrections were computed later in Refs.[22,23].The initial condition of the evolution is, by definition, free of mass logarithms, but it is affected by soft logarithms, that can be resummed to all orders too.[24][25][26][27] Virtues of the massive and massless schemes can be combined together by matching the fixed-order calculation performed in the massive scheme, with the all-order resummation of mass logarithms achieved by the fragmentation function approach [28][29][30][31][32][33][34].Since soft-gluon resummation is available for both the massive and the massless schemes, it is natural to explore the possibility of matching the two calculations at the resummed level, obtaining a theoretical prediction that fully accounts for all mass and soft logarithms.As it has been pointed out in the literature, see e.g.[12,[35][36][37][38][39] this is far from trivial, because the structure of soft logarithms significantly differs in the two approaches, hampering the construction of an all-order matching scheme.The goal of this paper is to overcome this difficulty and build a matching scheme that allows us to consistently resum both mass and soft logarithms in processes with heavy quarks.The crucial ingredient for the construction of a consistent resummation formula is a modification of the standard fragmentation-function calculation that, by exploiting the so-called quasi-collinear limit [40,41], correctly accounts for mass effects that originate from QCD radiation from all (hard) heavy-quarks in the process.Furthermore, we will also discuss in some details the inclusion of heavy-quark threshold effects that arise from the treatment of the QCD running coupling in a decoupling scheme, which is more natural than standard MS when considering processes with heavy flavours.In this work, we concentrate on parton-level results, leaving detailed phenomenological analyses to future work.
In our discussion, we are going to focus on the decay of an (off-shell) electroweak boson into two massive quarks, considering H → b b as a concrete example.However, because our construction is general, we will briefly discuss how to apply it to other processes involving heavy flavours, such as deep-inelastic scattering and heavy-quark decay.For the latter, we will also comment on similarities and differences between our calculation and an alternative approach developed in Refs.[38,39].
This paper is organised as follows.In Sect. 2 we review known results on soft-gluon resummation in the massive and massless schemes, highlighting the issues that we are set to solve.In Sect. 3 we revisit the problem in momentum space.This will allow us to better identify the relevant kinematic regions, leading to a more refined calculation.We will present our results for H → b b in Sect.4, comment on other processes and approaches in Sect.5, before drawing our conclusions in Sect.6.Details of the calculations and explicit results are collected in the appendices.

Heavy flavour pair production in weak decays
In this section we consider the decay of a colour-singlet massive system H, for instance an offshell photon, a Higgs boson or a Z boson, into a heavy quark-antiquark pair plus undetected radiation: (four-momenta are indicated in brackets.)In what follows, it is understood that the heavy quark is a b quark, but similar considerations can be extended to charm.We are interested in the differential decay rate dΓ dx with respect to the dimensionless variable which coincides with the fraction of the total available energy carried away by the b quark in the centre of mass frame.At lowest order in QCD perturbation theory, when no radiation is present, energy-momentum conservation gives x = 1: At higher orders the x dependence of the differential rate becomes non-trivial, and the limit x → 1 corresponds to kinematical configurations in which the emitted radiation is either soft, or collinear to the final-state antiquark.
The calculation of the differential rate can be carried out in two different schemes, sometimes referred to as the massive (or 4-flavour) and the massless (or 5-flavour) schemes.In the first approach, the spectrum dΓ dx is computed to some finite order in perturbation theory taking into account the finite value of the heavy quark mass m exactly.In this approach, collinear singularities are regularised by the heavy quark mass and, consequently, powers of log m 2 q 2 appear in the perturbative coefficients; such logarithms are large as q 2 ≫ m 2 , and may eventually spoil the convergence of the perturbative series.On the other hand, the kinematics of radiation is exactly taken into account in the whole range, up to a given order in α s .In the second approach, all flavours, including the b quark, are treated as massless, which is an accurate approximation as long as q 2 ≫ m 2 .Heavy-quark production is described in the fragmentation function formalism: factorisation is exploited in the sense that each final-state parton can fragment into a heavy quark with a suitable probability distribution, in close analogy to what happens for initial-state parton distribution functions.The mechanism of factorisation of collinear singularities is also very similar.In this approach, the heavy flavour mass is retained only as a regulator of collinear divergences, while contributions proportional to powers of m 2 q 2 are systematically neglected.Within this framework, the differential decay rate takes the factorized form where µ 2 F and µ 2 R are the factorisation and renormalisation scales, typically chosen of the order of the hard scale q 2 .The sum runs over all partons; the functions C i are processdependent partonic cross sections, which admit a perturbative expansion in powers of α s as long as q 2 is large enough.They are convoluted with the process-independent fragmentation functions D i .In this approach, only the dominant collinear region of radiation is included (as opposed to the massive scheme, where large-angle radiation is perturbatively taken into account).On the other hand, collinear logarithms of m 2 q 2 are resummed to all orders, up to a given logarithmic accuracy, through the solution of DGLAP equations for the fragmentation functions.Furthermore, the initial condition for the evolution equation of the b fragmentation function is given at an initial scale of the order of the b mass, and can therefore can be computed perturbatively [20,22,24].
The convolution product in Eq. (2.4) is turned into an ordinary product by Mellin transformation, defined as for a generic function f (x) defined in the range 0 ≤ x ≤ 1.We find where we have defined ξ = m 2 q 2 .The scale dependence of the fragmentation functions is governed by the DGLAP evolution equations The anomalous dimensions γ ij , Mellin transforms of the time-like splitting function, can be computed perturbatively; their explicit expressions to order α s can be found e.g. in [42].The solution of Eq. (2.7) can be schematically written in terms of a matrix of evolution kernels E ij and initial conditions, given at some reference scale µ 0F : (2.8) In the following, we will denote by Γ k (N, ξ) the (Mellin transformed) spectrum computed in the massive scheme up to order k, and by Γ (5) ℓ (N, ξ) the same quantity computed in the massless schemes, with evolution equations for the fragmentation functions solved at next ℓ -toleading logarithmic accuracy, and both the coefficient functions C i and the initial conditions for the evolution of fragmentation functions D 0 j computed up to order k.One may take advantage of both calculation schemes by combining the two results: where "double counting" stands for the perturbative expansion of Γ ℓ (N, ξ) to order k.In the following we will restrict ourselves to the case ℓ = 1, which is referred to as the FONLL scheme [28].
The perturbative coefficients of both quantities appearing in the rhs of Eq. (2.9) display a logarithmically divergent behaviour as N → ∞.Such logarithms arise as a remnant of the cancellation of soft singularities, which results in the presence of distributions in the physical spectrum, which dominate the perturbative coefficients in the vicinity of the threshold region x → 1. Mellin transformation maps the large-x region into the large-N region; it can be shown that at large N the Mellin transform of D k (x) is a polynomial of degree k in log N .These logarithmic contributions can be summed to all orders in the QCD expansion, up to a given logarithmic accuracy.The quantities appearing in Eq. (2.9) have different structures at large N .The order-α n s perturbative coefficient in the expansion of Γ k (N, ξ) is a polynomial of degree n in log N , plus terms that vanish as N → +∞.The coefficients of the polynomial carry the full ξ dependence: they include both powers of ξ, that vanish in the massless limit, and powers of log ξ, which regularize collinear singularities.The coefficient functions ℓ (N, ξ) are ordinary, subtracted partonic cross sections in massless QCD: the perturbative coefficients contain up to two powers of log N for each power of α s , corresponding to the remnants of soft singularity cancellation and collinear singularity subtraction.Finally, the initial conditions for the evolution of the fragmentation functions D 0 j (N, µ 2 0F , m 2 ) can be computed perturbatively; the corresponding order-α n s coefficients are polynomials of degree 2n in log N , whose coefficients depend on log ξ but not on powers of ξ.
As mentioned in the introduction, the resummation of soft-gluon contributions at all orders, up to a given logarithmic accuracy, for both Γ (4) k (N, ξ) and Γ (5) ℓ (N, ξ) can be performed, and has been studied in the in the past.In particular, for the massive-scheme case, we have where Γ (4,res) ℓ 1 (N, ξ) is the massive-scheme spectrum resummed to next ℓ 1 -to-leading log N , and the "double counting" term is its expansion to order k.Γ (4,res) ℓ 1 (N, ξ) was computed in [12] for ℓ 1 = 2. Similarly Γ (5) (N, ξ) = Γ (N, ξ) is the massless-scheme spectrum with both the coefficient functions and the initial conditions resummed to next ℓ 2 -to-leading log N , and the "double counting" term is its expansion to order k.
The combination of Eq. (2.11) and Eq.(2.12) in a single formula, which would incorporate the best of theoretical knowledge of the b energy spectrum in perturbative QCD, is not straightforward.Indeed, the sum of the results in Eq. (2.11) and Eq.(2.12) would require a further subtraction of doubly-counted contributions at all orders in α s ; such subtraction is far from trivial, because the two limits N → ∞ and ξ → 0 are known not to commute [12,[35][36][37][38][39]: one cannot simply subtract from the Γ (4,res) ℓ 1 (N, ξ) term in Eq. (2.11) its expression in the ξ → 0 limit, because it differs from the N → ∞ limit of Γ (5,res) ℓℓ 2 (N, ξ).In the following, after reviewing the expressions for the resummation of threshold logarithms in both the massive and the massless scheme, we will present a solution to this problem, valid to next-to-leading logarithmic accuracy.

Soft resummation in the massive scheme
The calculation of Γ (4,res) ℓ 1 (N, ξ) was performed in Ref. [43].It takes the exponentiated form where K(ξ, α s ) is a process-dependent factor, and γ soft , the so-called massive soft anomalous dimension, admits an expansion in powers of α s : where An explicit calculation gives Note that the resummed expression in Eq. (2.13) features, at most, single logarithms of N to any order in perturbation theory, i.e. α n s log m N , with m ≤ n.This is not surprising: collinear singularities are absent because of the finite quark mass and, consequently, collinear logarithms do not appear as logarithms of N but rather as logarithms of the mass.Thus, knowledge of the massive soft anomalous dimension to order ℓ 1 + 1 allows us to perform soft resummation at N ℓ 1 LL accuracy.The soft anomalous dimension is currently known to second order [44][45][46][47].However, in this work we restrict ourselves to resummation to singlelogarithmic accuracy, both in the massless and in the massive schemes, so we consider ℓ 1 = 0, i.e. we take the massive soft anomalous dimension at one loop, and the running of the strong coupling in Eq. (2.13) is taken into account at the leading logarithmic level.
At this accuracy, the Mellin transform Eq. (2.13) can be performed [48] by the replacement with N = N e γ E , γ E the Euler-Mascheroni constant and Θ the Heaviside step function.Using Eq. (2.17), after the change of integration variable z = 1 − x (which we will often adopt in the following) we obtain Γ (4,res) .18)In this expression, inverse powers of N are systematically neglected, while the dependence on m q is taken into account exactly.For this reason, Γ (4,res) ℓ 1 (N, ξ) is only accurate for values of N such that1 N < m q .As a consequence, the heavy quark threshold, corresponding to z = m q , lies within the integration range in the exponent of Eq. (2.18).In the variable-flavour number scheme for the running coupling, the number of active flavours that contribute to the QCD β function is increased by one unit at each quark mass threshold.Hence, the running coupling in Eq. (2.18) reads where It is interesting to study the small-ξ behaviour of Eq. (2.16).We find Thus, a class of collinear logarithms, specifically those that accompany a soft logarithm, are also resummed by Eq. (2.18).However, a full resummation of all mass logarithms can only be performed in a different calculation scheme, as detailed in the next section.For later convenience, we introduce a subtracted massive soft anomalous dimension: which vanishes as ξ → 0, β → 1.
Terms proportional to powers of log ξ with N -independent coefficients appear in the prefactor K.In the small-ξ limit we find 1 As pointed out in Ref. [12], this is not what one would expect: at order α k s , the leading collinear logarithm should be proportional to log k ξ.It was shown in Ref. [12] that the mass appearing in the double logarithm contribution at order α s is the mass of the undetected massive parton, b in the present case.This suggests that, conversely, spurious double logarithms of N from radiation off the undetected parton may arise in the 5-flavour scheme calculation.In the next section, we will confirm this expectation.

Soft resummation in the massless scheme
The resummation of soft logarithms for the b energy spectrum in the massless (5-flavour) approach was performed in Ref. [24] and the specific case of H → b b was considered in Ref. [49].Logarithms of the Mellin variable N appear both in the coefficient functions and in the initial condition for fragmentation functions into heavy quarks, which are typically given at an energy scale of the order of the heavy quark mass.The general expression for the resummed energy spectrum has the form where C is the coefficient function, E the DGLAP evolution kernel for the fragmentation function, and D 0 the initial condition.For future convenience, it will be useful to separate off the logarithmically enhanced contribution to E from the regular contribution at large N .
To this purpose, we define a subtracted evolution kernel where The renormalisation and factorisation scales µ 2 R , µ 2 F are chosen to be of the order of magnitude of q 2 , while the corresponding reference values µ 2 0R , µ 2 0F are of order m 2 .We remind the reader that knowledge of the DGLAP evolution kernel at the ℓ-loop order allows one to resum collinear logarithms to N ℓ LL.In this work we consider the case ℓ = 1.At NLL in N , i.e. ℓ 2 = 1, the resummed decay rate can be written The factor exp(J) in Eq. (2.26) describes soft radiation emitted collinearly to the observed b-quark.The jet function J has the following form: where Finally, the jet function J is related to radiation emitted collinearly to the b direction.One finds Note that beyond NLL the resummed exponent can no longer be written as the sum of two independent jet functions.The functions A, B and H have an expansion in powers of α s : with where n.The coefficients given in Eq. (2.32) are sufficient to achieve NLL accuracy.Explicit expressions for the evolution factor E, the coefficient function C and the initial condition for the fragmentation function D 0 at the relevant level of accuracy are given in App. A.
We now investigate the presence of double logarithmic terms, i.e. tems of order α k s log k+1 N , in the exponent of the resummed expression Eq. (2.26).Leading logarithmic contributions to the exponent J originate from the first term in the expansion of the function A(α s ) in Eqs.(2.25,2.28,2.29).We find where we have not indicated the dependence on α s and renormalisation and factorisation scales to simplify notations.Let us first evaluate the expression above at fixed coupling.We find which shows that here are no double logarithms of N at this level of accuracy.An analogous calculation for the jet function related to the unmeasured b leads to In this case double logarithms of N appear in the J contribution to the exponent already at order α s .We note that, while the calculation of D 0 , which contributes to the measured jet function J, is performed in the so-called quasi-collinear limit [40,41], the calculation of J is performed in the pure collinear limit.We will expand on this consideration in the following section, showing that the treatment of J and J on equal footing will allow us to consistently match Γ (4,res) ℓ=1,ℓ 2 =1 to all orders.The logarithmic structure of Eq. (2.33) is slightly modified when the running of the strong coupling is taken into account.The integration region in Eq. (2.33) is mapped into a rectangle The integration domain is divided in two regions by the curve q2 z 2 w = m 2 , or w = ξ z 2 ; in the lower region α s = α (4) s , while in the upper region α s = α (5) s (see Fig. 1).We must distinguish two cases.If ξ N 2 > 1 (green curve in Fig. 1), then 5)  s (q 2 z 2 w) + NLL. (2.37) The terms in the last line are N -independent.Using the leading logarithmic solution for α we get which does not contain any leading logarithmic term, i.e. tems of order α k s log k+1 N .
On the other hand, if ξ N 2 < 1 (red curve in Fig. 1), then dw w α (5)  s (q 2 z 2 w) + NLL.(2.40)In this case, terms of order α k s log k+1 N cancel at order α s , but not at higher orders.Indeed, using the leading logarithmic solution for α (n) s (µ 2 ) with the threshold condition and expanding in powers of α s (m 2 ) one finds Thus, we have found that the jet function J also exhibits leading logarithms of N , because of running coupling effects, but only in the region N < 1/ √ ξ.We note that the above results refer to the case of real and positive Mellin moments.When the resummed spectrum in Mellin space is continued analytically to the complex N plane, which is mandatory for Mellin inversion, their interpretation is not straightforward, especially as far as the hierarchy between different energy scales is concerned.A better understanding of this issue can be obtained by looking first at resummation in momentum space, as we show in the next section.

Resummation in momentum space
Let us briefly summarise how to perform NLL resummation in momentum space.[51] In this work we are interested in the energy spectrum of a final-state heavy quark, but the technique we are illustrating applies to a more general class of observables.We therefore consider of a generic infrared and collinear (IRC) safe observable V, a positive function of finalstate momenta that vanishes at Born level.It proves convenient to consider the normalised cumulative distribution i.e. the probability for the observable V to be smaller than some given value v, rather than the differential rate dΓ dv .Because we are ultimately interested in H → b b, we focus on a single quark-antiquark hard dipole with centre-of-mass momenta p 1 and p 2 respectively.We first address the case of massless quarks, p 2 1 = p 2 2 = 0, and then extend our results to massive quarks.
We begin by considering the O(α s ) contribution to Σ(v), corresponding to one-gluon emission (plus O(α s ) virtual corrections.)We denote by k the momentum of the emitted gluon.Neglecting for the moment the contribution of virtual corrections, the normalised distribution at order α s can be schematically written where dΓ (1) is the differential width for one-gluon emission.Therefore In the massless case, QCD emission probabilities, in the soft and collinear limits, behave as dkt kt dη i , where k t and η i are the transverse momentum and rapidity (assumed positive) of the emission with respect to emitting particle momentum direction: It is therefore convenient to express the kinematics in terms of these variables.The boundary of the phase space in terms of k t , η i is determined by the condition In our case, logarithmically enhanced contributions to the cumulative distribution originate from gluon emission either collinear to p 1 or to p 2 .We therefore define two regions of collinear emission: in terms of some large rapidity value η cut .The condition (3.7) implies as a consequence of Eq. (3.5).For large values of η cut we have (3.9) Restricting the integration region to the collinear-emission regions of the phase space, we have where M is the appropriate invariant amplitude, and Q is the maximum value allowed for the transverse momentum, It is shown in Ref. [51] that in most cases the observables of interest behave, in the soft and collinear limits, as where we have chosen Q as reference (hard) scale, and the index i labels the different collinear regions (in our case, the two regions collinear to either ⃗ p 1 or ⃗ p 2 .)The positive constants a i , b i , d i depend in general on the particular region of the phase space where the collinear limit is taken; in our case, we have two sets of constants.In the region collinear to the quark momentum p 1 we have, in the large rapidity limit, Similarly, in the region collinear to the antiquark momentum The squared amplitudes that appear in Eq. (3.10) can be computed, in the collinear limits, by exploiting factorization of collinear singularities.Let us first consider the region where the gluon is emitted in the direction of the quark momentum.It is convenient to adopt a Sudakov parametrisation of the gluon momentum: where p1 = (p 0 1 , −⃗ p 1 ) and κ is a spacelike vector such that p 1 • κ = p1 • κ = 0, and κ 2 = −k 2 t .The coefficient z1 is suppressed in the collinear limit.Indeed, the gluon mass-shell condition gives is much smaller than z 1 in the small-k t limit.In this limit, z 1 represents the fraction of the quark energy which is carried away by the gluon.The gluon rapidity with respect to the quark momentum direction is given by Conversely, when the emission is collinear to the antiquark, we may write the gluon momentum as and Thus, we have where is the appropriate timelike DGLAP splitting function.The −1 contributions inside square brackets in Eq. (3.22) account for virtual corrections, and IRC safety of the observable ensures the cancellation of the z i → 0 and k t → 0 singularities.Eq. (3.22) can therefore be written The choice of the reference hard scale of the process Q, and therefore of η cut , is immaterial at NLL [51].Following the literature we choose η cut = 0. Note that with this choice the upper bound for k t is larger than its kinematic limit q/2.However, this phase-space region is outside the jurisdiction of the resummed calculation and proper energy-momentum conservation is usually restored when matching to a fixed-order calculation.At LL accuracy, where each emission comes with a maximal number of logarithms, one can further assume that a single emission strongly dominates the value of the observable, leading to the exponentiation of the one-loop calculation also in momentum space (see for instance [52]): where the Sudakov form factor R(v), which represents the no-emission probability, is expressed as the sum of two jet functions, one for each collinear sector: with In order to achieve full NLL accuracy, Eq. (3.25) must be corrected in two ways.First, one should consider running of the strong coupling in Eq. (3.27) at two loops in the so-called Catani-Marchesini-Webber (CMW) scheme: [53] where α s (k 2 t ) is in the decoupling scheme, see Eq. (2.19).Second, we can no longer work in the strongly-ordered approximation, and the resummation must be performed either with numerical methods [51] or in a conjugate (e.g.Mellin) space, in order to factorise the observable definition.In the latter case, at the end of the calculation, the result should then be brought back to physical space.In some cases, this inversion can be done, to a given logarithmic accuracy, in closed-form, resulting in a correction factor, which only depends on derivatives of the Sudakov form factor. [54] However, we note that Eq. (3.25), supplemented by the prescription in Eq. (3.28), is enough to achieve NLL accuracy in Mellin space, provided that we identify v = N −1 , as shown in App.B. So finally where the replacement (3.28) in the jet functions j and j is understood.
The formalism described above can be generalised to the case of massive quarks.When m ̸ = 0, kinematics undergoes some minor modifications.The boundaries of the phase space gets slightly modified: the inequality (3.5) becomes Furthermore, the relationship between z i and zi , which originates from the mass-shell condition for the gluon momentum, is now given by zi Finally, the expressions Eqs.(3.18,3.21)for the gluon rapidity in the collinear regions still hold, up to corrections of higher order in m q and kt q .It was realised long ago that squared QCD matrix elements with massive partons factorise in the so-called quasi-collinear limit, [40,41] in which both k 2 t q 2 and m 2 q 2 go to zero, with their ratio kept constant.In this limit, the squared invariant amplitude for one-gluon emission takes the form The mass-dependent shift in the denominator of Eq. (3.33) acts as an effective lower bound of a logarithmic k t integration: Using Eqs.(3.18) and (3.21), we obtain η i < − 1 2 log ξ.This is the well-known dead-cone effect, [1,2] i.e. the fact that radiation off massive partons at angles below m/q is not logarithmic enhanced.Finally, we note that now the observable V can now explicitly depend on the heavy-flavour mass through ξ.

Lund plane interpretation
The procedure illustrated in the previous subsection is conveniently interpreted in terms of Lund diagrams.Lund diagrams offer a graphical representation of the kinematics described above.The Lund plane is usually represented in terms of pairs of the aforementioned logarithmic variables, so that, in the soft and collinear limit, equal areas correspond to equal emission probabilities.As an example, we show in Fig. 2, on the left, the Lund plane in case of a massless dipole in the dipole rest frame.
The variable in the vertical axis is the logarithm of The variable on the horizontal axis in the right half-plane is the rapidity η 2 of the emission with respect to the antiquark momentum, while on the left we report the rapidity η 1 with respect to the quark momentum.Note that where k ℓ is the longitudinal component of the gluon momentum along the direction of parton i, and θ i the angle between the gluon and the parton momenta.Hence, in the collinear limits cos θ i → 1 This is explicitly indicated on the horizontal axis in Fig. 2. In the massless case the boundaries of the phase space, namely the lines z 1 = 1 and z 2 = 1, are represented on the Lund plane by the lines log Lines of constant V in the soft and collinear limits, when Eq. (3.12) applies, are represented on the Lund plane by straight lines, as shown in Fig. 2.
The cumulative distribution also has a graphic interpretation on the Lund plane.The real emission cannot take place everywhere because it must give a contribution to the observable smaller than v. Thus, the emission probability corresponding to the phase-space region above the lines V 1 = v and V 2 = v, which is indicated with the shaded area in Fig. 2, is subtracted from 1 in Eq. (3.24).Thus, in order to compute the Sudakov form factor R, which represent the no-emission probability, we have to integrate the splitting functions (with running coupling at the appropriate accuracy, as previously discussed) over the shaded areas in in Fig. 2.
In the massive case, we have to employ the quasi-collinear limits of matrix elements squared, Eq. (3.33).This would seem to spoil the Lund plane interpretation.However, as noticed in Eq. (3.35), from the point of view of the logarithmic structure, we can still use the massless phase-space (but with quasi-collinear splitting functions) and the quark mass acts as a new phase-space boundary.On the Lund plane, this allowed region is bounded by the two vertical lines at η i = − 1 2 log ξ.For rapidities between these values, the heavy quark mass can be neglected altogether, while in the collinear regions η i > − 1 2 log ξ the quark mass acts as a cut-off for collinear singularities.This is shown in the right panel of Fig. 2. It is also useful to mark on the Lund plane the boundary between the regions of phase-space with 4 or 5 active flavours.This transition takes place when We have already noticed that, because the observable V can depend explicitly on the heavy-quark mass, we do not expect V 1 (k t , η 1 , ξ) and V 2 (k t , η 2 , ξ) to be of the same form as Eq.(3.12).Consequently, on the Lund plane, lines of constant 1−x are deformed with respect to the massless case.However, we are still allowed to represent them as solid-dashed straight lines on the Lund plane, provided that we take into account the boundaries at η i = − 1 2 log ξ, as shown in the right panel of Fig. 2.

Resummation of the cumulative distribution
In this section, we specialise the above discussion to the case of the observable V = 1 − x, with x defined in Eq. (2.2).In order obtain a resummed result for the cumulative distribution Σ(1 − x), we must discuss the parametrisations V i , i = 1, 2 of the observable V = 1 − x.Using energy-momentum conservation, we find which means that a measurement of 1−x is equivalent to a measurement of the invariant mass of the bg system recoiling against the b quark.The Sudakov parametrisations (3.15,3.19)give emission collinear to the quark, (3.39) emission collinear to the antiquark.
As we have anticipated, because of the presence of the mass term ξ, the above parametrisation has not the same form as Eq.(3.12).Consequently, on the Lund plane, lines of constant 1 − x would be deformed with respect to the massless case.However, to the accuracy we are working at, we are allowed to shift the transverse momentum as described in Eq. (3.35), provided that we take into account the boundaries imposed by the dead-cone.Thus, to NLL, we can use the following parametrisation: where we have used Eqs.(3.18) and (3.21).Note that this parametrisation coincides with the one one would find in the massless case.Therefore, to the accuracy we are working at, we can represent lines of constant 1 − x on the Lund plane as straight lines, parallel to the line z 1 = 1.
In order to obtain the resummed cumulative distribution, we compute the Sudakov form factor Eq. (3.27).In particular, we have to evaluate the measured jet function j and the unmeasured (or recoil) jet function j in momentum space.Using the shift Eq.(3.35) and the parametrisation (3.41), we have (3.43) Note that, in principle, we should also shift the argument of the running coupling.However, as shown in App.C, this effect only contributes beyond NLL accuracy.
In the following sections, we will express the jet functions as integrals over the running coupling.When this integrals are computed explicitly, both jets functions acquire an explicit dependence on the renormalisation scales µ 0R , µ R , for the 4-and 5-flavour components, respectively.In addition, the measured jet function also depends on the factorisation scales µ 0F , µ F .This dependence is omitted in Sects.3.3.1 and 3.3.2,but reinstated when we present our final results in App.D.
As a final remark, we note that, strictly speaking, the observable 1 − x is not IRC safe, because emissions that are arbitrarily collinear to the quark carry away some energy.However, because we are considering emissions off a massive quark, collinear singularities do not appear, and are replaced by large logarithmic mass corrections.With the formalism discussed in this section, we can only resum mass logarithms that have a non-vanishing coefficient when x → 1; the resummation of mass logarithms is achieved by DGLAP evolution, discussed in Sect.2.2.

Calculation of the measured jet function
In this section we compute the measured jet function in momentum space j(v, ξ) for v = 1−x, and we compare it to the analogous quantity J(N, ξ) computed in Mellin space in Sect.2.2.The Lund plane representation, depicted in Figs.3,4,5, is useful to identify the different integration regions.
The relevant integration region for the measured jet function is the area shaded in red in Figs.3,4,5, which refer to the cases 1 − x > √ ξ, √ ξ < 1 − x < ξ, 1 − x > ξ respectively.We see that the jet function receives contributions from both the regions below (n f = 4) and above (n f = 5) the red dashed lines which marks the heavy quark threshold.For values of 1 − x smaller than √ ξ the n f = 5 region freezes at 1 − x = √ ξ and the entire x dependence is given by the 4-flavour contribution.
We decompose j as and we consider first the case 1 − x > √ ξ.The k t integration range splits into two regions, Figure 3: Lund plane representation of the jet functions.The diagonal dashed lines in black represent the value of the observable 1 − x, which is taken to be larger than √ ξ.The shaded areas indicate the regions of phase space where emissions are vetoed, giving rise to the Sudakov form factor (in red the measured jet function j, in blue the recoil jet function j.)The corresponding Mellin-space jet functions, J and J are also indicated.separated by the red dashed line in Fig. 3, corresponding to the different numbers of active flavours.After translating the integration domain from the variables of the Lund representation to the pair z 1 , k 2 t , we obtain Note that in the first line of Eq. (3.45), which refers to the the n f = 5 region, we have neglected mass corrections.
Using the definition of the CMW strong coupling, Eq. (3.28), and of the relevant splitting functions, Eqs.(3.23,3.34), it can be checked that this result coincides with J(N, ξ) given in Eq. (2.27), provided that we identify 1/ N = 1 − x, µ 2 0 = m 2 , and µ 2 = q 2 .Thus, to NLL we have: In particular, the contribution from the n f = 4 region, i.e. the second line of Eq. (3.45), corresponds to the initial condition D 0 for the fragmentation function, while the n f = 5 contribution, i.e. the first line of Eq. (3.45), corresponds to the NLL contributions to E + ∆.As a side comment, it is interesting to note that the fragmentation function evolved up to the hard scale, i.e. the factor exp(E + D 0 ), can be obtained from Eq. (3.42), provided the integration is extended to the region bounded by the lines 0, thereby including not only the region where the emission is collinear to the antiquark momentum (shaded in blue in Fig. 3), but also part of the region outside the kinematic limit z 2 = 1.The role of the factor exp ∆ is therefore to subtract off this unphysical behaviour of the soft part of the fragmentation function. 3his calculation confirms the result of Eq. (2.40), namely the fact that the leading logarithms originating from the second and the third terms in Eq. (3.45) do cancel at first order but fail to do so at higher orders because the running coupling receives contributions from different flavour regions, i.e. from above and below the k 2 t q 2 = ξ line on the Lund plane.Furthermore, the analysis in momentum space allows us to identify the region in which this happens in terms of physical quantities, 1 − x > √ ξ, as opposed to ξ N 2 < 1 found in Eq. (2.40).We now turn to the case 1 − x < √ ξ, illustrated in Fig. 4. We find As anticipated, the n f = 5 contribution (first line) is now evaluated at the transition 1 − x = √ ξ and bears no x dependence.Furthermore, by inspecting the second line, we confirm that leading (double) logarithms in 1 − x cancel between the two contributions to all orders, as already found in the Mellin space computation, Eq. (2.37).
We also note that, in this case, double logarithms of the mass appear.Similarly to what happens in the region 1 − x > √ ξ for the double logarithms of 1 − x (or N ), this contribution cancels at O (α s ) but not at higher orders, because of the different evolution of the running coupling in different flavour regions.Finally, we note that there is no transition for j(1 − x, ξ) associated to 1 − x = ξ.As in the previous case, the result in Mellin space in the region 1 − x < √ ξ to NLL accuracy is given by These results are also in agreement with the fragmentation function analysis of Section 2.2.

Calculation of the unmeasured jet function
In this section we focus on the computation of j(1 − x, ξ) given in Eq (3.43).In this case there is a second transition at 1 − x = ξ, which was not present in the case of the measured jet function calculation.This transition is of kinematic origin and can be understood already at fixed coupling.With this assumption, the integrals in Eq. (3.43) are straightforward and Figure 5: Same as Fig. 3, for 1 − x < ξ.
the result reads Remarkably, after identifying 1−x = 1 N , the two expressions above reproduce the logarithmic terms present in the expansion at O (α s ) of the resummed formula in massless and massive schemes, Eqs.(2.23) and (2.18), respectively.Specifically, when 1 − x > ξ, j reproduces the double logarithmic behaviour of Eq. (2.23) at O (α s ).On the other hand, j captures both the single logarithms of Eq. (2.18) and the mass double logarithms of Eq. (2.22) when 1 − x < ξ.We also notice that at the transition point 1 − x = ξ, j is continuous.Therefore, the momentum space calculation of j allows us to better understand and resolve the issue of the non-commutativity of the x → 1 and ξ → 0 limits.Indeed, by computing the recoiling jet function in the quasi-collinear limit, rather than in the massless collinear limit (as done, for instance, in [24]), we are able to obtain an expression that interpolates between Eqs. (2.23) and (2.18).Furthermore, we are also able to confirm that the double logarithmic terms in Eq. (2.22) are entirely generated by j and the mass that appears is the one of the undetected antiquark b.Note that this effect is entirely driven by the behaviour of the observable in the hemisphere collinear to the antiquark, and it is absent on the other side of the Lund plane.
The above considerations can be generalised to NLL by including running coupling cor-rections (see also [55,56]).We have If 1 − x > √ ξ (see Fig. 3) the recoil jet function is computed with n f = 5 active flavours, and the quark mass neglected.Hence Thus, the calculation of j coincides with the one performed in Mellin space in the massless case, Eq. (2.30), after the replacement 1 where the first term is Eq.(3.51), evaluated at the transition 1 − x = √ ξ.In this case we found it convenient to swap the order of the integrations in z 2 and k 2 t .In this way, the z 2 integrals are all of the same form; the result, to NLL accuracy, has the form c 0 + c 1 log with c 0 and c 1 independent of k 2 t , but dependent on x, q 2 , ξ through the integration bounds.We note that the first and second integrals receive logarithmic contributions only by the softcollinear part of the splitting function, i.e. we can take in those integrals P gq = 2 z 2 , while in the third one, the upper limit of the z 2 integration is 1 and therefore the complete massless splitting function is needed.Finally, we must consider the region 1 − x < ξ.We find The first contribution does not depend on x and it is given by the result of the previous region, Eq. (3.53), evaluated at 1 − x = ξ.The x-dependent contributions arise from the n f = 4 region, as shown in Fig. 5.This term is sensitive to finite mass corrections, as discussed at the beginning of this section.In order to perform the integrals we adopted the same strategy explained above.Mellin space results in each region are obtained by the procedure adopted in the previous cases: Explicit results are given in App.D.

All-order matching and numerical results
We can now exploit the results of the previous section to combine, in a consistent way, the NLL resummed result differential rate for the process in Eq. (2.1) in the fragmentationfunction approach, Γ (5,res) ℓ=1,ℓ 2 =1 , with the resummed calculation Γ (4,res) in the massive scheme.To this purpose, we must take into account that constant (i.e.N independent) terms are not included in the momentum space results; they must therefore be recovered by matching with the calculations presented in Section 2. We have with Γ (1) (N, ξ) = Γ (5,res,sub) ℓ=1,ℓ 2 =1 exp j(1) (N, ξ) + j(1) (N ) , In the three expressions above, the suffix sub denotes the subtracted version of the 5-and 4-flavour calculations, that is, the factors that do not contain logarithms of N .In particular, Γ (5,res-sub) where E (sub) (N, µ 2 0 , µ 2 ) is defined in Eq. (2.24).Similarly, where By construction, Γ (1) (N, ξ) coincides with the 5-flavour result of [24].We note that Γ (3) (N, ξ) is fully consistent, to NLL accuracy, with Γ (4,res) ℓ 1 =0 given in Eq. (2.18), but differs from it in three respects.First, the argument of the running coupling in Γ (3) (N, ξ) was chosen to be µ 2 0 , which is of the order of the heavy quark mass, while it is chosen to be µ 2 ≃ q 2 in Eq. (2.18).The difference is subleading; the choice µ 2 0 appears to be more natural for a 4flavour calculation.Second, the momentum-space calculation of the jet functions is performed with the strong coupling at two loops, in the CMW scheme.Third, all logarithms of ξ are now exponentiated in the jet functions.
The contribution Γ (2) (N, ξ) performs the all-order matching between the two results.Note that, while our calculation determines all NLL contributions (both in ξ and N ), it does not determine the matching constant Γ (match) .Furthermore, because in this intermediate region 1 − x < √ ξ, we cannot simply match to the massless calculation, as done in the first region.However, we note that, in this region, illustrated in Fig. 4, the hard-collinear part of the jet function j is evaluated at the transition 1 − x = √ ξ.It seems therefore rather natural to choose Γ (match) = Γ (5,res,sub) The three functions Γ (1) , Γ (2) , Γ (3) are shown in Fig. 7 as functions of N , taken to be real and positive, for m = 5 GeV and two different values of q, namely q = 25 GeV, on the left, and q = 125 GeV, on the right.On each plot, we also indicate three vertical lines.Two of them mark the values of N that corresponds to the transitions in momentum space, i.e.Resummed formulae are not reliable beyond this value.Note that the vertical line corresponding to N = 1 ξ does not appear on the plots for q = 125 GeV, because in this case ξ −1 > NL , We also compare to the results obtained in the massless scheme, Eq. (2.26), with µ 2 = µ 2 F = q 2 , and µ 2 0 = µ 2 0F = m 2 (in blue), and in the massive scheme, Eq. (2.18), but with m 2 as the running coupling scale (in orange).
Let us comment the results in the different regions.We start by the functions to be used in the region 1 − x > √ ξ, which are shown in the top plots of Fig. 7.Because of the choice of the matching function, Eq. (4.5), we have that Γ (1) = Γ (5,res) ℓ=1,ℓ 2 =1 .Thus, in this region our result coincides with the one of Ref. [24].The massive result instead differs substantially and the difference because more pronounced at larger q 2 because of the presence of un-resummed logarithms.
The plots in the middle instead show our results for region 2 (in green), compared to the results present in literature.Here, in principle, our result differs from the massless calculation of Ref. [24] because of the treatment of the heavy-quark threshold in the running coupling.However, this difference turns out to be remarkably small, irrespectively of the value of ξ and N , as demonstrated by the fact that the blue and green curves largely overlap.Finally, in the third region (bottom plots) our result differs from both massless and massive calculations.While the former is not surprising because we have performed the calculation in the quasi-collinear limit, the latter is due to the exponentiation of the double logarithmic mass contributions discussed in Eq. (2.22).
It is also interesting to consider our Mellin space result Eq. (4.2) as a single function of N , on the real N axis.This way, we can interpret the different x-space regions Θ-functions in Nspace, identifying 1 − x → N −1 .Our results are shown in Fig. 7, again for q = 25 GeV on the left, and q = 125 GeV on the right.The colour-code for the different Γ (i) is the same as before, but now they are shown with solid lines within their N -space region of validity and transparent otherwise.We note that, while the exponential functions in Eq. (4.2) are continuous across these regions, the subtracted contributions introduce discontinuities.However, because of the choice of the matching condition, Eq. (4.5) the resulting function is actually continuous at N −1 = √ ξ and the only discontinuity appears at N −1 = ξ.This discontinuity has also been discussed in the literature, e.g.[39].
We conclude this discussion by noting that our calculation can be generalised to the case of the decay of a massive colour singlet into two quarks with different masses, which is relevant, for instance, for the decay of a W ± vector boson.If we consider the emission of a gluon, we have In order to work with an observable that vanishes on Born kinematics, we define a new scaling variable: Clearly χ → 1 when the gluon becomes either soft or collinear to p 2 .Thus, one should consider Mellin moments with respect to χ The calculation of the jet functions proceeds in the same way as described in the previous section, except that in this case the structure of integration regions is more complicated, due to the presence of two distinct thresholds.We have not performed the explicit calculations, but we have checked that logarithms of ξ 1 are resummed by DGLAP evolution, while double logarithms of ξ 2 are exponentiated in the recoil jet function, as expected.
Thus far, we have focussed our discussion on the decay of a colour singlet into a heavyquark pair.However, our calculation relies upon factorisation properties of QCD matrix elements in the quasi-collinear and soft limits, which are universal.Therefore, our formalism can be applied to different processes.In this section, we briefly discuss some applications, highlighting similarities and differences with other approaches.Detailed phenomenological studies at the level of physical processes are left to future work.
The first process we would like to discuss is heavy-quark decay where V is colour neutral boson, and X stands for unresolved QCD radiation.A very interesting case for collider physics is t → bW .In this context, one can combine a theoretical calculation that correctly treats the b-quark kinematics, with mass effects, together with a fragmentation function approach that is able to resum the large logarithms of m b /m t .In the limit where QCD radiation X becomes either soft or collinear to the fermions' directions, allorder resummation becomes relevant.Soft-gluon resummation for top decay has been studied in [35,37,57].In these studies, one is usually interested in the kinematic distributions of the b-quark (or B hadron).It follows that QCD evolution of the b-quark is always described by the measured jet function, and so there is no issue with non-commuting limits [12,35].Therefore, the matching of the massless and massive soft-gluon resummed calculations is straightforward.
The situation is rather different for kinematic distributions of the vector boson decay products: the W in the aforementioned t → bW decay, or the photon in b → sγ processes.Soft gluon resummation in these contexts has been widely studied, e.g [39,[58][59][60][61][62][63][64][65][66][67][68].In these examples, the final-state heavy quark is described by the recoil jet function and therefore the massless and soft limits no longer commute.This problem has been addressed by Aglietti et al., who were able to derive a resummed expression that is valid in both massive and massless regimes [38,39].The basic idea is the same as the one we have described in this paper, namely one has to consider the quasi-collinear limit with respect to the unmeasured quark.However, in contrast to our results, Aglietti et al. find a smooth function (in momentum space) that interpolates between the two regimes, rather than a function which is defined by cases.In our understanding, the main difference between the two approaches is that, while we compute the unmeasured jet function using the massive splitting function, Aglietti et al. exploit the dipole formalism with massive partons [41].Because the massive dipole reduces to the splitting function in the quasi-collinear limit, the two approaches agree at NLL level, but treat non-logarithmic corrections (related to recoil against soft and collinear emissions) differently.It would be interesting to numerically compare the two approaches.However, it is not clear to us how flavour thresholds are implemented in Refs.[38,39] and so we leave any detailed comparison to future work.
Another interesting process to be mentioned is deep-inelastic scattering (DIS) with heavy flavours.In particular, let us consider charged current (CC) DIS, either the scattering of a neutrino ν off a nucleus N , νN → lX, where l is a charged lepton, or lN → νX scattering.In both cases, we are interested in the situation in which charm flavour is identified in the final state.The parton-model contribution to these processes is given by where the initial-state quark is taken massless. 4It is interesting to note that CC-DIS experiments cover a rather large kinematical range in both the Bjorken scaling variable x B = Q 2 2P •q , where P is the momentum of the nucleon and Q 2 = −q 2 > 0 is the momentum transferred.Depending on the relative size of the charm quark mass m 2 c and Q 2 , a massless calculation or a massive one maybe more appropriate.Furthermore, if we want to investigate the large x B as done, for instance, in fixed-target experiments, then we face the issue of consistently combining the resummation of collinear logarithms and large-x logarithms.This problem was first investigated in Ref. [36], where two separate formulations of soft-gluon resummation were derived: one valid in massless limit m 2 c ≪ Q 2 and one that kept the full charm-mass dependence, appropriate for the region m 2 c ∼ Q 2 .The authors of Ref. [36] noted that the non-commutativity of the m 2 c /Q 2 → 0 and N → ∞ limits prevented them from obtaining a single resummed formula.Thanks to the formalism we have developed, we can bridge the gap between these two formulations, arriving at a unified resummed expression that interpolates between the two regimes.
Let us briefly describe the application of our formalism to this situation.In DIS experiments, a measurement of x B probes the dynamics of the incoming particle.Therefore, we can defined a measured jet function that describes radiation collinear to the initial-state quark.This is typically a massless jet function.However, note that, in contrast to the heavy-quark fragmentation case, where we had a perturbative initial condition, here we must supplement the description of radiation collinear to the initial state with non-perturbative PDFs, which also require specifying the flavour number scheme.The unmeasured jet function instead describes radiation collinear to the final-state charm.Our calculation of j allows us to interpolate between the massless m 2 c ≪ Q 2 and massive m 2 c ∼ Q 2 , thus consistently resumming both soft and mass logarithms.Phenomenological studies concerning the application of our formalism to CC DIS are work in progress.

Conclusions and outlook
In this paper, we have considered all-order calculations in processes with heavy flavours (mainly b quarks) in the final state.We have first reviewed the two different calculational frameworks that are usually employed for these processes, namely the massive and the massless schemes.The former allows us to fully take into account, at a given perturbative order, the dependence on the quark mass, while the latter allows us to resum, to all perturbative orders, large logarithms of the ratio of the quark mass over the hard scale of the process.Massive and massless calculations are usually combined together to obtain more reliable theoretical predictions.In particular, we have considered differential distributions, focussing on the energy fraction x of a heavy quark, produced in H → b b decays.The measurement of x introduces a further scale and logarithms of 1−x, which become large in the soft (or collinear) limit, should be resummed to all orders.We have reviewed the state-of-the-art results for soft-gluon resummation in both the massive and massless schemes.
Ideally, one would like to combine massive and massless calculations, both supplemented with their own soft-gluon resummation.However, standard matching procedures do not work, essentially because of the non-commutativity of the x → 1 and m 2 /q 2 → 0 limits.The main result of this paper is a modification of the massless calculation that makes this combination possible.This is achieved by performing the calculation of the jet function related to the unmeasured b quark in the quasi-collinear limit, rather then in the strictly collinear one, as it usually done in the massless calculation.A second achievement of this work is a detailed analysis of the role played by the heavy-quark flavour threshold on resummed expressions.
Our calculation is performed in momentum space and exploits the Lund representation in order to better identify the different kinematic regions.It is then transformed to Mellin space to next-to-leading logarithmic accuracy, so that it can be easily combined with results present in the literature.Thanks to our results, we can consistently resum soft and mass logarithms, while taking into account the effect of heavy-quark thresholds.While we leave detailed phenomenology to future work, we have performed some studies in Mellin space.We find that the effects we are studying become relevant only for rather low invariant masses, while are small for on-shell Higgs decay.Nevertheless, it would be interesting to combine our formalism with recent high-precision studies of heavy-quark fragmentation, see e.g.[26,27,69], and investigate our approach beyond NLL, where because of soft-gluon interference effects, it is no longer possible to factorise the resummation formula into the product of two independent jet functions.
We have also discussed other processes that present similar features, namely heavy-quark decay, commenting on similarities and differences of our approach and the one by Aglietti et al., and DIS with heavy flavours.We find the DIS case particularly interesting because experimental data span several orders of magnitude in Q 2 .We plan to focus on these processes for future phenomenological studies.
We conclude by noting that, as a by-product of our study, we have reached a detailed understanding of the use of the Lund plane to perform resummed calculations with massive quarks.In conjunction with recently developed IRC safe flavour-jet algorithms [70][71][72][73], this progress will allow us to extend existing resummed calculations for jet substructure to the case of heavy-flavour jets.These include, for instance, jet angularities [74][75][76], Soft Drop observables [77][78][79][80][81] and the primary Lund plane density [82], opening up exciting possibility to perform precision heavy-flavoured jet substructure studies at the LHC.
The resummed initial condition for the fragmentation function was computed in [24]: log N , n l = 4.The resummed coefficient function is given by: , (A.10) where both ∆ and J are computed in the strict collinear (massless) limit: ∆ N, µ 2 F q 2 , µ 2 R q 2 + J N, µ 2 R q 2 = c 1 log N + c 2 , (A.11) and λ = α s (µ 2 R )β (n f ) 0 log N , n f = 5.The only difference between photon, Higgs and Z boson is due to the overall constant C (1) (see [49] and [20] where we have changed one integration variable to ζ = k 2 t +z 2 2 m 2 q 2 z 2 2 .The strong coupling may now be expressed as a power series in log (z 2 2 (ζ − ξ)).We focus on the most singular term in the splitting function, and so we consider the integral: Let us consider first the situation 1 − x > ξ.In this case, up to power corrections in ξ, we can neglect the ξ dependence in the integrand and keep only the l = 0 contribution in the sum.This way, we recover the standard (massless) result.
The second case, 1 − x < ξ, requires a more careful treatment.We notice that

Figure 2 :
Figure 2: Lund plane representation of a soft-collinear emission off a hard quark-antiquark dipole, in the dipole rest frame for massless (left) and massive (right) quarks.The observable V under study takes a constant value v along solid dashed lines.The shaded area represents the region where the cumulative distribution differs from zero.The vertical lines η i = − 1 2 log ξ indicate the position of the the dead cone.The line k 2 t = ξq 2 = m 2 marks the separation between regions of phase-space with 4 or 5 active flavours.

. 52 )
Lowering the value of 1 − x, we enter the intermediate region ξ < 1 − x < √ ξ, where the x dependence arises both from the n f = 5 and from the n f = 4 contributions, as seen by inspection of the Lund diagram in Fig 4. Therefore, in the intermediate region we find j

N − 1
∼ 1 − x = √ ξ and N −1 ∼ 1 − x = ξ.The third line indicates the position of the Landau pole of the strong coupling,

Figure 6 := µ 2 0F = m 2 , and µ 2 = µ 2 F
Figure6: Our results, Eq. (4.2), plotted as a function of N , on the real axis, for two different values of the hard scale.Plots on the left are for q = 25 GeV, while the ones on the right for q = 125 GeV.In both cases we set have µ 2 0 = µ 2 0F = m 2 , and µ 2 = µ 2 F = q 2 .The vertical lines mark the values of N that corresponds to the transitions in momentum space, i.e.N −1 ∼ 1 − x = √ ξ, and N −1 ∼ 1 − x = ξ, and the position of the Landau pole.Note that for the 125 GeV case, ξ −1 > NL .

Figure 7 :
Figure 7: Our results, Eq. (4.2), plotted as a single function of N , on the real axis.