Hadron mass corrections in semi-inclusive deep-inelastic scattering

The spin-dependent cross sections for semi-inclusive lepton-nucleon scattering are derived in the framework of collinear factorization, including the effects of masses of the target and produced hadron at finite momentum transfer squared Q^2. At leading order the cross sections factorize into products of parton distribution and fragmentation functions evaluated in terms of new, mass-dependent scaling variables. The size of the hadron mass corrections is estimated at kinematics relevant for future semi-inclusive deep-inelastic scattering experiments.


I. INTRODUCTION
Determining the detailed flavor and spin structure of the nucleon remains a central challenge for hadronic physics into the 21st century. Considerable progress has been made over the past two decades in understanding the characteristics of the momentum and spin distributions of quarks and gluons (or partons) through precise measurements of the nucleon's parton distribution functions (PDFs) in various hard scattering reactions [1,2]. In addition to the traditional inclusive deep-inelastic scattering (DIS), Drell-Yan and other high-energy scattering processes, an increasingly important role in this quest has been played by semiinclusive production of hadrons in lepton-nucleon scattering.
Identification of specific hadrons, such as pions or kaons, in the current fragmentation region of a deep-inelastic collision serves as a tag of individual quark flavors, which in inclusive DIS are summed over. Important insights have been provided through semi-inclusive deepinelastic scattering (SIDIS) experiments on phenomena such as the SU(2) flavor asymmetry in the proton sea [3] and the ratio of strange to nonstrange quark distributions [4]. From experiments with polarized targets, SIDIS data have also provided fascinating glimpses of the possible flavor asymmetry in the polarized light-antiquark sea [5], while kaon production data has fueled the recent controversy concerning the sign of the polarized strange sea [6,7]. Furthermore, detection of forward baryons (in the center of mass frame) in the target fragmentation region of SIDIS is a potentially important avenue for extracting information on the pion cloud of the nucleon or the structure of the virtual pion itself [8][9][10][11].
In more recent developments, detection of both the longitudinal and transverse momentum distributions of hadrons produced in SIDIS measurements of various single-and doublespin asymmetries has opened up the largely unexplored realm of transverse momentum dependent parton distributions [12][13][14]. These reveal an even richer landscape of threedimensional momentum and spin distribution of partons in the nucleon, that will be the subject of increasingly greater attention at facilities such as Jefferson Lab [15] and COM-PASS [16], and a central component of the science program at the proposed Electron-Ion Collider [17].
The unambiguous interpretation of any SIDIS experiment in terms of leading twist PDFs or transverse momentum distributions requires control of various subleading 1/Q 2 corrections, such as target mass and higher twist effects, as well as knowledge of the fragmentation functions describing the parton hadronization. For inclusive DIS, the finite-Q 2 corrections are known to become important at low Q 2 values, particularly when the parton momentum fraction x is large [18][19][20]. Their effects on global fits of spin-averaged PDFs have been systematically studied in recent analyses by the CTEQ-Jefferson Lab (CJ) [21][22][23], ABM [24,25] and JR groups [26], and in spin-dependent PDF analyses by the JAM collaboration [27] (and to some extent also by the LSS [28], BB [29] and NNPDF [30] groups).
Typically, the effects of target mass corrections (TMCs) can be computed within a specific framework, while higher twist effects, which involve more complicated multi-parton correlations, are parametrized phenomenologically. The standard approach for computing TMCs has traditionally been within the operator product expansion, in which the mass corrections to inclusive DIS structure functions arise from twist-two quark bilinear operators with an arbitrary number of derivative insertions [18,31]. Extending this framework to processes involving particles in the final state is problematic, however, which has in practice limited the study of hadron mass corrections in SIDIS.
An alternative framework for TMCs was developed using techniques based on collinear factorization (CF) [32][33][34], in which the hard scattering is formulated in momentum space directly. The method has been applied to the computation of TMCs in inclusive scattering, both in unpolarized [35][36][37] and polarized [38] DIS, and in semi-inclusive hadron production in electron-proton anhilation [39,40]. For semi-inclusive hadron production in lepton-proton collisions, in contrast to inclusive DIS and e + e − anihilation, finite-Q 2 corrections can arise from both the effects of the target mass and the mass of the produced hadron. While earlier analyses [41,42] considered some of these corrections within the CF framework, the phenomenology of the combined effects of the target and produced hadron masses -which we refer collectively as "hadron mass corrections" (HMCs) -was systematically explored in Ref. [43] for unpolarized scattering.
In this work we extend the analysis of HMCs to the case of spin-dependent SIDIS at finite Q 2 . Because high energy spin-dependent data are generally more scarce than spinaveraged cross sections, a significantly larger fraction of the world's data set used to constrain spin-dependent PDFs lies in the low-Q 2 region (Q 2 ∼ 1 − 2 GeV 2 ). While target mass corrections have been incorporated in some global spin-PDF analyses [27][28][29][30], none of the analyses which have included polarized SIDIS data [28,44] have accounted for HMCs. With the increasing precision of new polarized measurements, and the consequently more accurate determination of spin-dependent PDFs, it is imperative to reliably account for subleading corrections which could impact the extraction of the leading twist distributions.
In Sec. II we outline the formalism used to compute the SIDIS cross sections at finite values of Q 2 within the collinear approximation in the presence of target and produced hadron masses. For completeness, we consider both polarized and unpolarized scattering, since the latter enters the calculation of the measured polarization asymmetries. In Sec. II B we review the collinear formalism and its application to hadron production in SIDIS. Expanding the hadronic tensor in terms of quark correlators, in Sec. II C we derive semi-inclusive cross sections, which at leading order are given by factorized products of PDFs and fragmentation functions expressed as functions of modified scaling variables. The relative importance of the HMCs is explored numerically in Sec. III, where we quantify the dependence of the finite-Q 2 cross sections on the kinematical variables and estimate the corrections for specific current and future experiments. Finally, we conclude by summarizing our results in Sec. IV.

