The twist-three distribution $e^q(x,k_\perp)$ in a light-front model

We discuss the twist-three, unpolarized, chiral-odd, transverse momentum dependent parton distribution (TMD) $e^q(x,k_\perp)$ within a light-front model. We review a model-independent decomposition of this TMD, which follows from the QCD equations of motion and is given in terms of a leading-twist mass term, a pure interaction-dependent contribution, and singular terms. The leading-twist and pure twist-three terms are represented in terms of overlap of light-front wave functions (LFWFs), taking into account the Fock states with three valence quark ($3q$) and three-quark plus one gluon ($3q+g$). The $3q$ and $3q+g$ LFWFs with total orbital angular momentum zero are modeled using a parametrization derived from the conformal expansion of the proton distribution amplitudes, with parameters fitted to reproduce available phenomenological information on the unpolarized leading-twist quark and gluon collinear parton distributions. Numerical predictions for both the quark TMD $e^q(x,k_\perp)$ and the collinear parton distribution $e^q(x)$ are presented, discussing the role of the quark-gluon correlations in the proton.


I. INTRODUCTION
Higher-twist parton distributions, especially when the dependence on the parton transverse momenta is taken into account, give access to a wealth of information about the nucleon parton structure [1][2][3]. They describe multiparton correlations inside the nucleon, corresponding to the interference between scattering from a coherent quark-gluon pair and from a single quark [4][5][6][7]. As such, they help understanding the quark-gluon dynamics inside the hadrons, and go beyond the probabilistic interpretation that applies to the leading-twist parton distribution functions.
Twist-three parton distributions functions (PDFs) and transverse momentum dependent parton distributions (TMDs) contribute to various observables in inclusive and semi-inclusive deep inelastic scattering (SIDIS), respectively. Although suppressed with respect to twist-two observables, twist-three structure functions are not small in the kinematics of fixed target experiments. One of the priority tasks of the future experimental program at JLab12 is the measurement of different higher-twist spin-azimuthal asymmetries in SIDIS [8][9][10]. A future electron ion collider would extend such experimental investigation by accessing different kinematical regions [11,12].
In this context, model studies have been shown to have important impact for the understanding of TMDs and the theoretical interpretation of related observables (see, e.g., Ref. [13]). Higher-twist quark PDFs and TMDs can in general be decomposed into contributions from leading-twist mass terms, singular terms and pure interactiondependent ("tilde") terms. This decomposition is obtained through the QCD equations of motion (EOM) and allows one to single the tilde term out as the contribution of quark-gluon correlation functions. Neglecting the tilde and mass terms is referred to as Wandzura-Wilczek approximation [14]. This approximation has been often used as starting point to simplify the description of twist-three SIDIS observables [15][16][17], and showed to be an useful numerical approximation [18]. However, there is no real experimental evidence of its validity, and it misses one of the main motivation to study sub-leading twist, i.e. the non-perturbative physics of quark-gluon correlations.
In this work we want to make a step forward with respect to quark-model calculations and take into account explicitly the contribution from intrinsic gluon degrees of freedom. To this aim, we use the language of light-front wave functions 1 (x, k ⊥ ) and the twist-three quark TMD e q (x, k ⊥ ), and we review the general decomposition of e q (x, k ⊥ ) derived from the QCD EOM. In Sec. III, we present the Fock-state expansion of the proton state, and describe a model-independent representation of the LFWFs for the 3q and 3q + g components with zero partons' orbital angular momentum, as derived originally in Refs. [41,43]. Then, we construct the corresponding LFWF overlap representation for the pure twist-three contribution to e q (x, k ⊥ ) and introduce our parametrization for the LFWFs, in terms of proton DAs. In Sec. IV, we fix the model parameters of the LFWFs by fitting the quark and gluon PDF f 1 (x) to the MMHT2014 parametrization [44] at the scale of 1 GeV 2 , and discuss our model predictions for both the TMD e q (x, k ⊥ ) and the PDF e q (x). We summarize the work in Sec. IV C. Technical details about the derivation of the model-independent decomposition of e q are given in App. A, and the expression of e q (x, k ⊥ ) in terms of our model LFWFs can be found in App. B. Quark TMDs are defined from the twist expansion of the following quark-quark correlators (see, e.g., Refs. [3,45]).
with k + = xP + 1 . The target state is characterized by its four-momentum P and the covariant spin vector S. In Eq. (1), ψ is the quark field and W is an appropriate Wilson line, that connects the bilocal quark operators and ensures gauge invariance [46]. It is defined as: where [a + , a − , a ⊥ ; b + , b − , b ⊥ ] denotes a gauge link connecting the points a µ = (a + , a − , a ⊥ ) and b µ = (b + , b − , b ⊥ ) along a straight line. At twist-three level, we find three unpolarized T-even quark TMDs, i.e. the twist-two TMD f q 1 (x, k ⊥ ), and the twist-three TMDs e q (x, k ⊥ ) and f ⊥ (x, k ⊥ ). In the following, we restrict ourselves to discussing the TMDs f 1 and e, which do not involve partons' orbital angular momentum transfer between the initial and final states. They are defined in terms of the quark-quark correlator as 1 We use light-front coordinates, with v ± = 1 √ 2 (v 0 ± v 3 ) and v ⊥ = (v 1 , v 2 ) for a generic four-vector v.
where P | . . . |P denotes the target spin-averaged matrix element and M is the nucleon mass. The quark collinear PDFs f q 1 (x) and e q (x) are obtained by integrating the corresponding TMDs over the transverse momentum k ⊥ . For later convenience, we also introduce the gluon unpolarized twist-two PDF f g 1 (x), defined in terms of the gluon correlator Φ g,ji (x, k ⊥ ) as where F µν a is the gluon field strength tensor and a is the color index. We now focus our attention on the quark TMD e q (x, k ⊥ ), and review its general decomposition in terms of leadingtwist quark-mass terms, singular terms and pure interaction-dependent contributions. The bilocal quark operator entering the definition in Eq. (4) can be rewritten as where we introduced the projection of the quark fields into the light-cone 'good' and 'bad' components, i.e., ψ + = P + ψ = ψ + and ψ − = P − ψ, respectively, with P ± = 1 2 γ ∓ γ ± . The bad component ψ − is a constrained field, as follows from the QCD EOM: where the covariant derivative is defined as with g s the strong coupling constant. If we assume that the plus component k + of the quark's momentum is strictly positive, we can invert the EOM (7) in a straightforward way using the Fourier expansion of the fields. Instead, problems arise when we include the contribution from zero modes corresponding to k + = 0. In this case, there appear singularities in the bad components of the field, and one needs a regularization prescription [47,48]. To this aim, we follow the procedure outlined in Refs. [6,49] and use the EOM (7) to derive the following operator identity where we defined . Using the expression in Eq. (9) for the bad component, Eq. (6) can be rewritten as with where the index j = 1, 2 labels the transverse component. The bad components ψ − contributes only to O s in Eq. (11). This operator, when inserted in the matrix element of Eq. (4) and integrated over z − , gives a singular contribution proportional to δ(x): This contribution is well known for the PDF e q (x), being related to the pion-nucleon-sigma term (see, e.g., Ref. [6]). The contribution to e q from the operator O m can be worked out using the Fourier expansion of the matrix element of the operator (12) and the definition of f q 1 in (3), as outlined in App. A, with the result where the singular term is a natural consequence of the divergences associated with the zero modes (x = 0). Limiting ourselves to the T-even sector and to the target-spin averaged matrix element, the contribution from the operator O tw3 can be rewritten as where G µν is the gluon field strength tensor. Using the results in App. A, Eq. (16) can be recast in the form whereẽ is a pure twist-three contribution defined in terms of the quark-gluon-quark correlation function [50] Φ Collecting the results in Eqs. (14)- (17), we end up with the following decomposition: The singular term beyond the contribution of e s is usually not discussed in literature and, to the best of our knowledge, has never been written explicitly in this form. The decomposition (20) is independent on the choice of the gauge. In the light-cone gauge A + = 0, with suitable boundary conditions at light-cone infinity for the transverse components of the gauge field, the gauge links in the correlators can be ignored and Φ A,j in Eq. (19) is replaced by the following correlator [50,51] In the framework of light-front quantization in the A + = 0 gauge, the effects of the final-state interactions associated with the gauge link can be reabsorbed in the LFWFs of the target, which acquire an imaginary phase [52]. As we are interested to T-even TMDs, these effects can be ignored, and in the following we will work only with real LFWFs.
Integrating Eq. (20) over k ⊥ , one obtains the corresponding decomposition for the PDF e q (x) from which one can easily infer well-known relations for the first Mellin moments of e q (x) [6].

