The Fragmentation Spectrum from Space-Time Reciprocity

Analyzing the single inclusive annihilation spectrum of charged hadrons in $e^+e^-$ collisions, I confront the hadronization hypothesis of local parton-hadron duality with a systematic resummation of the dependence on the small energy fraction. This resummation is based on the reciprocity between time-like and space-like splitting processes in $4-2\epsilon$-dimensions, which I extend to resum all the soft terms of the cross-section for inclusive jet production. Under the local-parton-hadron duality hypothesis, the resulting distribution of jets essentially determines the spectrum of hadrons as the jet radius goes to zero. Thus I take the resummed perturbative jet function as the non-perturbative fragmentation function with an effective infra-red coupling. I find excellent agreement with data, and comment on the mixed leading log approximation previously used to justify local parton-hadron duality.


Introduction
Local Parton-Hadron Duality (LPHD) posits that the momentum spectrum of partons, when calculated perturbatively, essentially determines the momentum spectrum of the resulting hadrons [1], so long as there is an underlying hard process at a scale well-above confinement. LPHD is consistent with the much more rigorous statement of factorization, which states that physics at different scales is not entangled beyond the anomalous dimensions of scaling quantities. Within QCD factorization, factorized sectors can be treated independently up to power corrections, and up to consistency conditions on their divergences. In a factorization theorem, any sector sensitive to scales where the running coupling is perturbative can be legitimately calculated in perturbation theory. In contrast, when taken in its strong form the LPHD hypothesis goes well beyond the statement of factorization, where one wishes to apply perturbation theory as close to the non-perturbative regime as possible. This is to say that the perturbative results for infra-red non-perturbative functions are to be the model for such functions. This can be illustrated with the classic observable of fragmentation: the total energy Q of the collision is broken up into various fractions, which are distributed amongst the asymptotic remnants of the hard interaction, labeled by h: (1.1) Within factorization, one constructs the fragmentation spectrum as the convolution of a non-perturbative boundary condition, the fragmentation function, with a perturbatively calculable renormalization group evolution, and a perturbatively calculable coefficient function. Within LPHD, up to questions about the scale of confinement which ends all colored interactions between fields, a purely perturbative calculation should also determine the shape of the fragmentation function. This is because the hadronization process itself is hypothesized to be unable to significantly change the pre-hadronization momentum distribution. Historically, this hypothesis when combined with a resummation of soft physics in the leading order Dokshitzer-Gribov-Lipatov-Alterelli-Parisi (DGLAP) equations led to the so-called "mixed leading log" approximation (MLLA), for review, ref. [2]. From the basic equations for the MLLA as applied to fragmentation, 1 two models may be derived that have seen phenomenological success: the distorted Gaussian approach of ref. [4,5], 2 and the so-called limiting spectrum of ref. [7] (and references therein). Both are intended to describe the hadron fragmentation spectrum in the region where the bulk of hadrons are produced, as the energy fraction vanishes. Both approaches however have difficulties in describing the full small energy fraction region while consistently extrapolating between different center-of-mass energies for the e + e − event. 3 Moreover, it is unclear how exactly the MLLA relates to genuine perturbative QCD calculations of the MS-splitting functions of QCD. It is clear that the leading order DGLAP kernel acts as the leading order kernel for the MLLA equations. None-the-less the MLLA derived anomalous dimension simply differs from the next-to-leading logarithmic MS anomalous dimension from small-z h resummed perturbation theory of refs. [10,11]. While this does not automatically invalidate the MLLA, one cannot simply match its results to known QCD perturbation theory beyond leading log.
Expanding upon a previous paper, ref. [12], I can construct a resummation of the soft terms in the cross-section that does not suffer from these limitations and derive a hadron fragmentation spectrum fully consistent with the standard perturbation theory results. I can precisely trace all quantities needed to define the hadron fragmentation spectrum to the space-like DGLAP anomalous dimensions in 4 − 2 -space-time dimensions. I do so by approaching the problem of hadrons as a problem of resumming inclusive jet production ( [13][14][15][16]) when the radius of the jet approaches 0, an idea initially explored in ref. [17]. Unlike the fragmentation function, the semi-inclusive jet function is perturbatively calculable for small but not extreme jet radii. However, one needs the jet function beyond the first few orders of perturbation theory, since DGLAP evolution drives all momentum fractions towards zero. Thus the jet function is dominated by its behavior at small energy fractions as the radius of the jet goes to zero. After I accomplish this resummation, I can invoke the LPHD hypothesis to derive the model I need for the non-perturbative fragmentation function: the jet function evaluated with a frozen coupling scale in the infrared. This is analogous to replacing the perturbative coupling with the effective coupling of ref. [18], see also ref. [19].
With this model of fragmentation, I can then fit the measured spectrum of charged hadrons with 3 parameters: the scale at which I freeze the coupling, the scale at which I end DGLAP evolution, and the basic normalization. I will find an excellent quantitative description of data from the Large Electron Positron Collider (LEP) and the Stanford Linear Accelerator at the Z-pole for essentially all energy fractions below ∼ 0.5. A good description is also obtained for older data at smaller center-of-mass energies, though hadron mass-corrections and issues of normalization limit the description of the extremely small-z h range.
Fits to the small-z h pion data have been made in ref. [20] using a systematic resummation for the time-like splitting functions and the coefficient function derived in [10,11,21], where the derived resummation for these quantities is equal to that obtained in ref. [12]. But the corresponding resummation for the fragmentation function was not attempted. Rather, a generic model for the fragmentation function was used, with success at describing the data.
The organization of the paper is as follows. First, I extend the formalism of ref. [12], and show how to resum the semi-inclusive jet function/fragmentation function. Given the form of the resummation, I can expect large cancellations between the resummation of the coefficient function and the jet function, which naturally suggests an organization of resummed spectrum. I then compare this arrangement to the MLLA approximation. Finally, I give the non-perturbative model for fragmentation suggested by LPHD and the resummed perturbation theory, and compare to data. Then I conclude.
I note that all results in this paper equally apply to semi-inclusive jet production, and in what follows I will always call the low-scale semi-inclusive jet function simply the "fragmentation function" to avoid the somewhat cumbersome "semi-inclusive jet function." Ultimately, the only difference between fragmenting to a jet and fragmenting to a hadron is the scale at which the evolution is stopped and the collinear function is evaluated, the functional form in lnz h being identical.

