The energy distribution of subjets and the jet shape

We present a framework that describes the energy distribution of subjets of radius r within a jet of radius R. We consider both an inclusive sample of subjets as well as subjets centered around a predetermined axis, from which the jet shape can be obtained. For r ≪ R we factorize the physics at angular scales r and R to resum the logarithms of r/R. For central subjets, we consider both the standard jet axis and the winner-take-all axis, which involve double and single logarithms of r/R, respectively. All relevant one-loop matching coefficients are given, and an inconsistency in some previous results for cone jets is resolved. Our results for the standard jet shape differ from previous calculations at next-to-leading logarithmic order, because we account for the recoil of the standard jet axis due to soft radiation. Numerical results are presented for an inclusive subjet sample for pp → jet + X at next-to-leading order plus leading logarithmic order.


Introduction
In this paper, we study the energy distribution of subjets with radius r inside a jet of radius R, as illustrated in fig. 1.We will consider jets defined through the anti-k T algorithm [1] or an infrared and collinear-safe cone algorithm.Subjets are obtained by reclustering the particles inside the reconstructed jet with radius parameter r < R. In addition, we consider the subjet of radius r centered around the standard jet axis or the winner-take-all R r Figure 1.Illustration of a subjet with radius r (red) inside a jet with a radius R (green).We focus on describing the energy fraction z r of the jet that is carried by the subjet.axis (WTA) [2].We mostly focus on an inclusive jet sample pp → jet + X, but briefly discuss how our framework can be extended when a veto on additional jets is imposed.
Specifically, we will develop the theoretical framework to calculate where z r is fraction of the jet energy contained in the subjet of radius r, and η and p T are the rapidity and transverse momentum of the jet with radius R. First proposed in ref. [3], the subjet distribution simultaneously provides information about the longitudinal and transverse energy distribution inside jets, through z r and r.The subjet observables considered in this work are connected to both the standard jet shape [4][5][6][7] and the jet fragmentation function [8][9][10][11].On the one hand, the jet shape is the average value of z r as function of r, for subjets centered on the jet axis.On the other hand, the jet fragmentation function describes the longitudinal momentum (or energy) fraction of hadrons in the jet.This is the r → 0 limit of the inclusive subjet energy fraction z r , where the collinear singularity is now cut off by hadronization instead of r.The jet shape [12][13][14][15] and the jet fragmentation function [16][17][18] have been measured by the LHC experimental collaborations in both proton-proton and heavy-ion collisions.We expect that similar measurements are feasible for the observables discussed in this work.There are several ways in which subjet distributions are valuable to present day collider phenomenology: Firstly, the various distributions of subjets discussed in this work provide a powerful test of our understanding of perturbative QCD at very high energies.Studying the energy distribution of subjets probes both the longitudinal and the transverse momentum distribution within jets at a more differential level, and may extend our current understanding of the underlying QCD dynamics.Secondly, the distribution of subjets can be used for discriminating QCD jets from boosted heavy objects, such as W bosons or top quarks.Their hadronic decays would produce two or three jets, which become the subjets of one fat jet due to their boost.This plays an important role in many searches for Beyond the Standard Model (BSM) physics [19][20][21][22].Many of the taggers used for identifying such a two or three prong decay are quite sensitive to soft radiation [23][24][25][26][27][28].For several of the subjet observables we consider, this soft sensitivity is power suppressed and collinear factorization is sufficient.In addition to being theoretically more robust, a reduced sensitivity to soft radiation is also advantageous experimentally due to the messy LHC environment.Our work on the inclusive subjet distribution and the distribution of subjets centered about a specified axis provides a first step on the direction of taggers that are less sensitive to soft radiation.An example of a more direct connection to BSM searches is given in ref. [29], where the authors proposed to use jet shapes ("jet energy profiles") to search for new physics at the LHC.Another application in the context of jet substructure is the discrimination of quark and gluon jets using e.g.fractal observables defined on subjets rather than hadrons [30].Thirdly, the subjet distribution is particularly suited for measuring the modification of jets in heavy-ion collisions, see e.g.[31][32][33].Jets that traverse the quark-gluon plasma get modified both in the longitudinal and transverse momentum direction, as can be collectively seen from the modification of longitudinal jet fragmentation function [17] and transverse jet energy profile [14].By identifying subjets inside a reconstructed jet, the correlations between these effects can be studied in a single measurement.An advantage of using subjets over hadrons is that the subjet distribution does not require the additional non-perturbative input of fragmentation functions.
Our setup relies on collinear factorization: first we exploit that the radius R is small, to factorize the dynamics of the jet from the rest of the cross section. 1 For pp collisions, we have Here f a,b denote the parton distribution functions and H c ab are hard functions describing the production of an energetic parton of flavor c with transverse momentum p T /z and rapidity η with respect to the beam axis.The subjet functions G jet c describe the subsequent conversion of that parton into a jet moving in (roughly) the same direction but with transverse momentum z × p T /z = p T , containing a subjet of radius r with fraction z r of the jet energy.The argument ω R = 2p T / cosh η of G jet c is the large light-cone component of the jet momentum, and the arguments r and R are suppressed.The symbols ⊗ denote convolution products associated with the variables x a,b and z, which are explicitly written out in eq.(3.41).Power corrections to the factorized cross sections are order R 2 suppressed.
We will consider both r R and r R. In the first case, only single logarithms of the form α n s ln n R need to be resummed to all orders.The subjet functions G jet c follow timelike DGLAP evolution equations allowing for the resummation of logarithms in the jet size parameter R [3,36,37] (see refs.[38,39] for a generating functional approach to jet radius resummation).For all subjet observables considered in this work, the resummation of the logarithms of R is the same.For r R, we encounter additional large logarithms of r/R.The structure and resummation for this class of logarithms depends on how the subjet of size r is identified.For an inclusive subjet sample, we perform an additional collinear factorization for the subjet, matching the subjet functions G jet c onto a semi-inclusive jet function [3,36] for the subjet.This enables us to resum single logarithms α n s ln n (r/R) using another DGLAP type evolution equation.
The refactorization for central subjets (i.e.those centered around a specific axis) in the limit r R differs from the inclusive subjet case, and crucially depends on the choice of axis.The standard jet axis is sensitive to soft radiation inside the jet, since the jet axis is aligned with the total jet momentum.By contrast, the winner-take-all axis is insensitive to soft radiation, but the location of the axis depends on the details of the collinear radiation.For the winner-take-all axis our factorization enables resummation to all-orders in perturbation theory using a (modified) DGLAP evolution [40], whereas for the standard jet axis this is complicated due to non-global logarithms.
We will also calculate the average z r value from eq. (1. 2) for the central subjet, which corresponds to the jet shape.Its cross section has a single logarithmic dependence on r/R for the winner-take-all axis and a double logarithmic dependence for the standard jet axis.Our factorization formula for the standard jet shape for r R differs from earlier approaches [5][6][7] 2 .Specifically, it involves a further refactorization to account for the soft radiation that recoils against the jet axis.The additional logarithms of r/R that we can resum, enter the cross section at next-to-leading logarithmic order.This is similar to the broadening event shape where the recoil of soft emissions in the direction of collinear particles (which only enters at NLL order) was initially overlooked [41] and only realized later [42].Like in ref. [43], the effect of recoil can be removed by using the winner-take-all axis.However, in this case this removes all soft sensitivity.
The outline of our paper is as follows: In sec. 2 we revisit the calculation of the semiinclusive jet function, addressing an inconsistency in the literature for cone algorithms, and presenting corrected analytical results.We discuss the inclusive production of subjets in sec.3 in terms of the subjet function, for both cone and anti-k T algorithms, and show numerical results for eq.(1.1) for pp → (jet j r ) + X at NLO+LL R +LL r/R .In secs. 4 and 5 we focus on subjets centered on the winner-take-all axis and the standard jet axis, respectively.In all sections, r R as well as r R are considered, and all matching coefficients are calculated at NLO.The jet shape is the second moment (average z r ) of the result in secs.4 and 5, and can be directly related to TMD fragmentation, as discussed in sec.6.We conclude in sec.7 and provide an outlook.