II. SEMI-INCLUSIVE SCATTERING WITH MASS CORRECTIONS
The semi-inclusive lepton-nucleon scattering process is illustrated in Fig. 1, where the incident lepton (with momentum ℓ) scatters from an initial state nucleon (p) to a recoil lepton (ℓ ′ ) via the exchange of a virtual photon (q), producing a final state hadron h (with momentum p h ). In this section we first review the external kinematics and choice of variables, before outlining the collinear factorization framework for describing the hard scattering process. After defining the hadronic tensor in terms of quark-nucleon and quark-hadron correlation functions, we derive expressions for the spin-averaged and spin-dependent cross sections in terms of parton distribution and fragmentation functions at leading order in the strong coupling constant and for finite values of Q 2 .

A. External kinematics
We expand the four-momenta of the external particles in terms of coplanar light-cone unit vectors n and n, satisfying n 2 = n 2 = 0 and n · n = 1 [32]. The "plus" and "minus" We work in a class of reference frames where the initial nucleon and virtual photon momenta are coplanar, so that where M is the nucleon mass, Q 2 = −q 2 , and the scaling variable ξ = −q + /p + . The (p, q) collinear frames are related to each other by a boost of parameter p + and contain, in particular, the target rest frame (where the nucleon p + component is p + = M/ √ 2) and the Breit frame (p + = Q/( √ 2ξ)). For other choices of reference frames see Appendix A of Ref. [43]. In the collinear (p, q) frame, ξ is identical to the Nachtmann scaling variable [31,45,46], where x B = Q 2 /2p · q is the Bjorken scaling variable, and in the Bjorken limit (Q 2 and q − → ∞ with x B fixed) one finds that ξ → x B . Momentum and baryon number conservation (and, for K production, strangeness conservation) impose an upper limit on with M b = M for π production and M b = M Λ (the lightest hyperon mass) for K production.
This value of x B corresponds to a final state consisting of the nucleon and hadron h at rest in the hadron's rest frame.
For spin-dependent scattering, the polarization vector S µ of the initial state nucleon can be parametrized as and satisfies the conditions p · S = 0 and S 2 = −1. In the case of a longitudinally polarized initial state nucleon (which we consider in this work), one has S L = ±1 and the transverse spin vector S µ T = 0. The incident and scattered lepton momenta can be decomposed as where ℓ µ ⊥ is the lepton transverse momentum four-vector, η = ℓ + /p + is the lepton momentum fraction, and we assume massless leptons, ℓ 2 = ℓ ′2 = 0. After some algebra one can show where y = q · p/ℓ · p and γ 2 = 4x 2 B M 2 /Q 2 . In the target rest frame y = ν/E is the fractional energy transfer from the lepton to the target, with ν and E the virtual photon and incident lepton energies, respectively. The magnitude of the lepton transverse momentum is set by four-momentum conservation, For the hadron produced in the final state, the momentum is parametrized as where ζ h = p − h /q − is the scaling fragmentation variable, and the hadron transverse momentum four-vector p µ h⊥ satisfies p h⊥ · n = p h⊥ · n = 0, with norm p 2 h⊥ = −p 2 h⊥ . The squared transverse mass of the produced hadron h is defined by m 2 h⊥ = m 2 h + p 2 h⊥ , where m h is the mass of the hadron. The variable ζ h can be related to the fragmentation invariant In the target rest frame it coincides with the ratio of the produced hadron to virtual photon energies, z h = E h /ν, which is frequently used in experimental analysis of SIDIS data. In the Bjorken limit the fragmentation variable ζ h → z h , while at finite Q 2 one has Since the produced hadron's energy is bounded from below by E h ≥ m h⊥ , one can show Combining Eqs. (12) and (13) with again M b = M for h = π and M b = M Λ for h = K. This limit corresponds to diffractive production of the observed hadron with maximal energy. As Q 2 → ∞, both the upper and lower limits become independent of Q 2 , z min h → 0 and z max h → 1.
As an alternative to the fragmentation invariant z h , one can define the invariant momentum fraction which is used in the study of hadron production in e + e − collisions. The choice of this variable avoids mixing inclusive (x B ) and semi-inclusive (z h ) variables as in Eq. (11). It also allows a clean separation of the current (z e > 0) and target (z e < 0) fragmentation regions in the Breit frame, in which z e = p z h /q z is the ratio of the longitudinal components of the hadron and photon momenta. In the current region, where the observed hadrons are produced with longitudinal momentum in the direction of the virtual photon, ζ h can be written in terms of the z e variable as Hadrons produced in the current region have ζ h > ζ h , which reflects the fact that a hadron produced at rest in the target rest frame belongs to the target region. Finally, in the Bjorken limit all three fragmentation variables become equivalent, ζ h → z h → z e , and the current region extends down to the The relation between the variable ζ h and z h is illustrated in Fig. 2 for several fixed values of x B and Q 2 . At finite Q 2 the kinematically allowed regions of z h are determined by Eqs. (13) and (14), and the boundaries between the current and target fragmentation regions occur h . For the production of pions, at low x B 0.3 and high Q 2 5 GeV 2 the differences between the two variables are almost negligible, and begin to be noticeable only for the highest z h values at lower Q 2 or higher x B . At these kinematics, pions are produced in the current region for z h 0.08 − 0.15. At high x B ≈ 0.6 and low Q 2 ≈ 1 GeV 2 , the effects are more pronounced, with deviations of ∼ 30% as z h → 1, and current fragmentation begins at a slightly higher z h . The effects of the kinematic lower limit in z h [Eq. (13)] is noticeable only at low Q 2 and high x B .
For kaons, the effects at high Q 2 and low x B are again negligible, although the larger K mass enhances the differences relative to the pion at the same kinematics. In particular, at low Q 2 = 1 GeV 2 the lower limit on ζ h is dramatically increased, and the current region is pushed to higher values of z h . At sufficiently large x B , the phase space for K production eventually vanishes; for Q 2 = 1 GeV 2 , for example, no kaons can be produced in the current fragmentation region with x B 0.66, and no kaons can be produced at all for x B 0.8.