Resumming the fragmentation function
In this section I will explain how to resum the fragmentation function.
I note the conventions for moment transforms as: The constant c is chosen to place all singularities off to the left of the contour. The soft region of z h → 0 corresponds to n → 0. The log counting I adopt is α s ∼ n 2 , and the order of the pole in n at 0 corresponds to the number of logarithms in momentum space: I start with the basic angular-ordered factorization posited in ref. [12], outline its connection to the standard formulation of the fragmentation spectrum in the factorization approach of perturbative QCD, focusing on the case of pure Yang-Mills theory. I generalize to full QCD in the following, as required to compare to data.

Defining the resummed fragmentation function
In the classic case of single-inclusive annihilation, one has a time-like hard momentum Q injected via the off-shell decay of a color-singlet current, and the objective is to know how many of the final state particles each carry the fraction z h of Q. This number density of hadrons carrying fraction z h factorizes as (ref. [22]):

4)
C is the perturbative coefficient function describing the production of a hard parton with momentum fraction z. Subsequently that parton decays into hadrons, including h. This hard parton is coupled to the color-singlet current via fluctuations at the scale Q. The fragmentation function d h describes how to convert this parton into hadrons, and is inherently non-perturbative, sensitive to the confinement scale Λ. Thus in perturbation theory D is infra-red divergent. However, when one calculates D in a regularization scheme, for instance using dimensional regularization, after the renormalization of ultra-violet divergences, one can take d h to be the infra-red divergences of the bare calculation of the number density, and thus acts as a counter-term for these IR divergences. Going to moment space, one has: Now γ T is the time-like DGLAP anomalous dimension. The exponentiation of anomalous dimension γ T reproduces all the IR divergences of the perturbative calculation, and is the MS fragmentation function, though not the physical one. Continuing the interpretation, Λ is the location of the Landau-pole for the renormalized coupling. Finally, µ is the arbitrary factorization scale, and renormalization group evolution in µ resums the logs of Q/Λ: Through this equation I resum the logs of Q/Λ, and γ T has been calculated to three loops in ref. [23], while partial results can be found in refs. [24,25]. In what follows I need to resum C, γ T , and d h , not just for the logs of Q/Λ, but also logs of z h as z h → 0. To accomplish these goals, I turn to a distinct factorization of D in terms of angles. Rather than dimensional regularization, I can regulate the infra-red divergences by cutting off the angle that any other emission may approach the observed hadron, defining the angular scale R ir . However I keep all quantities in 4 − 2 dimensions until after I have evolved in the angular cutoff R ir . This results in a number density D for the hadrons carrying z h that is a function of R ir and Q 2 . The number density satisfies the angularordered evolution equation: P is determined by the space-like anomalous dimension in 4 − 2 dimensions, introduced in ref. [26,27], and ρ is given by the beta function. For their detailed form, see ref. [12]. Solving this equation and evolving R ir to zero resums the logs of z h in both C and γ T of eqs. (2.4) and (2.7). The fact that the number density of hadrons obeys such an evolution equation implies that I can introduce an angular-ordered fragmentation function d containing all emissions above the angular sale R ir up to some angular factorization scale R uv : 4 R uv is the factorization angle. I use the notation R ir and R uv as R ir cuts off the angle for the infra-red splittings, and R uv cuts off the angle for the ultra-violet function. This factorization is independent of R uv , and I find that d again satisfies eq. (2.8). I can interpret the d in this equation as simply the solution to eq. (2.8) evolved between the angular scales R ir and R uv , with a delta function in δ(1 − z h ) initial condition. Put another 4 I call this scheme an angular factorization due to the argument of the splitting function in 4 − 2 dimensions appearing in eq. (2.8). When continuing to 4 − 2 dimensions, I continue the dimension of the transverse momentum of each loop integral. Thus the dimensional factor µ introduced in the regularization necessarily scales with the continued transverse momentum. Since the ratio µ 2 z 2 Q 2 R 2 ir appears in the splitting function of eq. (2.8), then zQRir must be a transverse momentum scale. Thus as zQ is the energy of the parent parton to the splitting function, I interpret Rir as the angle it makes to the collinear sector defined by the direction of the fragmented particle. way, d is now simply the green's function for the eq. (2.8), and D at R uv is the boundary condition. The independence on R uv follows from the semi-group property that any linear differential equation enjoys: I can evolve the initial condition straight to its end-point, or alternatively, I may stop at any intermediate time, take the solution at that time as a new initial condition, and solve again to the endpoint. The end result must be the same for any intermediate time.
The fragmentation function d resulting from the factorization and evolution in terms of angles is now the fundamental object I wish to consider.
The coefficient function C of (2.4) is resummed when R ir → 0 for d after solving eq. (2.8). Then taking µ = Q, Where R T (α s , n) is simply the small-z h resummed coefficient function C z, µ 2 , µ 2 in moment space. I introduce the notation R T rather than justC to emphasize that this is the resummed quantity for the small-z h logarithms. Likewise, I now posit that the fragmentation function d h of eq. (2.4) is resummed by formally taking the limit R uv → ∞: Where in an abuse of notation I set Λ = µ. This is in accord with standard effective field theory techniques regarding the treatment of cutoffs: a cutoff may be introduced to define a sector of the theory, but this cutoff must be sent to either to its IR or UV limit, depending on whether the modes of the function are UV or IR, respectively. 5 Of course, this relation of Eq. (2.12) cannot be literally true, since I are working in perturbation theory and the fragmentation function is non-perturbative. However, if I were to examine rather than fragmentation to hadrons, but rather semi-inclusive jet production, then such a relation should be true for the semi-inclusive jet function, which appears in the same type of the factorization theorem (eq. (2.4)) as the case of fragmentation, see refs. [15,16]. They share the same coefficient function and evolution equations, differing only in that the low scale function describes how a parton forms jets via some jet algorithm, rather than hadrons. When I claim that eq. (2.12) also resums the fragmentation function for hadrons, then I am invoking the LPHD hypothesis. This is to say that hadrons are but jets with a jet radius of the order of the mass of the hadron divided by the hard scale Q.