Inclusive cone jets revisited
In this section we review the calculation of the semi-inclusive jet functions (siJFs), which enter in the cross section for single inclusive jet production, pp → jet + X.Specifically, the cross section for inclusive jet production satisfies the factorization theorem in eq.(1.2), after replacing G jet c by the siJF J c [36].We first address an inconsistency in the literature for cone algorithms, before considering the calculation of subjet functions in the following sections (as we also present results for cone algorithms there).Our default notation will be for e + e − algorithms, where a jet is defined in terms of its energy E = ω R /2 and angle R.These results equally apply to pp algorithms, with the replacement ω R R → 2p T R in terms of a jet radius defined in (η, φ) coordinates, see e.g.ref. [35].For single inclusive jet production in proton-proton collisions at NLO, pp → jet + X, there are either one or two final-state partons inside the observed jet, whose possible assignments are illustrated in fig. 2. In refs.[11,34,36,[44][45][46] the cone algorithm was (effectively) implemented for two final-state partons as Partons in single jet: Partons in separate jets: where β 1 and β 2 are the angles of the final state partons with respect to the initiating parton.However, these regions of phase space are not complementary, and there are configurations with R < β < 2R that are double counted.The resolution depends on the specifics of the cone algorithm.For example, if only the particles themselves are used as seeds for the cone algorithm, the first criterion requires modification and the resulting algorithm happens to coincide with anti-k T (for two parton configurations).For the midpoint [47] and the SISCone [48] algorithms, the correct implementation is Partons in single jet: Partons in separate jets: Of course the midpoint and the SISCone algorithms will differ with additional particles.We now consider the calculation of the semi-inclusive quark jet function.The requirement that the partons are in separate jets in eq.(2.2), leads to the following expression for the case where the quark is inside the observed cone jet in fig.2(B), 3) The corresponding expression when the gluon makes the observed jet, fig. 2 The collinear phase space and (squared) matrix element in eq. ( 2.3) are given by where x is the momentum fraction and q ⊥ is the transverse momentum of (one of) the final partons with respect to the initiating quark.The angles β 1 and β 2 can be expressed in terms of x and q ⊥ as follows for the partons with momentum fraction x and (1 − x).The angle between the partons is After evaluating the integrals in eq. ( 2. 3) and combining it with the result when both partons are in the jet [11,34,36,49], we obtain the new results for the cone semi-inclusive jet function: where L R is defined as ( The result for gluon-initiated jets can be obtained in a similar way, In other words, [36]   q (z, ω R , µ) [36]   g (z, ω R , µ) (2.10) As the results for the semi-inclusive jet functions in ref. [36] are consistent with earlier analytical results for single inclusive jet production for cone algorithms in refs.[11,34,44,45], the above equations also provide a correction to these earlier cone jet results.As a consistency check, we note that only the updated cone jet results in eq.(2.10) satisfy the following momentum sum rule introduced in where the large momentum component of the initiating parton ω = ω R /z is held fixed.Similarly, the NLO matching coefficients for the jet fragmentation function as presented in [11,46] for cone jets are modified as follows For more details, we refer the interested reader to the earlier publications listed above.

Inclusive subjets
In this section, we study the (semi-inclusive) subjet function (SJF), which describes the energy distribution of all subjets inside a jet as in eq.(1.2).We gives its definition in sec.3.1, and calculate it to next-to-leading order (NLO) in sec.3.2.In sec.3.3 we derive the renormalization group equation (RGE) of the subjet function, which we use to resum the logarithms of the jet radius R. We subsequently consider r R in sec.3.4, performing the matching onto semi-inclusive jet functions (siJF) that describe the subjets of radius r, and use this to resum the large logarithms of r/R.The limit r → R is discussed in sec.3.5 and the fragmentation limit r → 0 is considered in sec.3.6.In sec.3.7 we discuss the subjet function for exclusive jet production.We will drop the adjective "semi-inclusive" in front of the SJF, since we restrict ourselves to inclusive jet samples everywhere else.In sec.3.8 we show numerical results for the momentum fraction of subjets in pp → (jet j r ) + X at NLO+LL R +LL r/R .

Definition of subjet function
We define the subjet function as a matrix element in Soft-Collinear Effective Theory (SCET) [50][51][52][53].In our definitions and calculations e + e − jet algorithms will be our default, which define a jet in terms of its energy E = ω R /2 and angle R, and similarly for the subjet.Our results directly apply to pp algorithms with the replacement ER → p T R, where the jet radius parameter R now refers to a distance in (η, φ) coordinates.The subjet functions for quark and gluon-initiated jets G jet q and G jet g are defined as suppressing the dependence on r and R in the arguments.Here n µ = (1, n) is a light-cone vector with its spatial component n along the jet axis, while nµ = (1, −n) is a conjugate light-cone vector such that n 2 = n2 = 0 and n • n = 2. χ n and B µ n⊥ are gauge invariant collinear quark and gluon fields in SCET.The sum over states |X runs over all final-state particles and includes their phase-space integrals.We sum over all reconstructed jets with radius parameter R in X and all subjets j r with radius r.The large light-cone momentum (approximately twice the energy) of the initiating parton, jet and subjet are denoted by ω, ω R and ω r , respectively.For the field producing the initiating parton this is encoded using the (label) momentum operator P. The variables z and z r describe the momentum fraction of the initiating parton carried by the jet, and that of the jet carried by the subjet, and are thus given by

NLO calculation
We now calculate the subjet function for quark-initiated jets, G jet q (z, z r , ω R , µ).For definiteness, we discuss the case where the anti-k T algorithm is used for reconstructing both the larger jet of size R and the subjet of size r, i.e. we consider "anti-k T -in-anti-k T ".However, at the end of this section we also present results for "cone-in-cone" and "cone-in-anti-k T ", see eq. (3.16).For "anti-k T -in-anti-k T " and "cone-in-cone" our calculations reveal that, at least at next-to-leading order, these results can be fully expressed in terms of known quantities by eq.(3.24).The calculation of the gluon subjet function G jet g follows the same steps.Note that the jet algorithms anti-k T , k T [54,55] and Cambridge/Aachen [56,57] yield the same results to the order that we are considering.
At leading order, the subjet functions are simply given by since the total energy of the initiating parton is transferred to the jet (of size R) and subjet (of size r).We perform the next-to-leading order calculation in pure dimensional regularization, where all virtual diagrams vanish.The collinear matrix element and phasespace were given in eq.(2.4).The five possible assignments of the partons over the (sub)jet are shown in fig.3, which we discuss in turn: Figure 3.The five configurations that enter for the quark subjet function at O(α s ): (A) the quark and gluon are inside the jet and subjet, (B) only the quark is inside the subjet but both partons are in the jet, (C) only the gluon is inside the subjet but both partons are in the jet, (D) only the quark is inside the jet and subjet, (E) only the gluon is inside the jet and subjet.
(A) The quark and gluon are inside the jet and subjet All the initial quark energy is transferred to the jet and subjet, so z = ω R /ω = 1 and where the θ-function encodes the constraint that both partons are inside the jet and subjet when using the anti-k T algorithm.Performing the q ⊥ and x integrals and expanding in powers of , we find where L r is defined as (3.6) We choose to write the results for all configurations (A) -(E) in terms of ω R .
(B) Only the quark is inside the subjet but both partons are inside the jet The energy of the quark initiating the jet is transferred entirely to the jet z = 1, but only the fraction z r < 1 is contained inside the subjet.This contribution is The θ-function encodes the constraint that the partons are sufficiently close together to be clustered into the jet but not so near that they are in the same subjet.Here L R is defined in eq.(2.8) and L r/R is given by (C) Only the gluon is inside the subjet but both partons are inside the jet This configuration is analogous to (B) but exchanging the quark and gluon as shown in fig.3(C).This amounts to replacing δ where the quark splitting functions we use are defined as (D) Only the quark is inside the jet and subjet For the configuration shown in fig.3(D), only a fraction z < 1 of the initiating quark's energy is transferred to the jet of size R.However, all the energy of the jet is inside the subjet, z r = 1.Thus its contribution to the SJF is given by The only constraint from the jet algorithm is that both partons are far enough apart that they are not clustered together into the jet.
(E) Only the gluon is inside the jet and subjet This is analogous to (D) but with the replacement δ(x−z) → δ(1−x−z) in eq.(3.11), Summing up the leading-order result as well as the five contributions at O(α s ), we obtain All 1/ 2 poles cancel in the sum as well as all double logarithms L 2 R and L 2 r/R .The result at NLO always involves δ(1 − z) and/or δ(1 − z r ), since the final state consists at most of two partons, but this structure does not generalize to higher orders.The remaining 1/ pole is a UV divergence that will be removed by renormalization and the resulting time-like DGLAP equation can be used to resum logarithms of R, as discussed in the next section.The result for the gluon SJF is: which involves the splitting functions These results agree with the general setup of ref. [3].However, there the contributions from (A), (B), (C) are treated separately (factorized) from (D) and (E) (see e.g.their eq.( 44)).This allows them to relate their expression to exclusive results.As these contributions are at the same scale, the factorization probably fails at higher order in α s .
For the semi-inclusive fragmenting jet function [46], a similar structure was obtained as in eqs.(3.13) and (3.14).However, this case involved an additional IR pole that cancelled in the matching onto the standard fragmentation functions.For the SJFs, this IR pole is regulated by the size of the subjet r leaving a single logarithmic dependence on the ratio r/R, i.e.L r/R .In sec.3.4, we are going to match the SJF onto a semi-inclusive jet function for the subjet.This will lead to another time-like DGLAP equation that can be used to resum logarithms of r/R.
We conclude this section by giving the renormalized one-loop results for the cone-incone and cone-in-anti-k T SJFs (3.16)

Renormalization and resummation of ln R
The renormalization and resulting RG equation of the subjet function is identical to that of the semi-inclusive jet function, because the additional measurement of the subjet does not modify the UV behavior.This also follows from consistency, since the SJFs and siJFs can be interchanged in factorization theorems.For completeness we still present the essential equations.The renormalization of the SJF is given by which leads to the following RG evolution equation with the anomalous dimension matrix γ ij .From our NLO calculation, we immediately obtain so the SJF satisfies the usual time-like DGLAP evolution equations.The one-loop renormalized quark SJF is given by This result for G can be used when the jet radii of inner and outer jets are comparable, r R. By evaluating it at the scale and evolving it with eq.(3.18) to the scale of the hard scattering the logarithms of µ R /µ H ∼ R are resummed.When r R, the logarithms L r/R in G jet i become large and require resummation as well.This is achieved by an additional factorization, which we discuss next.

Matching for r
R and resummation of ln(r/R) In the regime r R we have the following matching equation to all orders in α s (3.23) where J j are the semi-inclusive jet functions describing the subjet of size r.In this equation we have explicitly shown the dependence on the jet radius R and the subjet radius r in the arguments of the functions, to highlight its structure.This matching equation is very similar to that for the matching of the semi-inclusive fragmenting jet function onto fragmentation functions.In fact, the matching coefficients J ij are the same, since they are independent of r and we can therefore safely take the fragmentation limit r → 0. We have also verified this through a direct calculation, and therefore do not give these matching coefficients, but refer the reader to ref. [46] for anti-k T and sec. 2 for cone algorithms.
Interestingly, our calculations reveal that there are no O(r 2 /R 2 ) power corrections in eq.(3.23) at NLO.In fact, for anti-k T -in-anti-k T and cone-in-cone the NLO subjet function is fully determined by For cone-in-anti-k T , this equation does receive corrections when r > R/2.
Having performed the factorization in eq.(3.23), we can now resum the additional logarithms of r/R with the help of another DGLAP RG equation.The scale µ R in eq.(3.21) sets the large logarithms L R to zero in the matching coefficients J ij , whereas is the natural scale for semi-inclusive jet function J j describing the subjet in eq.(3.23).By using the RG equation to evolve the siJF from µ r to µ R , we resum the single logarithms of r/R.

The limit r → R
We will now start from the regime r R, for which the logarithms of r/R are resummed, and consider the limit when the inner jet radius becomes large.It is fairly straightforward to do this, because the absence of power corrections in eq. ( 3.23) at this order.The aim is to gain an analytical understanding of the r → R limit, for which we consider anti-k T -inanti-k T at leading-logarithmic accuracy.
We start by observing that for z r < 1, the only terms in eq.(3.13) that contribute are since all other terms are proportional to δ(1 − z r ).Note that the role of the variables z and z r are fundamentally different, in the sense that z is an integration variable for the convolution with a hard function and z r is the measured external variable.As can be seen from eq. (3.26), the subjet cross section at NLO is directly proportional ln(r/R) and, therefore, it can be considered as a direct probe of the ln(r/R) resummation.
In the limit r R, we refactorize the subjet function G jet i in terms of matching coefficients and evolved siJFs, see eq. (3.23).In order to recover the NLO result in eq.(3.26) from the resummed result in the limit r → R, we find that it is sufficient to consider the leading-order matching coefficients as well as the leading-order initial condition for the siJFs for the subjets Including O(α s ) corrections in the matching or initial condition would generate O(α 2 s ) terms in G jet i .Using the techniques of refs.[36,46], we solve the DGLAP equations associated with the resummation of both logarithms ln(r/R) and ln R in Mellin moment space.In order to perform the resummation, we take double Mellin moments of the subjet functions We only discuss the resummation of logarithms ln(r/R) associated with the variable z r and the Mellin variable N .(The DGLAP equation for resumming ln R was given in eq.(3.18) and is associated with the variable z and Mellin variable M .)The convolution of matching coefficients and the siJFs in eq.(3.23) turns into simple products The delta functions of the leading-order matching coefficients in eq.(3.27) integrate to 1 when taking moments, so the subjet function in Mellin space is given by the siJF evolved from µ r to µ R .For the quark siJF at LL accuracy, we find where P ji (N ) are the leading-order Altarelli-Parisi splitting functions in Mellin N moment space and β 0 is the first coefficient of the QCD beta function.Inserting the leading-order solution for the running strong coupling constant this becomes (3.33) In the limit r → R, Finally, we need to perform the double Mellin inverse transformation of the whole subjet function which is given by contour integrals in the complex M and N planes The first term in eq.(3.34) does not contribute to the cross section, since it does not contain poles in N .The second term in eq.(3.34) directly gives the NLO contribution for z r < 1 shown in eq.(3.26) above.Therefore, we have recovered the fixed-order cross section from the resummed result in the limit of r → R. We also verified this numerically in sec.3.8.Note that the inverse with respect to N in eq.(3.35) can be taken directly as it is associated with the observed external variable z r .However the resulting expression for z still needs to be convolved with the hard function.

3.6
The fragmentation limit r → 0 For sufficiently small r, the scale ω r r/2 of the semi-inclusive jet function with radius r becomes nonperturbative.At that point we can no longer speak of subjets but are really probing individual hadrons, and the subjet function should be replaced by a fragmenting jet function (inclusive in hadron species).This limit is continuous, since the matching coefficients are the same whether we match onto subjets or hadrons, as discussed in sec.3.4.We stress that our conclusions also apply to inclusive jet cross sections, as we will study the nonperturbative corrections to the semi-inclusive jet function.
We start by factorizing the siJF into a perturbative and nonperturbative component, In contrast to the rest of the paper, we work in terms of the large momentum component ω of the initiating parton i, instead of ω r = zω of the (sub)jet.J pert i is the perturbative result (which was calculated at NLO in refs.[3,36]) and J NP i captures the nonperturbative corrections.From boost invariance (or reparametrization invariance [58]) we infer that the arguments ω and r appear in the combination ωr, which was exploited in writing the arguments of J NP i .Since J pert i has the same anomalous dimension as the full J i , this implies that J NP i is independent of µ.Note that this crucially relies on including the nonperturbative corrections through a Mellin convolution in eq.(3.36), since the anomalous dimension is only diagonal in Mellin space.
Although J NP i is a two-dimensional nonperturbative function, we know its limits: ωr → ∞ : The first line is the perturbative limit, for which the nonperturative corrections (but not J NP i ) vanish.From the continuity of the r → 0 fragmentation limit of subjets, it follows that the semi-inclusive jet function turns into the fragmentation function (summed over hadron species), as shown on the second line.
We have extracted J NP q from Pythia [59], using parton-level and hadron-level inclusive jet spectra for e + e − collisions at a center-of-mass energy of 500 GeV, with the e + e − antik T algorithm.This implies that ω = 500 GeV (whereas ω r varies).Instead of considering the full z r dependence we take Mellin moments, such that J NP q is the ratio of hadron-level and parton-level cross sections.The result is shown in fig. 4 as function of 1/r.In the perturbative limit ωr → ∞, J NP q (N, ωr) → 1, so to visualize the nonperturbative effects we plot 1 − J NP q .Also, J NP q (N = 2, ωr) = 1, due to the momentum sum rule in eq.(2.11) (which relies on using ω rather than ω r ).In the left panel we limit ourselves to ωr/2 > 10 GeV, for which it is reasonable to carry out a series expansion in 2Λ QCD /(ωr).We find that the linear term vanishes and the quadratic term (fit shown as solid line in left panel) describes the points very well.There is a slight discrepancy for large values of r, but in this regime the O(r 2 ) corrections to the factorization theorem may no longer be negligible.In to the semi-inclusive quark jet function for Mellin moments N = 2 (black), 3 (blue), 4 (green), 5 (red) as function of 1/r.The points were extracted from parton-level and hadron-level Pythia data, using e + e − collisions at a center-of-mass energy of ω = 500 GeV, with the e + e − anti-k T algorithm.In the left panel we restrict to large values of r and show a fit to c/r 2 (solid lines).In the right panel, the asymptotic approach to the fragmentation limit is shown.The nonperturbative corrections for N = 2 vanish due to eq. (2.11).
the right panel of fig. 4 we focus on the nonperturbative regime, displaying the asymptotic behavior.For ωr/2 ∼ 5 GeV the nonperturbative corrections deviate from the quadratic fit by about 10%, whereas for ωr/2 1 GeV, the jets essentially consist of single hadrons.

Subjets in exclusive jets
Up to this point we have focussed entirely on inclusive jet production.However, one can also consider exclusive jet production, where additional jets are vetoed.The collinear (energetic) radiation is then forced to be inside the jet.Adding this additional restriction to the definition of the subjet function in sec.3.1, In this case there is no collinear radiation outside the jet, which is why we replaced X by J R .Consequently, there is no dependence on z, since z = ω R /ω = 1 always.
In the NLO calculation, the contributions in fig.3(D) and (E) are absent.This does not modify the dependence on z r and r, but removes the z dependence and introduces double logarithms of R.These logarithms of R can again be resummed using the RGE of the subjet function.However, instead of the convolution structure seen in eq.(3.18), the RGE of the subjet function is now multiplicative, The anomalous dimension γ i,excl is the same as for the unmeasured jet functions of ref. [49], and given in eq. ( 6.26) therein.The appearance of double logarithms of R indicate a sensitivity to soft radiation and so the factorization in eq. ( 1.2) must be modified to include a soft function.
The matching for r R in sec.3.4 onto semi-inclusive jet functions still holds, (3.40) but the matching coefficients J ij,excl are not the same as in the inclusive case.Rather, they are the same as those of the fragmenting jet functions for exclusive jet samples, which were calculated in ref. [60] for cone algorithms and in refs.[61,62] for anti-k T .

Phenomenology for pp → (jet j r ) + X
We present numerical results for the momentum fraction of subjets measured on an inclusive jet sample pp → jet + X.In analogy with the hadron-in-jet calculations in proton-proton collisions presented in refs.[11,46], we adopt the notation pp → (jet j r ) + X, where j r denotes a subjet of size r inside the larger jet of size R.The factorization formula for the subjet distribution in proton-proton collisions is given by dσ pp→(jet jr)X dp T dη dz r = a,b,c × where the sum on a, b, c runs over all relevant partonic channels.The PDFs are denoted by f a,b and the hard functions are given by H c ab , which have been calculated to NLO in refs.[44,63].For all numerical results presented in this section, we use the CT14 NLO set of PDFs [64].The variables s, p T and η correspond to the center-of-mass (CM) energy, the jet transverse momentum and the jet rapidity respectively.The hard functions depend on the corresponding partonic variables ŝ = x a x b s, pT = p T /z and η = η − ln(x a /x b )/2.The lower integration bounds x min a , x min b and z min can be written in terms of these variables and are listed for example in refs.[11,46].
The subjet function G jet c in eq.(3.41) is evolved to the hard scale µ ∼ p T by solving the DGLAP evolution equations associated with the logarithms ln R and ln(r/R).Numerically, we solve the DGLAP equations in Mellin moment space using the techniques developed in refs.[36,46], which in turn are based on the evolution packages of refs.[65,66].We jointly resum both single logarithms ln R and ln(r/R) with a combined accuracy of "NLO+LL R +LL r/R ".The evolved subjet function G jet c is divergent for z → 1.We can nevertheless perform the integrals in eq.(3.41) by adopting the prescription of ref. [67], as discussed in detail in ref. [36].Note that the factorized form of the cross section in eq.(3.41) is a purely collinear factorization, i.e. there is no soft function.All numerical results presented here are normalized by the total inclusive jet cross section, see eq. (1.1), for which we resum single logarithms of the jet size parameter ln R at NLO+LL R accuracy.
For our numerical results we choose representative LHC kinematics, taking a CM energy of √ s = 13 TeV and a rapidity range of |η| < 1.2.Both the jet and the subjets are identified using the anti-k T algorithm.We choose a jet radius of R = 0.6 for the outside jet.In fig. 5 we plot the momentum fraction z r for a subjet radius parameter of r = 0.2 for different bins of the transverse momentum p T of the jet.We multiply the results for the different p T bins by multiples of 10 for better visibility.One immediately notices that the plotted curves look like the QCD Altarelli-Parisi splitting functions.This behavior of the cross section can be most easily understood by looking at the fixed order results for the subjet function in eq.(3.20).Only the terms ln(r/R)P ji (z r ) have a non-trivial functional dependence on z r , as all other terms are proportional to δ(1 − z r ) and do not contribute at fixed order.The ln(r/R) resummation modifies the distribution slightly, so it is not exactly the splitting function.
We would like to point out an important difference of the results for the subjet distribution compared to the distribution of light charged hadrons inside jets, as presented in for example refs.[11,46].When measuring an identified hadron inside jets, the distribution falls continuously as z h increases.However, as can be seen from fig. 5, the distribution of subjets starts to rise again for sufficiently large z r .Whereas it becomes increasingly unlikely to find a hadron that carries a large fraction of the complete jet, a subjet with radius r < R may still contain most of the energy of the larger outside jet as long as r is not too small.In order to better see the dependence on the jet p T , we plot on the left panel of fig.6, the same curves as in fig. 5 but without multiplying them by multiples of 10.We observe only a relatively small dependence on the jet transverse momentum.
Next, we study the dependence of the subjet distribution on the subjet radius parameter r.In the right panel of fig.6, we show the momentum fraction of the subjet for four different values of r, ranging from 0.05 to 0.3.We choose the same kinematics as in fig. 5 and restrict the jet transverse momentum p T to the bin [50,100] GeV.We find a relatively strong dependence on r which is also due to the ln(r/R) resummation effects.
To make this point more clear, we show in fig.7 the dependence of the cross section for fixed values of z r as a function of r/R.Results are shown for four values of z r , and the two panels corresponds to different bins for the jet transverse momenta.One notices again the strong dependence on r which can span two orders of magnitude.For small z r the curves increase continuously as r decreases, since one finds more and more subjets.However, for sufficiently large z r the curves flatten out as r becomes small and can even turn over.This behavior is more pronounced for the smaller jet p T interval of [25,50] GeV, and arises because it is not possible to capture a very large energy fraction z r of the jet within only a narrow subjet.

Central subjets for the winner-take-all axis
In this section we focus on the energy distribution of the subjet centered about the winnertake-all (WTA) axis [2].In sec.4.1 we treat the case r R, which parallels the discussion in sec.3. We discuss the factorization for r R and the resummation of logarithms of r/R in sec.4.2.

Central subjet function for r R
The difference between the standard jet axis and WTA axis resides in the merging step of a clustering algorithm (and is thus not defined for cone jets).Specifically, it chooses the axis to be along the most energetic of the two particles (or pseudojets) that are being merged.For the configuration of at most two partons in the jet, the winner-take-all axis is along the most energetic one.This can simply be accounted for by an additional factor θ(z r − 1/2) compared to the O(α s ) calculation in sec.3.2.For example, for quark jets where the tilde for the SJF indicates that we are restricting to the central subjet about the winner-take-all axis.At higher orders there will be more partons inside the jets, leading to more significant differences between the calculation for the central subjet and an inclusive sample of subjets.
The renormalization of the central subjet function is the same as that of the semiinclusive jet function.This is immediate at O(α s ) from the above, but holds at higher orders because the rest of the factorization theorem does not depend on whether the energy fractions of the central subjet is measured or not.The central subjet function therefore satisfies the DGLAP evolution equation By evaluating Gi at its natural scale µ R ∼ ω R R/2 and evolving it to the scale µ H ∼ ω R /2 of the hard scattering, the logarithms of µ R /µ H ∼ R are resummed.

Matching for r R and resummation of ln(r/R)
In the regime r R, the central subjet function will contain large logarithms of r/R that require resummation.In direct analogy to eq. (3.23), this resummation is accomplished by the following matching equation to all orders in perturbation theory (4.3)We first describe this factorization formula and then explain why it holds for the winnertake-all axis.
The object Jj onto which we match is not the semi-inclusive jet function J j , which describes the distribution of energy fractions for all subjets produced by a parton.Rather, it only picks out the energy fraction of the subjet centered on the winner-take-all axis.It also differs from the central subjet function Gjet i because all partons are clustered together, since we have effectively taken R → ∞.Jj is defined as At order α s there are at most two partons and the winner-take-all axis is along the most energetic one, so once again This implies that the one-loop matching coefficients in eq. ( 4.3) are given by and thus directly related to those for the inclusive case in eq.(3.23).In eq.(4.3) the Gjet i on the left-hand side contains physics at angular scales R and r, that are factorized into the objects Jij and Jj on the right-hand side.The validity of this equation at next-to-leading order follows immediately from the above.However, to use it for resummation requires the factorization to hold to all orders in α s .In particular, the axis finding must factorize between the scales r and R, i.e. the axis cannot be sensitive to radiation at the jet boundary.This was shown in ref. [40] for the winner-take-all axis when using Cambridge/Aachen or anti-k T , in the context of transverse-momentum-dependent fragmentation.
Having performed the factorization in eq. ( 4.3), we can now resum the additional logarithms of r/R with the help of another RG equation.The scale µ R ∼ ω R R/2 is the natural scale for the matching coefficients J ij and the scale µ r = ω r r/2 is the natural scale for Jj .The evolution of Jj from µ r to µ R sums the logarithms of r/R, and is described by the following modified DGLAP equation, From the NLO expressions in sec.4.2 it follows that the one-loop anomalous dimensions γij are given by γ( 1) 5 Central subjets for the standard jet axis In this section, we discuss the energy distribution of the subjet of radius r centered about the standard jet axis.We start with r R in sec.5.1, which involves a similar calculation as in sec.3.2.In sec.5.2, we discuss the factorization for r R, which takes on a completely different form than in secs.3.4 and 4.2.In particular, the standard jet axis introduces a sensitivity to (the recoil of) soft radiation.We discuss how the double logarithms of r/R can be resummed in sec.5.3.This factorization suffers from non-global logarithms, obstructing an all-orders resummation.

Central subjet function for r R
We start by introducing the function describing subjets centered about the standard jet axis.To distinguish it from the subjet functions of secs. 3 and 4 we denote it by Ĝjet i .We remind the reader that our default notation is for e + e − algorithms, and that the central subjet thus corresponds to a cone of opening angle 2r.On switching to pp algorithms this of course becomes a "cone" in (η, φ) coordinates.It is defined by

Ĝjet
for quark jets, and analogously for gluon jets.There is no sum over subjets j r in the jet, because we now restrict ourselves to the momentum fraction of the central subjet.
The NLO calculation has the same ingredients as in sec.3.2, but the phase-space restrictions for configuration (A) through (C) are modified because we restrict to the subjet centered on the jet axis, The angles β 1 , β 2 and β were given in eqs.(2.5) and (2.6).The contributions (D) and (E) are not modified.There is also a new (and irrelevant) contribution from the configuration where neither parton is inside the central subjet.Performing the calculation, and carrying out the renormalization in the MS scheme, we find for the cone algorithm Ĝcone and for the anti-k T algorithm Ĝanti-k T q (z, z r , ω R , µ) (5.4) The renormalization of the central subjet function is again the same as that of the semiinclusive jet function which enables the resummation of logarithms of R.

Factorization for r R
In the regime r R, the central subjet function contains large logarithms of r/R that require resummation.This is achieved through a second factorization,

Ĝjet
where we made the dependence on r and R explicit in the arguments.We first describe the factorization formula, which differs significantly from eqs. (3.23) and (4.3), and then justify it.
The hard function H ij describes how the energetic parton i produces a jet initiated by parton j with longitudinal momentum ω R and jet radius R, carrying a momentum fraction z of parton i.The collinear function C j describes the fraction z r of collinear radiation produced by parton j, within an angle r of the standard jet axis.It takes into account that the initial collinear parton has a transverse momentum k ⊥ with respect to the jet axis, due to the recoil against the soft radiation, encoded in the soft function S j .The transverse momentum dependence causes the factorization in eq.(5.6) to suffer from rapidity divergences that require regularization.We will employ the η-regulator [68,69], for which ν denotes the corresponding rapidity renormalization scale.Other choices are possible too, see e.g.refs.[70][71][72][73][74]. Mode: Table 1.Modes and power counting in SCET that describe the momentum fraction of the central subjet in the regime r R.
The physical justification of eq.(5.6) is that the hard(-collinear) radiation cannot undergo a perturbative splitting inside the jet.Such a splitting would have a typical opening angle of order R and the contribution of such configurations to the central subjet of radius r R is power suppressed.(Generically, neither of the partons would lie within the central subjets.)Perturbative splittings outside the jet are of course allowed and encoded by the z dependence of the hard function.Collinear splittings inside the jet that affect the central subjet will have typical angle r, and are describe by the collinear function.The (collinear-)soft radiation is not energetic enough to influence the z r measurement, but its transverse momentum k ⊥ affects the jet axis, since the total transverse momentum with respect to the jet axis is zero, and must be taken into account.In the language of SCET, these correspond to distinct degrees of freedom with the parametric scaling of momenta summarized in table 1.
The collinear function has the following definition for j = q, and similarly for j = g.This describes the momentum fraction z r of the central subjet centered on the n axis.The recoil of the collinear radiation with respect to the jet axis due to soft radiation is taken into account through the δ 2 (P ⊥ −k ⊥ ).The definition of the soft function for j = q is given by The delta function sums the transverse momentum k i,⊥ of soft radiation inside the jet, β i < R. Y n is a soft Wilson line in the fundamental representation along the light-like direction n µ = (1, n) of the jet, and Y n is along the opposite direction nµ = (1, −n).For j = g the Wilson lines are in the adjoint representation and the overall normalization is modified 1/N c → 1/(N 2 c − 1).The hard function H ij does not have a direct matrix element definition in SCET, but instead it is defined by the matching relation in eq.(5.6).
The factorization for the standard jet axis in eq.(5.6) does not account for nonglobal logarithms (NGLs) [75,76].These arise because the transverse momentum of the (collinear-)soft radiation inside the jet is probed, but is unconstrained outside the jet. 4  The tree-level hard, collinear and soft functions are given by (5.10) We calculate the one-loop corrections in pure dimensional regularization, such that the virtual corrections are scaleless and vanish.The contributions to H ij come from perturbative splittings of the parton i where the jet consists solely of parton j.For H (1) qq and H (1) qg these can directly be read off from diagrams (D) and (E) in sec.3.2.For the anti-k T algorithm, we find (5.11) Similarly, the results for the cone algorithm can be written as (5.12) We next consider the soft function, which measures the transverse momentum of soft radiation in the jet.Performing the calculation using ref.[79], and noting that the jet region corresponds to rapidity y > − ln(R/2) with respect to the jet axis, we obtain (5.13) 4 Boosting to the frame where the jet becomes a hemisphere, the modes in table 1 become the standard SCETII hard, collinear and soft modes.Emissions into the other hemisphere are unconstrained and lead to additional (collinear-)soft Wilson lines [77,78].In the original frame these corresponds to emissions outside the jet described by Hij.
The result for S g follows by replacing C F → C A .The full collinear function is already complicated at NLO, because it involves two measurements. 5As our current approach is anyway limited to NLL order due to nonglobal logarithms, we simply consider k ⊥ = 0 (from the tree-level soft function).The calculation of the collinear function involves a slight modification to contributions (A), (B) and (C) to the central subjet function in sec.5.1.Since r R, the collinear radiation is close to the center of the jet and does not probe the jet boundary, removing the θ(β < R) and δ(1 − z) in eq.(5.2).For the quark case, the individual contributions are given by Here we needed to include the η-regulator for contribution (B).Adding up the various contributions and performing the renormalization, + 2 P qq (z r ) + P gq (z r ) ln z r . (5.15) A similar calculation yields the gluon collinear function + 2 P gg (z r ) + P qg (z r ) ln z r . (5.16) We have verified that these one-loop ingredients indeed satisfy eq.(5.6).This check involved a subtlety related to distributions: Since a naive expansion of Ĝcone i in r/R leads to improperly regulated plus distributions such as ln(1 − z r )/(1 − z r ) + instead of [ln(1 − z r )/(1 − z r )] + .Rather, we verify eq.(5.6) by first considering z r < 1 and then taking a suitable integral containing the point z r = 1, before expanding in r/R.Note that the difference between Ĝcone i and Ĝanti-k T i is completely captured by eq.(5.12) for r R.

Resummation of ln(r/R)
The resummation is achieved by evaluating each of the ingredients in eq. ( 5.6) at their natural scale and evolving them to common scales µ and ν using their RG equations To avoid a cumbersome complication for resummation in momentum space [81], it is much more convenient to carry out the rapidity resummation in impact parameter space. 6he one-loop anomalous dimensions directly follow from the expressions in sec.5.2,

Winner-take-all axis
The z r averaged results for the SJFs for central subjets along the WTA axis are given by: .
In this case we have a single logarithmic dependence on L r/R .The above expressions hold for both the cone and the anti-k T algorithm, as the algorithm-dependent pieces of Gjet i are contained in J i .
For r R, we have dz z Jj (z , ω r , µ) .(6.5) in direct analogy with eq.(6.2).However, the sum rules in eqs.(2.11) and (6.3) do not hold for Jj and Jij due to the additional theta functions, which is why the single logarithms of r/R persist.

Standard jet axis
Here we present the z r averaged results for subjets along the standard jet axis, showing separate results for the cone and anti-k T algorithm.We obtain (6.6) Note that here we have a double logarithms of L r/R and that these expressions contain power corrections of the form r/R. We have compared this with results available in the literature: Combining the in-jet calculation of ref. [7] (see also refs.[5,6,84] for earlier results obtained within standard QCD) with the out-of-jet contribution of ref. [36], we find agreement with eq.(6.6).Also, for r = R these results reduce to the semi-inclusive jet function in ref. [36].
Our refactorized cross section in the limit r R and, hence, the resummation of logarithms L r/R for the jet shape takes on a different form than in the literature.From eq. (5.6) it follows that where the refactorization of the soft (recoil) contribution is the new ingredient.This additional factorization is essential to resum the logarithms of r/R beyond LL accuracy.

Relation with TMD fragmentation
Because averaging is linear, the jet shape can directly be related to TMD fragmentation through the following sum rule dz r z r Gjet i (z, z r , ω R , µ) = and similarly for the standard jet axis (replacing tildes by hats).This formula describe the central subjet as the sum of the contributions of its hadron constituents.The TMD fragmentation function Gh i (z, ω R , k ⊥ , z h , µ) is the number density of hadrons of species h, momentum fraction z h and transverse momentum z h k ⊥ with respect to the winner-take-all axis [40].The restriction to the central subjet of radius r is encoded by We have verified that eq.(6.8) holds for anti-k T with the winner-take-all axis as well as the standard jet axis using the one-loop results in refs.[40] and [85,86].

Conclusions
In this paper we considered the energy fraction z r of subjets of size r inside a jet of size R.
We presented analytical results for the following three cases: inclusive subjets obtained by reclustering all particles in the jet with jet radius parameter r, as well as central subjets along the winner-take-all axis and along the standard jet axis.The single logarithms of the form α n s ln n R are the same in each case and can be resummed to all orders by solving the associated DGLAP evolution equation.
We also considered the logarithms of the ratio of the jet size parameters ln(r/R), whose structure depends on the particular subjet observable.For each case, we performed an additional refactorization of the cross section in the limit r R, enabling the resummation for this class of logarithms.For central subjets along the WTA axis, this refactorization is known to all-orders in α s but we are currently restricted to leading logarithmic resummation because of our knowledge of the anomalous dimensions.For central subjets along the standard jet axis, an all-orders factorization formula is hindered by non-global logarithms.We presented numerical results for the z r -distribution of inclusive subjets measured on an inclusive jet sample pp → (jet j r )X, leaving numerical results for central subjets to future work [87].In addition, we considered the average energy fraction of these results, which is known as the jet shape for subjets centered on the standard jet axis.For the jet shape, our factorization formula in the limit r R involves an additional refactorization compared to the literature, to account for the recoil effect of soft radiation on the jet axis.Along the way, we also pointed out an inconsistency in the literature for analytical results for inclusive cone jets and their substructure.
There are various possible applications of our work in the future.First of all, it will be very interesting to perform numerical calculations for central subjets along the two axes we considered in this work.This will be particularly relevant since our factorization for the standard jet shape in the limit r R differs from that in the available literature at next-to-leading logarithmic order, and it is possible to compare to experimental data in this case.In addition, for the more differential case (the z r -dependent case), experimental measurements will be feasible which can shed new light on the substructure of jets at the LHC.In addition, it will be interesting to further explore the possibility of how these "relatively inclusive" jet substructure observables can be used to discriminate between QCD jets and boosted objects.Finally, we expect that the different energy distributions of subjets considered in this work can be very relevant to better understand the properties of the quark-gluon plasma created in heavy-ion collisions.

Figure 2 .
Figure 2. The three configurations that enter for the quark semi-inclusive jet function at O(α s ): (A) the quark and gluon are inside the jet, (B) only the quark is inside the jet, (C) only the gluon is inside the jet.

3 Figure 6 .
Figure 6.Left: Same as fig.5but without multiplying the result for the different p T intervals by multiples of 10.Right: Subjet distribution measured on an inclusive jet sample, using the same kinematics as in fig.5and the p T bin of[50, 100]  GeV.The distribution is shown for different values of the subjet radius: r = 0.05 (green dotted), r = 0.1 (black dashed), r = 0.2 (blue dot-dashed) and r = 0.3 (red solid).

Figure 7 .
Figure 7.The subjet distribution measured on inclusive jets as a function of r/R, using the same kinematics as in figs.5 and 6, and the jet p T bin [25, 50] GeV (left) and [100, 200] GeV (right).Four representative values of the ratio z r are shown: 0.01 (green dotted), 0.1 (black dashed), 0.8 (blue dot-dashed) and 0.95 (red solid).