Jet fragmentation functions in proton-proton collisions using soft-collinear effective theory

The jet fragmentation function describes the longitudinal momentum distribution of hadrons inside a reconstructed jet. We study the jet fragmentation function in proton-proton collisions in the framework of soft-collinear effective theory (SCET). We find that, up to power corrections, the jet fragmentation function can be expressed as the ratio of the fragmenting jet function and the unmeasured jet function. Using renormalization group techniques, we are able to resum large logarithms of jet radii R in the perturbative expansion of the cross section. We use our theoretical formalism to describe the jet fragmentation functions for light hadron and heavy meson production measured at the Large Hadron Collider (LHC). Our calculations agree very well with the experimental data for the light hadron production. On the other hand, although our calculations for the heavy meson production inside jets are consistent with the PYTHIA simulation, they fail to describe the LHC data. We find that the jet fragmentation function for heavy meson production is very sensitive to the gluon-to-heavy-meson fragmentation function.


Introduction
Collimated jets of hadrons are a dominant feature of high energy particle interactions, especially at the current highest energy hadron collider, the Large Hadron Collider (LHC), where jets are abundantly produced. The internal structure of these jets has become an important tool to test the fundamental properties of Quantum Chromodynamics (QCD), and to search for new physics beyond the Standard Model [1,2]. Needless to say, a good understanding of jet substructure allows deeper insights into QCD dynamics and serves as a prerequisite for further progress.
One of the jet substructure observables proposed and explored in more detail recently is the jet fragmentation function, which describes the longitudinal momentum distribution of hadrons inside a reconstructed jet [3][4][5][6][7][8][9][10][11][12][13]. Experimental studies on hadron distribution inside jets have been pioneered at the Tevatron [14] in the 1990s. More recently, both the ATLAS and the CMS collaborations have measured the distributions of light hadron [15][16][17][18] and heavy meson [19] production inside jets at the LHC. The jet fragmentation function is an interesting and important observable: since it probes the hadron fragmentation at a more differential level, it can reveal detailed information about the jet dynamics involved in producing the identified hadron. At the same time, it can provide further information about the non-perturbative hadronization encoded in the standard fragmentation functions. One might even gain insight into the nontrivial spin correlation through the study of azimuthal distribution of the hadron inside jets [20][21][22][23][24].
Since gluon jets are much more abundant in proton-proton collisions at high energy hadron colliders, jet fragmentation functions should be more sensitive to gluon fragmentation. We will show that this is the case especially for heavy meson production inside jets. This situation is very different from the e + e − → h X and e p → e h X processes, where the gluon fragmentation function does not enter at leading-order in the perturbative calculation and, thus, can only be probed through QCD evolution or higher-order radiative corrections.
There is also strong motivation to study the jet fragmentation function in heavy ion collisions at high energies, where hot and dense QCD medium -the quark-gluon plasma -is produced. By comparing the jet fragmentation function measured in ultra-relativistic nucleus collisions and the one in proton-proton collisions, one can understand how the presence of the strongly interacting medium produced in heavy ion collisions modifies the hadron distributions inside jets. Understanding the light and heavy flavor dynamics in the medium will help further determine the precise properties of the QGP. For recent experimental measurements of the jet fragmentation function in heavy ion collisions at the LHC, see [16][17][18]. For some theoretical work along this direction, see [25][26][27].
In this paper, we study the jet fragmentation function in proton-proton collisions using soft-collinear effective theory (SCET) [28][29][30][31][32]. Previously, in [10,13] a full next-to-leading order (NLO) calculation was performed. Closely related work with emphasis on heavy flavor was also recently presented in [33,34]. As we will show below, within SCET the hadron distribution inside jets is governed by the ratio of two quantities: the fragmenting jet function (FJF) G h i (ω, R, z, µ) introduced and studied in [3][4][5][6][7][8], and the unmeasured jet function J i (ω, R, µ) introduced in [35]. Here, i is the parton that initiates the jet with energy ω and radius R, while z is the fraction of the jet momentum carried by the identified hadron h. The FJF G h i (ω, R, z, µ) can be further written as a convolution of perturbatively calculable Wilson coefficients J ij and the fragmentation functions D h j (z, µ). Using the renormalization group techniques, we are able to simultaneously resum logarithms of the form ln R and ln(1 − z), which have a significant numerical impact. Such resummations were not addressed previously in the fixed NLO calculation of [13]. We use the formalism to describe the experimental data at the LHC for the distribution of light hadron and heavy meson production inside jets. The study of the jet fragmentation function in heavy ion collisions using SCET will be performed in a forthcoming paper [36]. Some of the input for this calculation, such as the final-state in-medium splitting functions [37] and mediummodified fragmentation functions applied to leading hadron production [38,39], are already available.
Here, we would like to remind the readers that, although the jet fragmentation function and the fragmenting jet function look very similar, they have different meanings. It is important to understand their differences and relations since they appear throughout the entire paper. The jet fragmentation function is an experimental observable describing the distribution of hadrons inside jets. On the other hand, the fragmenting jet function is a theoretical quantity which enters the factorized expression in the calculation of the jet fragmentation function. See Sec. 2 for more details.
The rest of the paper is organized as follows. In Sec. 2, we first provide the defini-tion of the jet fragmentation function. We then derive a factorized expression for the jet fragmentation function, which involves the FJF and the unmeasured jet function. We give the matching coefficients for the FJF to be convolved with the standard fragmentation functions, and in particular for jets reconstructed using the anti-k T jet algorithm, which is used in almost all jet reconstruction at the LHC. We collect the detailed derivations of the matching coefficients in the Appendix Sec. A. In Sec. 3, we present the numerical results of our calculations for light hadron and heavy meson production inside jets and compare with the experimental data at the LHC. We also explore the theoretical uncertainty, the sensitivity of the observable to the jet algorithm (either cone or anti-k T ), and the radius dependence. We summarize our paper in Sec. 4.