The fragmentation function to all orders
Going to moment space, it is straightforward to derive the resummed fragmentation function in terms of the resummed coefficient function to all orders. I then recall the expansion of evolution kernel P in eq. (2.8) within dimensional regularization as worked out in ref. [12]: (2.14) Setting µ = Q, then this equation can be iteratively solved by introducing the functions: Then the iterative expansion is obtained to be: Where all sums are from i = 1 to ∞. Critically, since all functions are defined in dimensional regularization, I have the result: That is, all integrals are scaleless and set to zero. Then at any intermediate angle θ, since d also satisfies eq. (2.9), I have the result: But comparing to eqs. (2.10) and (2.12) this implies the chief result: . (2.20) The resummed fragmentation function is but the inverse in mellin space of the coefficient function. It is straightforward to check against refs. [15,16] that this resummed formula indeed gives the correct logarithms for the n f = 0 semi-inclusive gluon jet function to the known orders. Given this result, the coefficient function and the fragmentation function would cancel exactly to a delta function in eq. (2.4), even after evolution, except for the fact that the running coupling is evaluated at different scales. For moderate evolution one would expect still a large cancellation, only seeing large effects when Q 2 Λ 2 .

Full QCD
I discuss briefly the extension to full QCD for the angular-ordered evolution equations. For now I restrict myself to leading log accuracy, and I write the evolution equation for the evolution matrix: Where we are summing over the flavors in the flavor singlet sector of the theory. Rather than the angular-ordered fragmentation function, I now have the angular-ordered fragmentation matrix. Then given an initial condition f i , I can form a solution for the angularordered fragmentation function: At leading log, the splitting matrix is (in a (g, q) basis for the singlet flavor space, and representing matrices in bold face): Since C 2 = C, this implies that I can solve eq. (2.21) with: whered(n, R 2 ir , R 2 uv ) solves eq. (2.14) truncated to leading log with µ = Q. Then the projector for the coefficient function and jet function respectively is given as: Introducing the standard factorization formula for fragmentation in full QCD: then when coupled with the correct initial conditions (f ir and f uv ), after dropping (canceling) the IR or UV divergences, these matrices of eq. (2.26) will resum the coefficient and fragmentation functions respectively: I then introduce the initial and final conditions in the (g, q) basis for the resummed e + e − cross-section relevant for hadron production: For the order to which I are currently working, R T will given by eq. (3.6).Ū n, Q 2 , Λ 2 is the matrix solution of eq. (2.7) when generalized to multiple flavors in moment space, using the small-z h resummed time-like splitting functions. I have written the hard projected state as: One may wonder why the initial state is not: Strictly speaking, if I expand eq. (2.36), I reproduce perturbation theory. If I expand eq. (2.35), I will get a term proportional to δ(1 − z) that does not appear in the explicit calculation of the e + e − coefficient function. The answer lies in the numerical evaluation of the contour integral performing the mellin inversion, in say eq. (2.34). As long as I stay strictly away from z = 1 exactly, I never pick up the delta function that appears if I were to perform the perturbative expansion and then mellin invert. Thus I automatically implement the delta function subtraction that eq. (2.36) explicitly implements.