B. Collinear factorization
At leading order in the strong coupling α s , the SIDIS reaction proceeds through the hard scattering of the virtual photon from an initial state quark with momentum k to a quark with momentum k ′ = k + q, which then fragments to a hadron h, as illustrated in Fig. 3. Higher order processes involving gluon radiation and scattering from qq pairs can be considered, but for clarity of the derivation of the finite-Q 2 corrections we restrict ourselves to the leading order calculation.
The parton four-momenta can be parametrized, in analogy with the external variables in Sec. II A, in terms of the light-cone vectors n and n as where x = k + /p + is the light-cone momentum fraction of the nucleon carried by the struck quark, and z = p − h /k ′− is the light-cone fraction of the fragmenting quark carried by the hadron h. The parton transverse momentum four-vectors k ⊥ and k ′ ⊥ are orthogonal to n and n, k ⊥ · n = k ⊥ ·n = 0 and k ′ ⊥ · n = k ′ ⊥ ·n = 0, with norms k 2 The essence of the collinear factorization approach is to expand the hard scattering term around "collinear" momenta k and k ′ , where the initial and final parton virtualities k 2 and k ′2 are kept for generality. In this approximation the transverse momentum of the initial parton is neglected, and the transverse momentum of the final parton is taken along the direction of the hadron h, making these collinear to the proton target and produced hadron, respectively.
Defining the partonic analog of the Bjorken scaling variable x B at finite Q 2 bŷ one can show thatx lies within the rangex min ≤x ≤x max , where [37,43] Here the lower limitx min corresponds to the minimal mass of collinear spectator partons, while the upper limitx max arises from the minimum value of the current jet invariant mass.
These limits are consistent with the limit on x B in Eq. (4) for any k 2 ≥ x(ζ h −1)Q 2 /ξ, and in the Bjorken limit the range ofx reduces to being between x B and 1, as required. Analogous limits can be found for the partonic fragmentation variable z, z min ≤ z ≤ 1, where For the practical implementation of the finite-Q 2 kinematical limits, a choice of the initial and scattered parton virtualities k 2 and k ′2 must be made. For the initial parton the requirement that the collinear parton virtuality is independent of x leads to the restriction k 2 ≥ 0. For light, bound initial-state partons ( k 2 ≤ 0) this constrains therefore the collinear expansion to be around k 2 = 0. Determining the virtuality k ′2 of the scattered parton is generally less clear, on the other hand. For the leading order hard scattering in Fig. 3, conservation of four-momentum and the condition k 2 = 0 constrain the parton momentum fraction x to equal ξ h , where In order for x to respect the limits in Eq. (23), the choice k ′2 = m 2 h /ζ h was proposed in Ref. [43], giving . Larger k ′2 values would also allowx to fall within the bounds in Eq. (23); however, the above choice is the closest to the physical quark mass, and is the one adopted in our numerical analysis here. The dependence of the calculated cross sections on the choice of k ′2 is explored further in Sec. III.

C. SIDIS at finite Q 2
The cross section for the SIDIS process can be written as a product of hadronic (W µν ) and leptonic (L µν ) tensors [13], where α is the electromagnetic fine structure constant. The leptonic tensor can be computed straightforwardly from QED, where λ is the lepton helicity.
The semi-inclusive hadronic tensor is defined in terms of matrix elements of the electromagnetic current operator J µ between the initial state nucleon with spin S and the final state with a hadron h and unobserved hadrons X, where we use the shorthand notation d 3 p X /2E X = i∈X d 3 p i /2E i , with p X the total momentum of teh unobserved hadrons. The hadron tensor can be expressed at leading order in α s in terms of quark-quark correlators Φ q and ∆ h q , associated with the quark distribution and fragmentation functions, respectively [33,34,42], where the sum is taken over quark and antiquark flavors q. Note that since the parton momenta in the δ-function have been approximated by their collinear components, the integrations over dk − d 2 k ⊥ and dk ′+ d 2 k ′ ⊥ act directly on the correlators Φ q and ∆ h q . The correlator associated with the parton distribution function is defined, in the lightcone gauge, as [47] Φ where ψ q is the quark field operator, and E X ′ and p X ′ are the energy and momenta of the intermediate state X ′ corresponding to a nucleon with a quark removed. Similarly, for the quark fragmentation correlator one has [47] ∆ where E Y ′ and p Y ′ are the energy and momenta of the state Y ′ resulting from the fragmenting quark with the hadron h removed.
The leading twist contributions to the SIDIS cross section can be extracted by expanding the integrated correlators and selecting the terms linear in the light-cone vectors. For the parton distribution function, one finds where the first and second terms define the spin-averaged, q(x), and spin-dependent, ∆q(x), distribution functions, and the ellipsis indicates higher twist contributions [13]. The fragmentation function D h q (z) is analogously defined from the ∆ h q correlator by Inverting Eqs. (32)- (33), one can write the quark distribution and fragmentation functions explicitly by projecting with the appropriate Dirac operators, where ω ± are light-cone coordinates. (For ease of notation we omit the Q 2 dependence in the arguments of the quark distribution and fragmentation functions.) The fragmentation function D h q (z) here is defined with the standard normalization, h 1 0 dz z D h q (z) = 1 [42]. To compute the hadronic tensor in Eq. (29), we can decompose the δ (4) ( k + q − k ′ ) function along the +, − and transverse components of the momenta. The δ-functions for the + and − components constrain the partonic variables to x = ξ h and z = ζ h , respectively. The δ-function for the transverse component forces the transverse momentum of the produced hadron h to vanish, p h⊥ = z k ′ ⊥ = 0. Nonzero transverse momentum hadrons can be produced via higher order perturbative QCD processes, or from intrinsic transverse momentum in the parton distribution functions themselves [13]. With these constraints, the hadron tensor can then be factorized into products of parton distribution and fragmentation functions evaluated at ξ h and ζ h , respectively, The main effect of the hadron masses at finite kinematics is therefore a replacement of where the reduced unpolarized and polarized cross sections are defined as In Eqs. (38) the arrows denote the spins of the lepton and target nucleon, and the dependence of the functions on the scale Q 2 is made explicit. In Eqs. (38) the variable is the ratio of longitudinal to transverse photon flux, and in Eqs.
Note that at the maximum value of x B allowed for SIDIS [see Eq. (4)] the finite-Q 2 variable ξ h satisfies ξ h < ξ h (x B = x max B ) < 1. As in the case of inclusive DIS [37], the SIDIS cross section therefore does not vanish as x B → x max B , which reflects the well-known threshold problem in which the leading twist structure function is nonzero for x B ≥ 1 [48][49][50][51]. Analogously, for the finite-Q 2 fragmentation varialble one has ζ h < ζ h (z h = z max h ) < 1, and since the fragmentation function does not vanish as z h → z max h , the perturbatively calculated SIDIS cross section can also exceed the fragmentation threshold.
Before exploring the dependence of the SIDIS cross sections on the finite-Q 2 scaling variables in the next section, it is useful to first establish the purely kinematic corrections, independent of the PDFs and fragmentation functions, which augment the finite-Q 2 results from their scaling limit. In Fig. 4  ≈ 25% below unity at large z h . This behavior should be kept in mind when assessing the numerical effects of the HMCs in the cross sections in the next section. One should also note that the upturn in J h is almost entirely confined to the target fragmentation region (small z h ), where the validity of calculations based on the perturbative handbag diagram in Fig. 3 is more questionable, and factorization in terms of fracture functions [52][53][54] may be more appropriate.