Jet fragmentation function
In this section we give the definition of the jet fragmentation function and calculate it using the factorized expression in SCET. The evaluation involves the fragmenting jet function G h i (ω, z, R, µ), and we provide the Wilson coefficients J ij to be convolved with the fragmentation function D h j (z, µ). We give the results for jets reconstructed using cone and anti-k T algorithms, as J ij depends on the jet algorithm. The results for cone jets are available in [7], while those for anti-k T jets were first written down in the appendix of [40]. We provide the detailed derivations of J ij for anti-k T jets in the appendix, and the results are consistent with [40].

Observable and factorized expression
The jet fragmentation function F (z, p T ) describes the longitudinal momentum distribution of hadrons inside a reconstructed jet. We will compare our calculations with the jet fragmentation functions measured in proton-proton collisions, p + p → (jet with h) + X. Here, F (z, p T ) is defined as follows, where dσ h /dydp T dz and dσ/dydp T are the differential cross sections of jets with and without the reconstruction of the hadron h in the jet. Here, y and p T are the jet rapidity and transverse momentum. z is the fraction of the jet transverse momentum carried by the hadron, z ≡ p h T /p T , with p h T the transverse momentum of the hadron. Jets are reconstructed using either the cone or the anti-k T algorithm with the jet radius R, and the R-dependence is suppressed in the expression for F (z, p T ). As we will see, jet fragmentation functions will be different for jets reconstructed using different jet algorithms.
Because the contribution from the soft radiation to the longitudinal momentum is power suppressed [41], it suffices to illustrate the SCET factorized expression for the jet fragmentation function in e + e − collisions (Fig. 1). Following [3,7,12,35,41,42], the differential cross section for N -jet production with the jet p T i and y i , the hadron h inside Figure 1. Illustration of the N -jet production in e + e − collisions, where a hadron is measured in the jet labeled by J 1 with rapidity y and transverse momentum p T . z is the fraction of the jet momentum carried by the hadron. Jets are reconstructed using a jet algorithm with radius R. We impose an energy cutoff Λ outside the jets to ensure the N -jet configuration. Λ is a low energy scale constraining the soft radiation (red lines). The green lines represent the collinear splittings. one jet (labeled by 1), and the energy cutoff Λ outside all the jets can be written as follows, where H(y i , p T i , µ) is the hard function describing the short-distance production of the N jets with rapidities y i and momenta p T i . S n 1 n 2 ···n N (Λ, µ) is the soft function with N soft Wilson lines along the jet directions. The energy cutoff Λ outside the jets is imposed to ensure the N -jet configuration. The hadron h measured inside jet 1 is described by the FJF G h ω 1 (z, µ), with the jet radius R suppressed. J ω i (µ) (for i = 2, · · · , N ) are the unmeasured jet functions introduced in [35] 1 , with ω i representing the large light-cone component of the jet momentum and ω i = 2p T i in the frame where the jet is in the transverse direction. The factorized expression is valid for collimated jets up to power corrections of the type Λ/Q or R.
On the other hand, the differential cross section for N -jet production is given by with the same hard function H(y i , p T i , µ), soft function S n 1 n 2 ···n N (Λ, µ), and unmeasured jet functions J ω i (µ) with i = 2, · · · , N . The only difference is that G h ω 1 (z, µ) in Eq. (2.2) is replaced by the unmeasured jet function J ω 1 (µ) in Eq. (2.3) since we do not measure the hadron. The distribution of the hadron h inside jet 1 then becomes, All the hard, soft and unmeasured jet functions (except for jet 1) cancel in the ratio. Taking the average over the jet production cross section, with proper phase space (PS) cuts on both jet rapidity y and transverse momentum p T , e.g. the rapidity interval and the width of the p T bin, the jet fragmentation function F (z, p T ) becomes where dσ i /dy dp T is the cross section to produce the jet initiated by parton i, and we have written out explicitly the arguments for both the FJF G h i (ω, R, z, µ) and the unmeasured jet function J i (ω, R, µ).
In the next subsection we will provide explicit expressions for the fragmenting jet function G h i (ω, R, z, µ) and the unmeasured jet function J i (ω, R, µ). Here it is instructive to point out that G h i (ω, R, z, µ) and J i (ω, R, µ) have the same renormalization group (RG) evolution [5,7,35] and the ratio G h i (ω, R, z, µ)/J i (ω, R, µ) is renormalization group invariant, with possibly different characteristic scales for G h i and J i .