NLL versus MLLA
Going back to pure Yang-Mills, I am now in a position to write down the resummed spectrum to all orders. Once I derive the resummed anomalous dimension γ T and coefficient function R T , I have: All singular n-dependence has been resummed into the functions γ T and R T , and the factors B uv and B ir are simply the UV and IR boundary conditions to the evolution. In this sense, they cannot be deduced from the evolution eq. (2.8) or the function d in eq. (2.9), but must be obtained from a matching calculation. It formally holds in pure-Yang Mills theory to at least N 3 LL accuracy, counting α s ∼ n 2 , see ref. [12], and to N 2 LL accuracy I have B uv = B ir = 1.
To compare to the MLLA, I content myself to LL accuracy in the coefficient and fragmentation function. This results in an effective NLL accuracy for the resummed crosssection, so that I call eq. (3.1) NLL accurate when R T is truncated to leading (in small z h log counting) accuracy due to the ratio form the resummed coefficient and jet functions take. The ratio form implies a cancellation at fixed order, such that the expanded ratio is down by a power of α s than it would naively be when α s is taken at a common scale. Thus when the anomalous dimension is taken at NLL level for the small-z h logs, and the coefficient and fragmentation functions are taken at LL in the small-z h counting, I call the whole formula NLL. Of course, I could take all resummed quantities to be NLL in the small-z h counting, and this is perhaps the more consistent thing to do, but I leave such investigations to later work.
The MLLA result is written as: Where γ M LLA is the mixed leading log anomalous dimension, given as (refs. [5,7]): However, the DGLAP splitting kernel to NLL is: Clearly, these are in-equivalent. However, the MLLA formula assumes that all logs are explicitly exponentiated via the anomalous dimension, while the NLL cross-section has implicit logs of Q 2 and Λ 2 in the fragmentation and coefficient functions via the running of α s . This suggests that to compare to MLLA, I ought to also force the exponentiation of the fragmentation and coefficient functions, so I write using the leading log resummation of R T : I have used the renormalization group equation for α s to rewrite: Note that I truncated the solution to the RG equation to the first non-trivial order in α s , and made a similar approximation in R T (α s (Λ 2 ), n) once I redefined the running coupling. This is an unreliable approximation to the leading log running when Λ 2 Q 2 , never mind subleading logs. Pressing on none-the-less, I have: Where in the formula, I assume that α s is evaluated at the scale Q 2 . Note that the difference with the MLLA is entirely due to running coupling effects being truncated in an unjustified fashion. Note that since the difference is independent of n to the order I truncated, the difference will be just an adjustment to the scale evolution of the normalization of the spectrum, assuming I can truncate the running coupling as I do. That is, the collinear log ln Q 2 Λ 2 and its small-z h resummed coefficient, is handled correctly only for moderate evolution.