III. PHENOMENOLOGICAL IMPLICATIONS
In this section we examine numerically the phenomenological consequences of the finite-Q 2 rescaling of the SIDIS cross section derived in Sec. II, and explore their impact on future hadron production experiments.

A. Kinematical dependence of HMCs
To disentangle the separate HMC effects in the SIDIS cross sections arising from the PDFs and fragmentation functions, in Fig. 5 we illustrate the ratios of PDFs at finite Q 2 to the corresponding massless limit distributions. For a systematic comparison we consider The most dramatic feature in the ratios is the steep rise at large x B , which sets in at smaller x B values for decreasing Q 2 . The results are qualitatively similar for the unpolarized and polarized distributions, with the rise delayed to slightly larger x B for the latter. The differences between the unpolarized and polarized PDF ratios for the most part reflect the differences in the shapes of the respective input PDFs. However, the qualitative features of the results do not change when using other leading order distributions, such as the unpolarized GJR [57] and polarized BB [58] parametrizations. Generally, the behavior of the ratios observed in Fig. 5 is reminiscent of that found in previous studies of TMCs for inclusive DIS [19].
For the fragmentation functions, ratios of the finite-Q 2 isoscalar functions for π + + π − production to those in the massless limit D (0) are displayed in Fig. 6 for fixed values of x B and Q 2 , using the leading order HKNS parametrization [59] for the fragmentation functions. for z h 0.7, before rising rapidly as z h → 1. The effect is less pronounced with increasing Q 2 , with a smaller enhancement of the ratio and a delayed (though even more dramatic) rise at large z h . The fragmentation function ratios at x B = 0.8 in Fig. 6(b), on the other hand, display a significantly stronger enhancement, particularly at the lowest Q 2 value. Here the effect is about an order of magnitude larger, and features a striking upturn at z h 0.2, where the finite-Q 2 fragmentation function becomes several times larger than the high-Q 2 limit. As outlined in Ref. [43], this arises from the shape of the fragmentation function at finite-Q 2 kinematics. The general features of the results, however, remain unchanged if one uses the KKP parametrization [60], for instance.
In particular, expanding the hadron mass corrected fragmentation function in a Taylor series about the massless limit, one can write the corrected to uncorrected ratio as The z h dependence of the correction is thus determined by the negative shift in the fragmentation variable (ζ h − z h ) and by the z h slope of D(z h ). Since the pion fragmentation function is generally a decreasing function of z h at small z h , the ratio is driven upward as z h → z min h , where |ζ h − z h | is maximum. Note that for kaons and protons, in contrast, the slope at small z h can be positive, which would lead to a suppression of the mass corrected function at z h ∼ z min h . At the other extreme, in the exclusive production limit the fragmentation The open circles denote the boundaries between the target and current fragmentation regions.
function ratio becomes divergent for the same reason as the PDFs; namely, the functions in the scaling limit vanish as x B → 1 or z h → 1, whereas the finite-Q 2 scaling variable and the corresponding rescaled functions remain finite.
Combining the effects of the HMCs in the parton distribution and fragmentation functions, in Fig. 7 we show ratios of the SIDIS spin-averaged and spin-dependent cross sections with and without HMCs as a function of z h , for several fixed values of Q 2 and x B . Specifically, we consider scattering from a proton target, with the production of π + + π − mesons in the final state. The massless limit cross sections σ Overall, the z h dependence of the cross section ratios follows the trends indicated in Figs. 5 and 6 for the PDF and fragmentation function ratios, with strong enhancement of the finite-Q 2 cross sections at large z h , and decreasing effects at higher Q 2 . At the lower x B value (x B = 0.3 in Fig. 7(a)), the HMC effects in the spin-dependent (thick lines) and spin- polarized ratios are almost indistinguishable. The differences are more striking at larger x B (x B = 0.8 in Fig. 7(b)), where the effects on the spin-averaged cross section are somewhat larger than on the spin-dependent cross section. This stems directly from the delayed rise above unity of the ∆q/∆q (0) PDF ratio in Fig. 5(b) at high x B values compared with the corresponding q/q (0) ratio in Fig. 5(a).
While most of the existing SIDIS data have involved the production of charged pions, the detection of heavier mesons and baryons can provide complementary information on the flavor and spin structure of PDFs, as well as on the dynamics of hadronization. The production of kaons, for instance, tags strange or antistrange quarks, and has been used with polarized targets as an independent means of determining the ∆s distribution in the nucleon [6], and for unpolarized scattering to determine the magnitude of the s distribution at small x B [4]. In Fig. 8 the spin-averaged and spin-dependent cross section ratios with and without HMCs are shown for the production of charged pions (π + + π − ) and kaons (K + + K − ) at x B = 0.3 and Q 2 = 5 GeV 2 . The effects are enhanced significantly with increasing hadron mass, particularly at low z h values, mostly because of the (1 + m 2 h /ζ h Q 2 ) factor in the ξ h variable in Eq. (25). Increasing values of m 2 h /ζ h shift the argument of the PDF to higher x B , where the smaller magnitude of the distributions effectively suppresses the cross section. This phenomenon inherently arises from the choice of invariant mass squared k ′2 for the scattered quark, discussed in Sec. II B. While the choice of the mass becomes irrelevant at high Q 2 , the dependence on k ′2 can be appreciable at low Q 2 values, as Fig. 9 illustrates.
Here the spin-averaged cross section ratios computed with k ′2 = m 2 h /ζ h are compared with those for massless partons, k ′2 = 0, as used in Ref. [40]. For the production of pions, the dependence on the quark virtuality is negligible at Q 2 = 5 GeV 2 , but becomes evident at lower Q 2 for small z h , z h 0.5. Overall, the σ h /σ

B. Mass corrections for specific experiments
The relevance of the HMCs to future SIDIS experiments is illustrated in Fig. 10, where the difference between the finite-Q 2 and Bjorken limit cross sections, is evaluated relative to the finite-Q 2 cross section at kinematics typical of planned 11 GeV Jefferson Lab experiments [61][62][63]. The cross sections are computed using the same spinaveraged PDF (CT [55]) and fragmentation function (HKNS [59]) parametrizations discussed in Sec. III A above.
The effects are pronounced mostly at large z h , where the ratio of uncorrected to corrected cross sections σ The effect of HMCs on the actual (reduced) spin-averagedσ h as well as spin-dependent ∆σ h cross sections is shown in Fig. 11, where the cross sections are calculated at the same kinematics as in Fig. 10. The impact of HMCs on the cross sections are more pronounced for increasing x B and decreasing Q 2 .
For spin-dependent scattering, the effects on the semi-inclusive polarization asymmetry A h 1 can also be quantified by defining the difference with respect to the massless limit asymmetry A where (neglecting the transverse g 2 contribution) A h 1 = A h /D. Here A h = ∆σ h /2σ h is the parallel asymmetry, and D = (1 − (1 − y)ε)/(1 + εR) is the photon depolarization factor, with R the ratio of longitudinal to transverse photoproduction cross sections. At leading order, the polarization asymmetry is then given in terms of the ratios of sums of polarized and unpolarized PDFs, where the kinematic prefactor √ 1 − ε 2 /D = 1 + γ 2 /(1 + γ 2 y/2). In the Bjorken limit The corrections to the asymmetry δ (HMC) A h 1 as well as to the actual A h 1 asymmetry are shown in Fig. 12 as a function of z h for the production of π + mesons from a proton target, at fixed x B and Q 2 values corresponding to planned 11 GeV Jefferson Lab experiments [63].
At low x B values the differences are very small except at very small z h , where the effects increase as z h → 0. As for the spin-averaged cross sections in Fig. 10, the effects increase with increasing x B and with decreasing Q 2 . At the highest x B value, x B = 0.725, the asymmetry is reduced by ≈ 0.03 for z h 0.3. If the asymmetry were to approach unity in the x B → 1 limit, this would amount to a correction of ≈ 3% − 4%. Somewhat larger corrections are obtained using the GJR [57] and BB [58] PDF parametrizations, although for this combination the A h 1 asymmetry is not guaranteed to respect the positivity constraint at large x B . On the other hand, any dependence of the HMCs on the input PDFs and fragmentation functions is in principle artificial, since in the actual global analyses of SIDIS data the distributions would be determined uniquely and self-consistently by implementing the hadron mass corrected expressions in Eqs. (38) directly.