Unmeasured jet function
For convenience, we provide all the relevant results for the unmeasured jet function J i (ω, R, µ).
At O(α s ) [35], where L = ln ω tan (R/2) µ , (2.8) and d q/g,alg J represents the algorithm-dependent pieces, The unmeasured jet function J i (ω, R, µ) satisfies the RG equation with the anomalous dimension given as follows: (2.14) Here, Γ i cusp and γ i are the cusp and non-cusp anomalous dimensions, with the perturbative with T F = 1 2 , n f the number of active quark flavors, and The solution of the RG equation for the unmeasured jet function is where µ J is the characteristic scale of J i (ω, R, µ), which eliminates the large logarithms in the fixed-order calculation. From Eqs. (2.6) and (2.7), the choice of µ J ∼ ω tan (R/2) ≡ p T R eliminates the logarithm L. We denote this scale as "p T R " for later convenience.

Fragmenting jet function
The fragmenting jet functions G h i (ω, R, z, µ) [5,7,11] are closely related to the fragmentation functions D h j through matching coefficients J ij where D h j (z, µ) is the fragmentation function of a parton j fragmenting into a hadron h. Eq. (2.19) is valid for a light hadron h up to power corrections of order Λ 2 QCD /ω 2 tan 2 (R/2). Thus, to avoid large non-perturbative power corrections, R should not be too small. On the other hand, for heavy meson fragmenting jet junction Λ QCD should be replaced by the heavy quark mass m Q in the above equation [12].
The Wilson coefficients J ij depend on the jet algorithm. The results for cone jets were given in [7], while those for anti-k T jets were first written down in the appendix of [40]. We provide the detailed derivations of J ij for anti-k T jets in the appendix, and the results are consistent with [40]. Here we only list the final results: where the functionsP ji have the following expressions [6] J alg ij (z) represent pieces that depend on the jet algorithm. For cone jets [7], The fragmenting jet function G h i (ω, R, z, µ) satisfies the following RG equation where the anomalous dimension γ i G (µ) = γ i J (µ) is the same as that of the unmeasured jet function J i (ω, R, µ) [5,7,35] in Eq. (2.14). The solution to the RG equation is where the scale µ G should be the characteristic scale that eliminates the large logarithms in the fixed-order perturbative calculations. In the large z region, the scale choice µ G = ω tan (R/2) (1 − z) ≡ p T RZ resums [7] both ln R and ln (1 − z). However, for consistency, this would require extracted fragmentation functions D h j with a built-in resummation of logarithms in (1 − z), which is currently not available. It might be instructive to point out that with such a scale, the power corrections in Eq. (2.19) will be of the order of Λ 2 QCD / ω 2 tan 2 (R/2)(1 − z) 2 , similar to the usual threshold resummation, see, e.g. Ref. [44]. For the numerical calculations presented in the next section, we will choose µ G = ω tan (R/2) to resum ln R and comment on the effect of ln (1 − z) resummation.
Let us make a few comments about our resummation formalism. As we have pointed out already at the end of Sec. 2.1, since G h i (ω, R, z, µ) and J i (ω, R, µ) follow the same RG evolution equations, as given in Eqs. (2.37) and (2.18), respectively, the ratio as given in the factorized formalism Eq. (2.5) is thus RG invariant. In other words, this ratio does not depend on the scale µ. Choosing µ G = µ J = ω tan(R/2), the whole RG exponential forms cancel in the ratio. However, this does not mean that resummation effects disappear in our framework. On the contrary, the resummation effect is shifted entirely into the scale µ G -dependence of the standard fragmentation function D h i (z, µ G ) through Eq. (2.19). In other words, we are resumming ln(R) logarithms in this case through the DGLAP evolution equations of the fragmentation functions. This type of resummation was not achieved previously in the fixed NLO calculation of [13]. It will be very interesting to explore the exact relation between our work and the previous NLO calculation [13], which we are going to address in a future publication.

Phenomenology
In this section, we present the numerical results of our theoretical formalism and we compare our calculations with the experimental data for both light hadron and heavy meson production at the LHC. We will also explore the theoretical uncertainties of our formalism.

Light hadron jet fragmentation function
We first study the distribution of light hadrons inside jets in proton-proton collisions. Both ATLAS and CMS collaborations at the LHC have measured the distribution of light, charged hadrons h = h + + h − inside jets. We perform the numerical calculations using the CT14 NLO parton distribution functions [48] and the DSS07 NLO fragmentation functions [49,50]. We keep the Γ i 0,1 and γ i 0 terms in the series expansion of the anomalous dimension γ i J,G with i = q, g. Therefore the calculation is at next-to-leading logarithmic accuracy.
In Fig. 2, we compare our calculations with the experimental data from ATLAS [15] in proton-proton collisions at the center-of-mass (CM) energy of √ s = 7 TeV. Jets are reconstructed using the anti-k T algorithm with R = 0.6 within the rapidity range |y| < 1.2. The transverse momenta p T of jets are measured across a wide range, from 25 GeV to 500 GeV. The numbers in square brackets correspond to different jet transverse momentum bins, e.g. [25,40]  while the solid blue curves are the "nominal" theoretical calculations, where we make the scale choices µ = p T , µ G = µ J = p T R defined in the last section. The green bands are the estimated uncertainties of our theoretical calculations from the scale variations for all three scales µ, µ G , µ J by a factor of 2 around the above central values. See detailed discussions in Sec. 3.2 below. Note that the DSS07 fragmentation function parameterizations for D h i (z, µ) are only valid for 0.05 < z < 1 and 1 < µ 2 < 10 5 GeV 2 . Thus, all the calculations outside these regions are based on the extrapolations of the DSS07 parameterizations provided by the distributed package from the authors [49,50]. As we have expected, the theoretical uncertainties from the scale variations are relatively small, due to the fact that G h i and J i follow the same RG running as discussed in Sec. 2.1. At the same time, as one can see, there is good agreement between our theoretical calculations and the ATLAS data. Our calculations slightly overshoot the experimental data at large z for jets with low p T . Since there are large uncertainties for fragmentation functions in the large z region [51,52], jet fragmentation function measurements in proton-proton collisions can help constrain them in this region. The solid red circles are the ATLAS preliminary data [53], while the magenta triangles are the CMS data [17]. The blue solid curves are the "nominal" theoretical calculations, with the green bands representing the theoretical uncertainties estimated from scale variations.
In Fig. 3, we compare our calculations with the preliminary ATLAS data [53], as well as the CMS measurements [17] in proton-proton collisions at the CM energy √ s = 2.76 TeV.
Here, jets are reconstructed using the anti-k T algorithm with R = 0.4 within the rapidity range |y| < 1.6 for ATLAS, whereas for CMS R = 0.3 and 0.3 < |y| < 2. The solid red circles are the ATLAS data, while the magenta solid triangles are the CMS data. As one can see, our calculations agree with the data rather well. Note that the CMS data has a very different trend for low z 0.05 compared to the ATLAS data. Our theoretical predictions in Figs. 2 and 3 also agree with the results in [13] that use the full NLO calculation.

Algorithm and radius dependence, and theoretical uncertainty
Here, we study the dependence of the jet fragmentation function on the jet algorithm and the jet radius. We will also estimate the theoretical uncertainty by varying the characteristic scales in our formalism.  We choose the scales µ = p T and µ G = µ J = p T R . The solid red curve is for anti-k T jets, while the dashed blue curve is for cone jets. Lower panel: the ratio of the jet fragmentation functions F (z, p T ) cone /F (z, p T ) kT for cone and anti-k T jets.
We will first explore the jet algorithm dependence. In the upper panel of Fig. 4, we plot the jet fragmentation function F (z, p T ) for light charge hadrons as a function of z inside jets with 60 < p T < 80 GeV, |y| < 1.2, R = 0.6 at √ s = 7 TeV as an example. We choose the scales µ = p T and µ G = µ J = p T R . The solid red curve is for anti-k T jets, while the dashed blue curve is for cone jets. As we can see from this plot, F (z, p T ) for cone jets is smaller (larger) than that for anti-k T jets at large (small) z. This is a consequence of two combined effects: in the low z region, the FJF G h i for cone jets is larger than that for an anti-k T jet. As z gets closer to 1, the FJF G h i for cone and anti-k T jets approach the same value because there is little radiation left in the jet to distinguish between jet algorithms. Also, the unmeasured jet function J i for cone jets is larger than that for anti-k T jets. To see the difference more clearly, we plot the ratio F (z, p T ) cone /F (z, p T ) k T between the jet fragmentation functions for cone and anti-k T jets in the lower panel of Fig. 4. We now study the jet radius R dependence. We choose the scales µ = p T and µ G = µ J = p T R . In Fig. 5, we plot as an example the jet fragmentation functions F (z, p T ) as a function of z for four different jet radii R = 0.2 (solid red), R = 0.4 (dashed blue), R = 0.6 (dotted black), and R = 0.8 (dash-dotted magenta) for jets with 60 < p T < 80 GeV, |y| < 1.2 at √ s = 7 TeV. We find that in the large z 0.1 region F (z, p T ) gets smaller as R increases. In the small z 0.1 region, F (z, p T ) becomes larger as R increases because of the normalization of F (z, p T ). This is related to the scale dependence of D h i (z, µ G ), which is governed by the DGLAP evolution equations: D h i (z, µ G ) increases (decreases) as µ G increases for small (large) z [54]. Since µ G = p T R = 2p T tan(R/2), increasing R will increase µ G .
Finally, we estimate the uncertainty of our theoretical calculations by varying the scales µ, µ J , and µ G around We independently vary the scales by a factor of 2 around their central values, i.e.,   Fig. 6, while the blue dashed curve represents the central value with µ = p T , µ G = µ J = p T R . The uncertainty of our calculations is generally small for the moderate z region, and it is compatible with the results based on the full NLO calculation in [13], where only the variation of the scale µ is implemented. This gives us confidence that the RG evolutions for both the FJF G h i and the unmeasured jet function J i indeed improve the convergence of the theoretical calculation.
When z gets closer to 1, one can see that the scale uncertainty band becomes larger. As we have shown in the last section, there is an explicit dependence in the FJF G h i on ln (1 − z). These logarithms become large as z approaches 1, i.e. in the hadronic threshold limit. We may [7] simultaneously resum logarithms of the jet radius R and (1 − z) by choosing the scale µ G ∼ 2p T tan (R/2) (1 − z) ≡ p T RZ . We plot this by independently varying the scales as follows, Such scale variations correspond to the green band in the upper panel of Fig. 6, while the red dashed central curve represents the calculation with µ = p T , µ J = p T R , µ G = p T RZ . As one can clearly see, the uncertainty of the calculation with ln(1 − z) resummation is largely reduced in the large z region. In order to see the effect of ln(1 − z) resummation more clearly, in the lower panel of Fig. 6, we plot the ratio R (p T RZ , as a function of z and we set µ = p T and µ J = p T R . As one can see, resumming ln (1 − z) leads to an enhancement of the jet fragmentation function F (z, p T ) in the large z region. For z 0.8, the enhancement is about a factor of 2. Even though the theoretical uncertainty is reduced with the scale choice µ G = p T RZ , we do not use this scale when comparing to data in Figs. 2, 3 above and Fig. 7 below. This is due to the fact that the fragmentation functions that we use in our numerical studies are extracted using fixed-order calculations [49][50][51][52]. In order to be consistent, we have to adopt the conventional scale choice µ G = p T R . However, we want to make an important point. If one performs a fit for fragmentation functions using the F (z, p T ) data, the extracted functions would differ significantly in the large z region when the more accurate calculation with ln (1 − z) resummation is used. Our conclusions here are similar to the observations made in [55] in the context of threshold resummation.

Heavy meson jet fragmentation function
Our theoretical result in Eq. (2.5) was derived for light hadron production inside jets. However, it can also be applied to describe heavy meson production inside jets using the Zero Mass Variable-Flavor Number Scheme (ZMVFNS) [56,57]. Such a scheme applies when the perturbative scales Q are much larger than the heavy quark mass m Q : Q 2 m 2 Q . In this kinematic regime, the heavy quarks are expected to behave like massless partons. One can, thus, treat heavy quarks as the other light partons, and logarithms associated with m Q are resummed using the DGLAP evolution. Power corrections of O(m 2 Q /Q 2 ) are neglected in this formalism. In our case, the ZMVFNS applies in the kinematic regime where µ, µ J , µ G m Q . The ATLAS collaboration has recently measured the distribution of D * ± -mesons in jets with p T > 25 GeV and R = 0.6 [19]. Given the fact that the charm mass is relatively small m c ∼ 1.3 GeV [58], the jet transverse momentum is large and the radius is moderate, this satisfies the requirement for using the ZMVFNS.
Within the ZMVFNS, the only change in our theoretical formalism is to also include the charm production in Eq. (2.5): i=q,g,c with q and c representing light and charm flavor, respectively. Like in light hadron calculations, we make the scale choices µ = p T , µ G = µ J = p T R for the "nominal" calculations. We follow Sec. 3.2 to calculate the theoretical uncertainties from the scale variations. We use the charm-meson fragmentation functions extracted from the inclusive production of a single charm-meson D in e + e − collisions: e + e − → D X. The parameterizations for D h i (z, µ) with i = q, g, c and h = D are available in [59], which yield a good description of the inclusive D-meson production in proton-proton collisions at the LHC [60]. Thus, we will use this parametrization in our calculations.  . The calculation of jet fragmentation functions for D * ± meson production compared to the experimental data from the ATLAS collaboration at √ s = 7 TeV [19]. Jets are reconstructed using the anti-k T algorithm with R = 0.6, and the jet rapidity is within |y| < 2.5. We show 6 different panels which correspond to different jet p T ranges. The solid blue circles are the experimental data measured by ATLAS [19], while the empty red circles are the PYTHIA simulations provided in the ATLAS paper [19]. The solid red curves are our default theoretical calculations using the ZMVFNS. The green bands are the estimated theoretical uncertainties from the scale variations. The dashed blue curves are our calculations using an enhanced gluon-to-D meson fragmentation function: In Fig. 7, we compare our calculations for the D * ± jet fragmentation function with the ATLAS experimental data at CM energy √ s = 7 TeV [19]. Jets are reconstructed using the anti-k T algorithm with R = 0.6, and the jet rapidity is within |y| < 2.5. We show 6 different panels which correspond to different jet p T ranges covering 25 < p T < 70 GeV. The solid blue circles are the experimental data measured by ATLAS [19] and the empty red circles are the PYTHIA simulations provided in the ATLAS paper [19].
The solid red curves are our default theoretical calculations, which use the central values of the D-meson fragmentation functions D h i (z, µ) from [59]. The green bands are the theoretical uncertainties estimated from the scale variations. As one can clearly see, our theoretical calculations are consistent with the PYTHIA simulations for all different jet p T bins. However, they are significantly below the experimental measurements from the ATLAS collaboration.
As we have mentioned, the D-meson fragmentation functions are extracted in e + e − collisions, where the gluon fragmentation function D D g (z, µ) does not enter at leadingorder in the theoretical formalism. Thus, gluon fragmentation is only indirectly probed through QCD evolution and/or higher-order corrections. This leads to a large uncertainty of the extracted gluon-to-D meson fragmentation function. Note that Ref. [59] does not provide the uncertainty of the extracted charmed-meson fragmentation functions. However, comparing different extractions from the same sets of e + e − data [59,61,62], we find that the gluon-to-D meson fragmentation function D D g (z, µ) can differ by a factor of 3, while quark-to-D meson fragmentation functions D D q,c (z, µ) do not vary so dramatically [59]. Other than that, the various extractions [59,61,62] differ only by the initial scales for the QCD evolution or by the treatment of the heavy quark mass. This provides a strong hint that the current extraction of the gluon-to-D meson fragmentation function could have a very large uncertainty.
To explore the uncertainty of the gluon-to-D meson fragmentation function, we reperform our calculations of the jet fragmentation functions for D * ± meson with the gluonto-D meson fragmentation function enhanced by a factor of 2, i.e. D D g (z, µ) → 2 D D g (z, µ). These calculations are shown by the dashed blue curves in Fig. 7. They lead to much better agreement with the ATLAS data. We have also tried enhancing other quark-to-D meson fragmentation functions D D q,c (z, µ) by a similar factor, but none of them could lead to such an efficient enhancement in the jet fragmentation function. We conclude that jet fragmentation functions of heavy mesons in proton-proton collisions have great potential to constrain the gluon-to-heavy meson fragmentation functions.

Summary
In this paper, we studied jet fragmentation functions for light hadrons and heavy mesons inside reconstructed jets. We wrote down a factorized expression in SCET for the jet fragmentation function in proton-proton collisions. We found that, up to power corrections, the jet fragmentation function can be expressed as the ratio of the fragmenting jet function and the unmeasured jet function. These two functions satisfy the same renormalization group equation, and the fragmenting jet function can be further expressed as a convolution between the fragmentation functions and the matching coefficients. Using SCET, we were able to simultaneously resum large logarithms of the jet radius R and (1 − z), which has a significant impact on the phenomenology considered in this work. We used the theoretical formalism to describe the jet fragmentation functions for light hadron and heavy meson production measured at the LHC. We found that our calculations agree very well with the experimental data for light hadron production. We explored the jet algorithm and the R dependence of the jet fragmentation functions, and we estimated the theoretical uncertainty by scale variation. For heavy meson production inside jets, although our calculations are consistent with PYTHIA simulations, they fail to describe the corresponding LHC data. We found that enhancing the gluon-to-heavy meson fragmentation function leads to much better agreement with the experimental data. We emphasize that the jet fragmentation function for heavy meson production in proton-proton collisions is very sensitive to the gluon-to-heavy meson fragmentation function. In the future, we plan to extend our calculations to describe jet fragmentation functions in heavy ion collisions in order to understand nuclear modifications of hadron production inside jets.
with z = q − /ω. The FJF G h i (ω, R, z, µ) is therefore related to the fragmenting jet functions G h i (m 2 J , z, µ) [5] with the extra measurement of jet mass through where δ alg = δ cone or δ anti-k T are the constraints given in Eqs. (A.3) and (A.4). The FJF G h i (ω, R, z, µ) can be matched onto the fragmentation function D h i (z, µ): and J ij are the matching coefficients. The FJF G j i (m 2 J , z, µ) with i, j = q, g has been extensively studied in [5,11]. Using pure dimensional regularization with 4 − 2 dimensions in the MS scheme, the bare results at O(α s ) can be written in the following compact form [11,63]: where the functions P ji (z, ) are Substituting Eq. (A.7) into Eq. (A.5) and performing the integration over m 2 J with the constraints imposed by the jet algorithm δ alg , one obtains the bare FJF G j i,bare (ω, R, z). We present the results for anti-k T jets here, as the explicit expressions are not available in the literature: +P qq (z) (L + ln z) (A.14) +P gg (z) (L + ln z) where, as a reminder, β 0 and L are given by 16) andP ji have the expressions [6] given in Eq.  15) are infrared (IR), and they match exactly those in the fragmentation functions D j i as we will show below. The renormalization of G h i,bare (ω, R, z) is given by where i is not summed over on the right hand side. The corresponding renormalization group (RG) equation is given by where the anomalous dimension γ i G (µ) is The solution to Eq. (A.18) is then where the scale µ G should be the characteristic scale chosen such that large logarithms in the fixed-order calculation are eliminated. The counter terms Z i G (µ) are given by From these results we obtain the anomalous dimension γ i G (µ) with the following form, and higher-order results can be found in [5,[44][45][46][47].