Parton-Hadron Duality
In this section I extend the small-z h resummed fragmentation function into a genuine model for charged hadron production. First I illustrate the basic features of the model, built on the assumptions that the coupling freezes in the IR, and that the shape of the fragmentation function is given by the small z h -spectrum from resummed perturbation theory. Moreover, I assume that the fragmentation spectrum is dominated by its behavior at small z h , so that naively forcing momentum conservation (a power correction) coupled with iso-spin symmetry dictates the evolution of the normalization. This enforcing of the normalization could be dispensed with if I matched the small-z h resummation to the full DGLAP equations, which would then have momentum conservation built in from the beginning, but I leave that for later work.

Fragmentation model from lnz h resummation
Here I propose that according to the LPHD hypothesis, the shape of the fragmentation function is simply the resummed jet function with an effective α s coupling in the infra-red: 1 + anomalous dimension terms . I have omitted the terms generated by the anomalous dimension that vanish when µ 2 = Λ 2 . For my purposes, I adopt the simple effective coupling by shifting the argument of α s by a constant mass scale, thus freezing the running coupling at a quasi-perturbative scale: I have also checked that using a hard cutoff in the running coupling, stopping it at a fixed mass scale, also works as a non-perturbative model, yielding decent fits to the data. In total, this implies I have 3 degrees of freedom in the model for the fragmentation function: the scale I freeze the coupling given by c, the normalization of the fragmentation function N h ± , 6 and finally the scale I stop DGLAP evolution, the scale of confinement Λ 2 . While this effective coupling is strictly a guess that at the form of the non-perturbative dynamics, in the spirit of LPHD, the form of the function is completely determined by the small-z h resummation, and was calculated in perturbation theory. There will be limitations to such a fragmentation function, particularly in the extremely small-z h limit, z h ∼ m 2 h Q 2 , when the mass of the hadron m h causes significant changes to the dispersion relation between momentum and energy. Further it is expected to be wrong in the threshold region z h → 1.
In pure Yang-Mills, I write the prediction for the fragmentation spectrum as: Note that for consistency sake I shift the argument of the coupling everywhere, but in practice this does not matter for the high-scale function or anomalous dimension. N (Q 2 , Λ 2 ) is a function dictated by the demand the fragmentation spectrum should obey the sum-rule: That is: As far as the expansion in z h is concerned, this is a perfectly valid modification to the spectrum, since formally these are power corrections that go beyond the small-z h approximation, and does not therefore invalidate the logarithmic counting or changes the shape of the spectrum. This is a phenomenological motivated modification that is based on the following assumptions: • That I will predict the shape of the spectrum correctly for all Q.
• Momentum is conserved, and the shape of the spectrum of neutral hadrons is identical to charged hadrons. In full QCD, this assumption hinges on the extended iso-spin symmetry governing up, down and strange quarks, since light hadron production is dominated by nucleons, pions, and kaons, and thus can be viewed as combining the LPHD hypothesis with iso-spin symmetry.
• The spectrum is suppressed at z h → 1, and the small z h region starts at z h ∼ 0.5, so the fact that the fragmentation function is wrong above this point is of little consequence, as long as it is small.