IV. CONCLUSION
With the imminent completion of the 12 GeV upgrade of Jefferson Lab, and with ongoing programs at existing facilities, a new generation of SIDIS experiments will vastly improve our understanding of the spin and flavor decomposition of parton distributions in the nucleon, and explore the relatively new domain of transverse momentum dependent parton distributions. A full utilization of the new data will require that effects associated with kinematical constraints at finite energy are properly taken into account. Following earlier work which studied the dependence of unpolarized SIDIS cross sections on the masses of hadrons in the initial and final states [41][42][43], in this work we have presented a comprehensive analysis of hadron mass corrections to both spin-averaged and spin-dependent cross sections and asymmetries at finite values of Q 2 .
Using the framework of collinear factorization, we have derived formulas for SIDIS cross sections in the presence of HMCs, which at leading order in α s result in a rescaling of the PDFs in terms of the modified Nachtmann variable ξ h and of the fragmentation functions in terms of the finite-Q 2 fragmentation variable ζ h . Our results respect all kinematical limits at finite Q 2 , and reproduce the standard expressions in the Bjorken limit. An interesting feature of the modified formulas is that, in the presence of HMCs, the parametric dependence on the scattering and fragmentation variables in the parton distribution and fragmentation functions becomes kinematically intertwined. While this complicates the analysis of SIDIS data in certain regions of kinematics, the effects are calculable within our framework.
We have quantified the hadron mass effects numerically as a function of the kinematic variables x B , z h and Q 2 in order to determine the regions where the corrections are most relevant. Generally, the HMCs are strongest at large x B and low Q 2 (as for target mass corrections in inclusive DIS), and for large as well as very low values of the fragmentation variable z h . The effects are also more dramatic for heavier hadrons such as kaons than for pions at the same kinematics. Extending the previous analysis of mass corrections in unpolarized SIDIS [43], we have also examined the dependence of the HMCs on the choice of scattered parton virtuality, k ′2 . In some cases the uncertainty due to this choice is quite significant, particularly at small values of Q 2 and low z h for pions, and over a larger z h range for kaons, where the correction is observed to change sign.
To illustrate the importance of HMCs in practical applications, we have computed the corrections to SIDIS cross sections and polarization asymmetries that would need to be applied at kinematics relevant to upcoming experiments at Jefferson Lab [61][62][63]. Here the x B and Q 2 are necessarily correlated, so that usually the data bins at small x B correspond to lower Q 2 values, while at large x B the Q 2 is typically higher. For unpolarized pion production, the HMCs are strongest at large z h , for all kinematics. At intermediate z h values the corrections at low x B are relatively small, 10%, but increase to ∼ 40% − 50% at higher x B (x B = 0.725), even at moderately large Q 2 (Q 2 ≈ 9 GeV 2 ). Qualitatively similar behavior is observed for the semi-inclusive polarization asymmetry A h 1 , which receives larger HMCs at higher x B values, although there is stronger sensitivity to the specific behavior of the input PDFs. Overall, our analysis suggests that mass corrections may be an important ingredient in future analysis of SIDIS data from facilities such as Jefferson Lab, especially at high values of x B , and particularly for hadrons heavier than the pion.
An immediate application of the results derived here will be in upcoming global spin PDF analyses, such as by the JAM Collaboration [27], which aims to fit an expanded set of high-energy scattering data, including SIDIS, down to Q 2 = 1 GeV 2 . Future theoretical development of this work should include extending the calculation to next-to-leading order in α s , which will necessitate consideration of hadron production at nonzero transverse momenta, p h⊥ . The work can also be extended to other types of distributions measured in various single-spin asymmetries in SIDIS reactions, which provide information not only on the x B and z h distributions but also on the transverse momentum of the partons.