III. LIGHT-FRONT FOCK-STATE EXPANSION
In this section, we derive the LFWF overlap representation for the TMD e q (x, k ⊥ ). Since we will work under the assumption k + > 0, in the following we will not consider the singular terms in Eqs. (20) and (22).
In the framework of light-front quantization in the A + = 0 gauge, and restricting ourselves to the contributions from the 3q and 3q + g Fock-states, the light-front Fock-state expansion of a proton state with momentum P and light-front helicity Λ reads with In Eqs. (24) and (25), Ψ Λ 3q and Ψ Λ 3q+g are, respectively, the LFWF for the N = 3 and N = 4 parton Fock state with λ i the parton light-front helicity, q i = u, d and g the quark and gluon flavor index, c i the parton color index, and k i the parton momentum. For the argument of the LFWFs, we used a collective notation, . Furthermore, the sum over the color indexes is understood, and, then, using also the sum over the flavor index, the color matrix in Eq. (25) can be saturated with the color index of the first quark only.
The integration measures are defined as: The partial contribution of the N -parton Fock state is defined as where P N is the probability to find the N -parton state in the proton. Using the expressions of Eqs. (23)-(24) for the proton state in the definition (18), the LFWF overlap representation Eq. (30) involves the overlap of LFWFs for the 3q and the 3q + g Fock states, giving direct information on the quark-gluon correlations inside the proton. The momentum of the active quark in the 3q LFWF is set to the external momentum (k + = xP + , k ⊥ ) that is also the sum of the gluon momentum and one of the quark momentum in the 3q + g LFWF. This makes evident the partonic interpretation of theẽ term as an interference between scattering from a coherent quark-gluon pair and from a single quark [4][5][6].
The corresponding results for the LFWF overlap representation of the twist-two contribution to e q , related to the unpolarized distribution f q 1 , are given in terms of the sum of the squares of the N -parton LFWFs, according to the density interpretation of the twist-two distributions (see, e.g., Ref. [53]).