Implementation
I evolve the DGLAP equations in momentum space, using a Runge-Kutta based discretization scheme detailed in App. B. I thus avoid any diagonalization of matrices as must be performed for mellin space evolution. The appropriate spectrum for e + e − is calculated as: The function F is given by eq. (2.34), and U ij is calculated by starting the Runge-Kutta evolution (in the notation of App. B) with the initial conditions H g (z, t = 0) = δ(1 − z), H q (z, t = 0) = 0, and evolving with the full NLL small-z h resummed singlet splitting matrix given in ref. [11]. At the b-quark threshold I switch to a 4 flavor scheme for the splitting functions and α s , though remaining in a 5 flavor scheme makes little difference in the fits. To compute F , I deform the integration contour into a staple as follows: first I where is a small number. Then ±in 0 is approximately where the cut furthest from the origin begins in the integrand. The cuts corresponding to n ∼ √ α s at hard scale are closer to the origin and bounded by the staple automatically.
I deform the mellin inversion contour into a staple shape surrounding the cuts, stretching from −∞ − in 0 , to − in 0 , then up to + in 0 , and from there back to −∞ + in 0 . All integrals on each leg of the staple are then nicely numerically convergent. I then numerically compute the convolution, and derive N (Q 2 , Λ 2 ) by enforcing the sum rule in Eq. (4.4). Finally, I use the two-loop running of α s , with α s = 0.1187 at the Z-pole.
To gain some intuition about the predicted shape of the fragmentation function I derived in this paper, in Fig. 1 I plot the numerical form of the fragmentation function derived in Sec. 2.2 where the resummation is performed to leading log accuracy for small energy fraction logs. Thus R T is defined with eqs. (3.6) and (3.3). I also present the
The fragmentation spectrum for e + e − → h ± + X at center of mass energies 58 and 91 GeV. The gray region is excluded from the χ 2 /d.o.f. determination based on hadron mass corrections or z h ∼ O(1), and the quoted χ 2 /d.o.f. is the goodness-of-fit for that data set alone using the best fit parameters. The fragmentation spectrum for e + e − → h ± + X at center of mass energies 189 and 206 GeV. The gray region is excluded from the χ 2 /d.o.f. determination based on hadron mass corrections or z h ∼ O(1), and the quoted χ 2 /d.o.f. is the goodness-of-fit for that data set alone using the global best fit parameters.
The upper end of these regions is determined by when I expect hadron mass corrections to become significant. To estimate this, I compare the energy fraction (z h ) versus the 3momentum fraction of the hadron (x p ), whose relation is given by: with m h the hadron mass (ref. [35]). Since I am fitting for generic light hadron production, I take m h ∼ 0.5 GeV, the mass of the kaon, and demand that |x p −z h | < 0.2. More stringent criteria tends to improve the χ 2 /d.o.f. for the fit. For comparison to the limiting spectrum of the MLLA approximation, I also compare to OPAL data at Q = 91GeV, which was fitted to the limiting spectrum, and which is largely consistent with the SLD data, though with larger error. This data I do not include in any of the χ 2 determinations. The total χ 2 per degree of freedom (χ 2 /d.o.f. ) is 4.6, however, the fit is completely dominated by the SLD data due to the extremely small error bars quoted. Restricted to SLD, ALEPH, and L3 data alone the χ 2 /d.o.f. is 1.1. Finally, I note that the SLD data is consistent with the earlier OPAL data of ref. [8], and ALEPH data at Q = 206 GeV is consistent with the L3 data at the same energy.
Though the above quoted regions contain the points I use for the χ 2 determination, I note that qualitatively the fits extrapolate essentially for all z h , and this is even true quantitatively for the data at the Z-pole and above. For instance, if I compute the χ 2 /d.o.f. on the L3 data at 206 GeV alone, I obtain 0.3 for all data points between 0.5 < −lnz h < 6.6 and similarly for the ALEPH data where I would obtain a χ 2 /d.o.f. of 2.5 in the same range, without changing the fit.
I find that the best fit values are obtained with N h ± = 1.188, c = 1.0, Λ = 1.0 GeV. I note that the TASSO data at lower Q 2 data prefers a larger value of Λ, at 2 GeV, but the same normalization. I speculate as to reasons for this in sec. 5.4. I note that I can obtain good fits to the data with c = 2.0 or c = 0.5 if I also change the end point of the splitting evolution, and one could also explore attempting to fit the value of α s as well.

Comparison to MLLA
I can compare to the two common approximations for the MLLA, the so-called limiting spectrum and the Gaussian saddle-point. For completeness sake, I give the formula for the limiting spectrum and the Gaussian spectrum in App. A, as well as the fit values. The spectra are supposed to be valid for Q → ∞ and describe the peak region of the distribution, and good fits were obtained to Z-pole data in ref. [8] using 3 active quark flavors. In Fig. 5, I present the original OPAL data, assuming a 1% systematic error along with their quoted statistical uncertainty. I plot their MLLA limiting spectrum fit, as well as the Gaussian 5 active flavor fit, and quote χ 2 /d.o.f. for the fit region 1 < −lnz h < 4.1. Note that this fit region is smaller than the one used for the SLD data above, due to the Gaussian's more limited range of applicability. I see that all schemes appear capable of describing the peak of the OPAL data, 7 however, when evolved to the L3 scale, the predicted MLLA spectra fail to describe the data, except at moderate to small values of −lnz h , while the NLL is uniformly valid throughout the fit region, and even beyond. I note also that the NLL result uses a physical number of active quark flavors at LEP energies.

Low Q Data
The global fit described at the beginning of this section is dominated by the very precise SLD data at the Z-pole, which evolves nicely to describe the higher Q LEP data. However it is at tension with the lower Q TASSO data. I can refit to the lower Q data, and in Fig.   7 I also refit the NLL spectrum to the OPAL data, this results in a slightly larger normalization.  1), and the quoted χ 2 /d.o.f. is the goodness-of-fit for that data set. Compared to the NLL curve is the evolved 3 flavor scheme MLLA limiting spectrum (purple LS) from ref. [8], and a 5 flavor scheme evolving the moments of a Gaussian according to the MLLA anomalous dimension, as in ref. [5] (orange G).  Figure 6. The fragmentation spectrum for e + e − → h ± + X at center of mass energies Q = 35 and 44 GeV. The gray region is excluded from the χ 2 /d.o.f. determination based on hadron mass corrections or z h ∼ O(1), and the quoted χ 2 /d.o.f. is the goodness-of-fit for that data set. I compare the global fit prediction to a fit to the lower Q data alone. The improved fit is achieved by only changing the scale at which splitting evolution is ended. 6 I plot the new fit to the lower Q data of TASSO, while plotting for comparison the global fit curve. The new fit differs only in that it prefers a higher value of Λ, at 2.0 GeV. The normalization and coupling freezing scale are the same. That is, since Λ TASSO = 2Λ Global , c TASSO = 1 4 c Global , so that the mass scale used to freeze the coupling remains the same, c TASSO Λ 2 TASSO = c Global Λ 2 Global . I give the χ 2 /d.o.f. for both curves as a metric to see how the fit improves. This also nicely illustrates the impact of varying the end-point of just the evolution, as the other fit parameters are the same. Finally, I note I can drop the χ 2 /d.o.f. to 5.4 for the TOPAZ data with a Λ = 1 GeV like the global fit, but taking the normalization to 1.1.
At present, I can only speculate as to why ending the DGLAP evolution at a higher energy improves the fit. Clearly, there could be an issue of the missing higher order terms at NNLL or beyond in the resummation, which are becoming more important at lower vales of the hard initial scale. This should be investigated. However, the global fit and the specific TASSO fit are largely consistent at all z h values leading up to the peak, differing mainly in the peak itself and the region I omit due to potentially large hadron mass correction. Interestingly, the new fit interpolates into this region as well. This suggests that as hadron mass corrections become critical, they can be accounted for by tailoring the end-point of the DGLAP evolution to each hadron mass. Thus if I were to fit each charged species separately, varying the end point of evolution yet keeping the same basic fragmentation model proposed in this paper, I could by and large account for hadron mass effects. Then as Q → ∞, the hadron mass effects become suppressed, and the endpoint of the DGLAP evolution tends to a universal value at the confinement scale. However, such investigations are beyond the scope of this paper.