A. Relation to nucleon distribution amplitudes
In this section, we present a parametrization of the proton LFWFs in terms of leading-twist (twist-three) and nextto-leading-twist (twist-four) proton DAs. To this aim, we consider the component of the proton state corresponding to vanishing total orbital angular momentum of the partons [41,43,54] where g ↑ (g ↓ ) denotes the gluon state with positive (negative) light-front helicity.
With respect to the expansion in terms of LFWFs of Eqs. (24) and (25), here the flavour and helicity structure of the parton composition is made explicit, and the light-front wave amplitudes (LFWAs) ψ (j) are scalar functions, which depend only on the parton momenta, i.e. the argument i = 1, 2, 3, 4 stands fork i . The dependence of the LFWAs on the factorization scale is implicit. The proton state with negative helicity is obtained from Eqs. (31)-(33) by applying the light-front parity transformation Y , which corresponds to a parity operation followed by a 180 o rotation around the y axis. It acts on a state of momentum P and light-front helicity Λ as where s is the total spin of the state and η is the intrinsic parity of the hadron (η = +1 for a proton and a quark state, η = −1 for a gluon state). The representation ofẽ q (x, k ⊥ ) in terms of overlap of the LFWAs in Eqs. (31)-(33) is reported in App. B. Following Ref. [41], we can write the LFWAs using the following factorized form where the k ⊥ -dependence of the functions Ω N is assumed to be Gaussian, i.e.
with the normalization: The pure x i -dependent part of the LFWAs in Eqs. (35)- (38) can be directly related to the proton DAs [38,41,55].
In particular, the LFWA which enters the 3q component of the proton state coincides with the leading twist-three proton DA, i.e.
The LFWAs for 3q + g Fock-state are related to the next-to-leading-oder twist-four DAs [56] by 3 The distribution amplitudes can be expanded into a basis of orthogonal polynomials, leaving all the factorization-scale dependence in the coefficients of the expansion. For the three-quark DA, we keep the first few terms of the conformal expansion, corresponding to [56] For the twist-four DAs, we adopt the asymptotic form, given by The normalization constant f N (µ) in Eq. (45) and the parameters λ q i (µ) of the expansion in Eqs. (46)- (48) have been evaluated using QCD sum rules [41] at the scale µ = 1 GeV and in the chiral limit of vanishing quark masses, while the "shape" parameters A, B in Eq. (45) have been calculated on the lattice from the DA moments [57,58].
In our approach, we introduce an explicit dependence on the mass of the "constituent" partons. This amounts to effectively taking into account the non-perturbative nature of the proton state modelled in terms of few Fock-state components. Accordingly, we use the Brodsky-Huang-Lepage prescription [42], and modify the k ⊥ -dependent part of the LFWA by replacing k 2 ⊥ with k 2 ⊥ + m 2 , i.e., we use For the x-dependent part of the LFWAs, we keep the same expressions in terms of proton DAs as in Eqs (45)- (48).
Since we moved away from the chiral limit, we do not use the results of QCD sum rule and lattice QCD for the coefficients of the DAs. Instead, we treat them as free parameters, as specified in the following.

A. Fixing the parameters from twist-two parton distributions
The parameters of the LFWFs are fixed to reproduce available phenomenological parametrizations for the unpolarized PDF f 1 . In particular, we used the results of the MMHT2014 parametrization [44] at Q 2 = 1 GeV 2 , i.e. the same scale of the results from lattice and QCD sum rules. The valence-quark contribution to f 1 is fitted in the range x ∈ (0, 1), while we limit to x ∈ (0.4, 1) for the fit of the gluon unpolarized PDF. In the last case, we did not include lower values of x, because in that region our nonperturbative constituent parton model is not able to catch the low-x dynamics of the gluons and, at the same time, the phenomenological parametrizations are not so well constrained at low x and low scales.
In the fit procedure, the renormalization constant f N of the leading-twist DA in Eq. (45) and the parameters λ g 1 , λ g 2 and λ g 3 of the next-to-leading twist DAs in Eqs. (46)- (48) are considered to be distributed as Gaussian variables with means and standard deviations given by the estimates of the QCD sum rules and lattice calculations. The shape parameters A and B in Eq. (45) are sampled as uniform variables within a larger range than the error band of the lattice calculations. This allows us to having more flexibility to take into account the effects of finite quark and gluon masses. The parton masses along with the width of the k ⊥ distributions enter only in the parametrization of the functions Ω N in Eq. (49). The parton masses are fitted with the constraint Furthermore, we assume the same width for the k ⊥ -distributions of the gluon with positive and negative helicity, i.e.
and we take a 4 along with the width a 3 for the 3q state as additional free parameters. The fitted values of the coefficients of the leading-twist and next-to-leading-twist DAs are shown in Table I, in comparison with the values obtained from QCD sum rules in the chiral limit and lattice calculations. For the quark and gluon mass we obtain m q = 0.161 GeV, m g = 0.050 GeV.
For the widths, the fit results are which correspond to have a quark-gluon state slightly more compact in the transverse-momentum space as compared to the valence three-quark configuration. We note that the two sets of parameters in Table I are consistent, within the uncertainty bands, except for the coefficients A and B. However, we found very small sensitivity of the fitted PDFs to the A and B parameters.  Table I and in Eqs. (52) and (53). Dashed curve: results from the parametrization of Ref. [44], with the corresponding uncertainty bands.
The fit results for the unpolarized parton distributions for the valence up and down quarks, and for the gluon are shown in Fig. 1, in comparison with the MMHT2014 phenomenological parametrizations [44] at Q 2 = 1 GeV 2 . The results for the quark PDFs reproduce very well the MMHT2014 parametrization at larger x (from x ≥ 0.2, in the case of the up quark, and x ≥ 0.4, for the down quark), while they have a much faster fall-off for x → 0 than the Regge-motivated behavior of the parametrization. The results for the gluon PDF, in the fitted range of x ≥ 0.4 reproduce well the MMHT2014 parametrization. Furthermore, for the probability of the 3q and 3q + g components of the proton state we find the following results which are consistent with the values of Ref. [41].

B. Twist-three distributions
Using the explicit parametrization for the LFWAs given in Sec. IV A, and the LFWF overlap representation in App. B, we obtain the results shown in Figs. 2 for the PDF xe q (x) of up and down quarks. The red short-dashed curves show the contribution from the twist-two term, while the blue long-dashed curves correspond to the genuine twist-three contributions. The relative size of twist-two and genuine twist-three contributions depends on the quarkmass parameter, which enters as proportionality constant that weights the twist-two term in Eq. (22). In our model calculation, the partons' mass also enter the functions Ω N in Eq. (49). However, in this case, the dependence on the partons' mass is such that it does not affect the relative size of the the twist-two and twist-three contributions to xe q (x), and slightly changes their behaviour as function of x.
We also note that the twist-two and genuine twist-three terms have a quite different x-dependence. The twist-two contribution is peaked at x ≈ 0.2, with a fast fall-off at larger x, more pronounced in the case of up quarks than down   [57,58] for the coefficients of the twist-three DA in Eq. (45), parametrizing the 3q LFWAs, and the QCD sum rules estimates [41] of the coefficients of the twist-four DAs in Eqs. (46)- (48), parametrizing the 3q + g LFWAs. quarks. On the other side, the xẽ q (x) contribution is peaked at x ≈ 0.5, and becomes the dominant contribution at larger x, especially for down quarks. We also note that the pure twist-3 contributions for the up and down quarks have very similar size, whereas the twist-2 contribution for the up quark is approximately twice as large as the twist-2 contribution for the down quark (note the different scales on the vertical axis).
Recently, the CLAS collaboration has reported preliminary results of a measurement of the beam asymmetry in dihadron SIDIS, using a longitudinally polarized 6 GeV electron beam off an unpolarized proton target [10,59]. These data have been analyzed to extract the following flavor combination of the valence-quark contributions to e q (x): In Fig. 3, the preliminary CLAS data points at the scale Q 2 = 1.5 GeV 2 are compared with our model predictions at the scale Q 2 = 1 GeV 2 . The model results are also split in the contributions from twist-two and genuine twist 3-terms, corresponding to the short-dashed blue curve and the long-dashed cyan curve, respectively. Our results are in quite good agreement with the experimental extraction at the two higher-x bins, but they are not able to reproduce the observed fast rising at lower x. This could be due to a lack of our model, that, according to the fit results for the unpolarized PDF f 1 , shown in Fig. 1, is less reliable in the lower-x region. We also notice that the genuine twist-three contribution in the considered x-range is very small, supporting the results within the light-front constituent-quark picture that was used in Ref. [21] and showed to be able to reproduce the results of the CLAS data at higher x. However, one should bear in mind that these data are still preliminary and have unestimated systematic uncertainties.

C. Transverse momentum
Next, we turn our attention to the k ⊥ dependence of the TMDs. We define the x-dependent mean squared transverse momentum of a generic TMD j (x, k ⊥ ) as follows The corresponding results for the quark and gluon unpolarized distribution f 1 (x, k ⊥ ) and for the quark TMDẽ (x, k ⊥ ) are shown in the left and right panel of Fig. 4, respectively. The quark results refer to the valence-quark contribution, since in our model we neglect the sea quarks. The results for the gluon are obtained from the 3q + g component of the proton LFWF, and therefore correspond to the intrinsic non-perturbative gluon contribution. All the results refer to the model scale of Q 2 = 1 GeV 2 . We find that in the case of f 1 , the size and the x-dependence is very similar for up and down quarks, and for gluons. They are all peaked at x ≈ 0.3, and fall down rapidly at higher x, with a very similar slope in the case of the down quarks and the gluon. The similar behavior for down quark and gluons is an indication that the results for the down quark contribution to the f 1 TMD are more dominated by the 3q + g component of the LFWF than in the case of up quark. For the quark contribution, there exists a recent extraction, based on a fit of the unpolarized TMD f 1 (x, k ⊥ ) to available experimental data measured in SIDIS, Drell-Yan and Z boson production [60]. This extraction assumes no quark-flavour dependence and uses data in the range of 5 · 10 −2 ≤ x ≤ 0.4. Therefore, the corresponding results for the mean squared transverse momentum get a sizeable contribution from sea-quarks, which results in a x-dependence quite different from the bell-shaped structure we find in our calculation. Nevertheless, we find that our results are in the same range of values as the phenomenological extraction.
In the case of the mean squared transverse momenta ofẽ q (x, k ⊥ ) we observe a more pronounced quark flavor dependence w.r.t. to the case of the unpolarized twist-two TMD. Furthermore, they are peaked at higher-x values, with the contribution from up quarks slightly shifted at larger x w.r.t. to the one from down quarks.
By integrating Eq. (57) over x, we obtain the results for the k 2 ⊥ width of the TMDs. The results for f 1 andẽ at the model scale of Q 2 = 1 GeV 2 are presented in Table II. They are very similar for all the partons in the case of f 1 , while they have an appreciable quark-flavor dependence in the case ofẽ. We also notice the broader widths for theẽ distribution with respect to f 1 .  Finally, in Fig. 5 we show the three-dimensional plots ofẽ for the up and down quark, as function of x and k 2 ⊥ . For both quark flavors, the k 2 ⊥ -dependence of the distributions slightly moves away from a Gaussian shape, as expected from our model Ansatz for the k 2 ⊥ -dependence of the LFWFs in Eq. (49).

Conclusions
In this work, we studied the sub-leading structure function e q . First, we reviewed its model-independent decomposition, which follows from the QCD EOM. In particular, this decomposition contains a genuine twist-three part, i.e. the "tilde" term encoding the quark-gluon correlations, a pure twist-two contribution and a singular (δ-like) term. The singular term is in turn given by a well-know contribution that can be related to the pion-nucleon-sigma term, and an additional term that, instead, is poorly discussed in literature and has been written down explicitly here for the first time.
Then, we focused on the modeling of the tilde term. To this aim, we constructed a model-independent representation in terms of overlap of LFWFs for the 3q and 3q + g Fock-components of the proton state. For the calculation, we restricted ourselves to the LFWFs corresponding to zero orbital angular momentum of the partons. This allowed us to relate the 3q and the 3q + g LFWFs in the limit of zero-transverse separation to the leading-twist and nextto-leading-twist DAs of the proton, respectively. The transverse-momentum dependence of the LFWFs was built up by assuming a modified Gaussian Ansatz, with explicit dependence on the mass of the partons. The DAs and the transverse-momentum dependent part were parametrized in such a way to reproduce available phenomenological parametrizations for the unpolarized PDF f 1 of quarks and gluons at the scale of Q 2 = 1 GeV 2 . The fit gave values for the parametrization of the DAs that are consistent with lattice calculations [57,58] and predictions from QCD sum rules [41].
With these ingredients, we provided predictions for both the quark PDFẽ q (x) and the quark TMDẽ q (x, k ⊥ ), in comparison with the corresponding twist-two contribution given in terms of f q 1 . In the case ofẽ q (x), we found that the pure twist-three contribution has a distinctive x-dependence from the twist-two contribution, and is about 20-30% of the twist-two contribution at the peak position, using our model-results for the quark masses. We also considered preliminary results from a phenomenological extraction of a particular flavor combination of the valencequark contribution to e q (x) [59]. For this quark-flavor combination, our model predictions are almost saturated by the twist-two contribution to e q (x), and showed a quite good agreement with the extracted results in the large x-region.
The results of this work are encouraging for an extension to other subleading-twist TMDs. However, twist-three TMDs other than e q involve a transfer of orbital angular momentum between the initial and final proton states. This requires the modeling of the LFWF components with non-zero orbital angular momentum. Work in this direction is in progress.