Conclusions
The concept of local parton-hadron-duality gives important insight into the process of hadronization, and helps explain the success of QCD perturbation theory, beyond the traditional estimates of QCD factorization. It is important to note that using the traditional power counting regions for the modes informing the factorization for fragmentation, local parton-hadron-duality could be considered surprising, see refs. [36,37]. And indeed, treating the fragmentation spectrum with the standard momentum regions, tacitly or explicitly assuming z h ∼ O(1), leads to using the fixed-order expansion for all quantities. This does not achieve the fragmentation spectrum given here. A soft resummation is needed.
However, the original formulation of LPHD rested on unsystematic approximations to the QCD perturbation theory, making it unclear whether the success of the LPHD hypothesis rested on these unsystematic approximations, and whether a more grounded approach to QCD perturbation theory would destroy the LPHD concept. I find that this is not the case: systematically improving the small-z h resummation to NLL accuracy, in the coefficient function, anomalous dimension, and jet/fragmentation function gives a uniformly valid description of the fragmentation data, by extrapolating semi-inclusive jet production to hadrons and freezing the infra-red coupling at a quasi-perturbative scale. Even though the small-z h resummed ansatz worked well, there is still room for improvement. First one can extend the small-z h resummation to NNLL and even N 3 LL order. Further, one can match to the threshold region at z h = 1 (see refs. [38,39]), hopefully enforcing the sum-rule for the spectrum "naturally." Moreover, one would also like to use the small-z h resummation to improve the description of jet production rather than just hadron fragmentation, constructing a jet function that is uniformly valid in all regions of the energy fraction.
From a more formal point of view, the results are very pleasing: the fragmentation spectrum of hadrons is truly determined by reciprocity with the space-like branching pro-cess, when expressed in angular ordered evolution equations. Initially this reciprocity was postulated for the anomalous dimensions alone in refs. [40,41], and found a critical role in the scaling properties of energy-correlation functions (refs. [42,43]). But to capture the fragmentation spectrum, this reciprocity must be extended to the full cross-section as done here.

Acknowledgments
D.N. wishes to thank Felix Ringer and Nobou Sato for discussion of fragmentation and Wouter Waalewijn for reading the manuscript. D.N. was supported by the Department of Energy under Contract DE-AC52-06NA25396 at LANL and through the LANL/LDRD Program.
be handled more carefully. In particular, I use the formula: