Saturation effects in exclusive rhoT, rhoL meson electroproduction

We use recent results for the gamma*L ->rhoL and gamma*T ->rhoT impact factors, computed in the impact parameter representation within the collinear factorization scheme, to get predictions for the polarized cross-sections sigmaT and sigmaL of the diffractive leptoproduction of the rho meson at high energy. In this approach the helicity amplitude is a convolution of the scattering amplitude of a color dipole with a target, together with the virtual gamma wave function and with the first moments of the rho meson wave function (in the transverse momentum space), given by the distribution amplitudes up to twist 3 for the gamma*T ->rhoT impact factor and up to twist 2 for the gamma*L ->rhoL impact factor. Combining these results with recent dipole models fitted to DIS data, which include saturation effects, we show that the predictions are in good agreement with HERA data for photon virtuality (Q**2) larger than typically 5 GeV**2, without free parameters and with a weak dependence on the choice of the factorization scale, i.e. the shape of the DAs, for both longitudinally and transversely polarized rho meson. For lower values of Q**2, the inclusion of saturation effects is not enough to provide a good description of HERA data. We believe that it is a signal of a need for higher twist contributions in the rho meson DAs. We also analyze the radial distributions of dipoles between the initial gamma* and the final rho meson states.

Abstract: We use recent results for the γ * L → ρ L and γ * T → ρ T impact factors, computed in the impact parameter representation within the collinear factorization scheme, to get predictions for the polarized cross-sections σ T and σ L of the diffractive leptoproduction of the ρ meson at high energy. In this approach the helicity amplitude is a convolution of the scattering amplitude of a color dipole with a target, together with the virtual gamma wave function and with the first moments of the ρ meson wave function (in the transverse momentum space), given by the distribution amplitudes up to twist 3 for the γ * T → ρ T impact factor and up to twist 2 for the γ * L → ρ L impact factor. Combining these results with recent dipole models fitted to DIS data, which include saturation effects, we show that the predictions are in good agreement with HERA data for photon virtuality (Q 2 ) larger than typically 5 GeV 2 , without free parameters and with a weak dependence on the choice of the factorization scale, i.e. the shape of the DAs, for both longitudinally and transversely polarized ρ meson. For lower values of Q 2 , the inclusion of saturation effects is not enough to provide a good description of HERA data. We believe that it is a signal of a need for higher twist contributions in the ρ meson DAs. We also analyze the radial distributions of dipoles between the initial γ * and the final ρ meson states.

Introduction
In the high energy limit, exclusive processes and more particularly the diffractive leptoproduction of vector mesons, provide a nice probe to study the hadronic properties. From the experimental side, data have been extracted in a wide range of center-of-mass energies, from a few GeV at JLab to hundreds of GeV at the HERA collider. The kinematical range which is at the heart of the present paper is the large energy in the center-of-mass of γ * p, denoted W , for which the HERA collaborations H1 and ZEUS measured ρ-meson electroproduction starting from the early period of HERA activity [1,2] till now, with an increasing precision leading to a complete analysis [3,4] of spin density matrix elements, polarized and total cross-sections describing the hard exclusive productions of the ρ and the φ vector mesons V in the process (1.1) These matrix elements and polarized cross-sections can be expressed in terms of helicity amplitudes T λ V λγ (λ γ , λ V : polarizations of the virtual photon and of the vector meson). The ZEUS collaboration [3] has provided data for different photon virtualities Q 2 , i.e. for 2 < Q 2 < 160 GeV 2 , 32 < W < 180 GeV ( |t| < 1 GeV 2 ), while the H1 collaboration [4] has analyzed data in the range 2.5 < Q 2 < 60 GeV 2 , 35 < W < 180 GeV (|t| < 3 GeV 2 ) . The high virtuality of the exchanged photon allows the factorization of the amplitude into a hard subprocess described within the perturbative QCD approach and suitably defined hadronic objects, such as the dipole-nucleon scattering amplitude or the vector meson wave functions and distribution amplitudes (DAs) [5][6][7]. HERA data are thus interesting observables to test the properties of these non-perturbative objects such as the saturation dynamics of the nucleon or the transverse momentum dependence of the vector meson wave functions.
On the theoretical side, three main approaches have been developed. The first two, a k T -factorization approach and a dipole approach, are applicable at high energy, W ≫ Q ≫ Λ QCD . They are both related to a Regge inspired k T -factorization scheme [8][9][10][11][12][13][14], which basically writes the scattering amplitude in terms of two impact factors: one, in our case, for the γ * − ρ transition and the other one for the nucleon to nucleon transition, with, at leading order, a two "Reggeized" gluon exchange in the t-channel. The Balitsky-Fadin-Kuraev-Lipatov (BFKL) evolution, known at leading order (LLx) [15][16][17][18] and next-to-leading (NLLx) order [19][20][21][22], can then be applied to account for a specific large energy QCD resummation. The dipole approach is based on the formulation of similar ideas althought not in k T but in transverse coordinate space [23,24]; this scheme is especially suitable to account for nonlinear evolution and gluon saturation effects. The third approach, valid down to W ∼ Q, was initiated in [25] and [26]. It is based on the collinear QCD factorization scheme [27,28]; the amplitude is given as a convolution of quark or gluon generalized parton distributions (GPDs) in the nucleon, the ρ-meson DA, and a perturbatively calculable hard scattering amplitude. GPD evolution equations resum the collinear quark and gluon effects. The DAs are also subject to specific QCD evolution equations [29][30][31].
Though the collinear factorization approach allows us to calculate perturbative corrections to the leading twist longitudinal amplitude (see [32] for NLO), when dealing with transversely polarized vector mesons related to higher twist contributions, one faces endpoint singularity problems. Consequently, this does not allow us to study polarization effects in diffractive ρ-meson electroproduction in a model-independent way within the collinear factorization approach. An improved collinear approximation scheme [33] has been proposed, which allows us to overcome end-point singularity problems, and which has been applied to ρ-electroproduction [34][35][36][37].
In this study, we consider polarization effects for the reaction (1.1) in the high energy region, s = W 2 ≫ Q 2 ≫ Λ 2 QCD , working within the k T -factorization approach, where the helicity amplitudes can be expressed 1 as the convolution of the γ * → ρ impact factor Φ γ * (λγ )→ρ(λρ) (k 2 , Q 2 ) with the unintegrated gluon density F(x, k 2 ) density, which at the Born order is simply related to the nucleon impact factor Φ N →N (k 2 , Q 2 ), where k is the transverse momentum of the t−channel exchanged gluons. For the transverse amplitude the end-point singularities are naturally regularized by the transverse momenta of the t−channel gluons [38][39][40]. At large photon virtuality, the γ * − ρ impact factors can be calculated in a model-independent way using QCD twist expansion in the region k 2 ≫ Λ 2 QCD . Such a calculation involves the ρ-meson DAs as nonperturbative inputs. The calculation of the impact factors for Φ γ * L →ρ L , Φ γ * T →ρ L is standard at the twist-2 level [41], the next term of the expansion being of twist 4, while Φ γ * T →ρ T was only recently computed [39,40] (for the forward case t = t min ), up to twist-3, including two-and three-parton correlators, which contribute here on an equal footing.
In a previous study [42], we used the results [40,41] for the Φ γ * L →ρ L and Φ γ * T →ρ T impact factors and a phenomenological model [43] for the proton impact factor. It was pointed out that the region k 2 ≫ Λ 2 QCD gives the dominant contribution to the helicity amplitudes while the soft gluon (k 2 < 1 GeV 2 ) contribution cannot be neglected. The soft gluon contribution to the amplitudes involves the interaction of large size color dipole configuration |r| (|r| ≡ 1/|k|) in the fluctuations of the probe and saturation effects could then play an important role.
Following this idea, we have shown in ref. [44] that the helicity amplitudes, expressed in impact parameter space and then computed in the collinear factorization scheme, factorize into the dipole cross-section and the wave functions of the virtual photon combined with the first moments of the ρ meson wave functions parameterized by the DAs, given by the twist expansion up to twist 3 for the production of a ρ T and to twist 2 when producing a ρ L . Note that in [38] a very similar approach was followed, where the DAs were replaced by the Taylor expansion of models of the ρ meson wave function at small qq pair transverse size. The results of ref. [44] link the k T −factorization approach, in particular the results of [40], with the calculations performed in refs. [45][46][47][48][49] within the dipole approach. The main difference between our present approach and the one of refs. [45][46][47][48][49] is that instead of using the light-cone wave functions φ(z, r) which in practice need to be modeled, the amplitude of our approach involves the DAs which parameterize the first moments of the wave functions.
Our main point here is that one can assume the dominant physical mechanism for production of both longitudinal and transversely polarized mesons to be the scattering of small transverse-size quark-antiquark and quark-antiquark-gluon colorless states on the target. This allows to calculate corresponding helicity amplitudes in a model-independent way, using the natural light-cone QCD language -twist-2 and twist-3 DAs.
This paper is organized as follows. In section 2, we introduce the impact factor representation of the helicity amplitudes as well as the kinematics of the process. In section 3, we first recall some results for the γ * L → ρ L and γ * T → ρ T impact factors computed in momentum space respectively up to twist 2 [41] and twist 3 [40] accuracies using the collinear approximation to parameterize the soft part associated to the production of the ρ meson by distribution amplitudes (DAs), calculated in ref. [50,51]. We recall then the expression of these impact factors in the impact parameter space according to the results of ref. [44] allowing to decouple the dipole-nucleon scattering amplitude from the amplitude of production of dipoles in the initial (γ * (λ γ )) and final (ρ(λ ρ )) states. We terminate section 3 by expressing helicity amplitudes and polarized cross-sections in terms of the dipole crosssection. In section 4, we give a brief review of the main properties of the models for the dipole cross-section [52][53][54] or the proton impact factor [43] that we use in our study. We compare our predictions with the data of HERA [3,4] in section 5 and we obtain a good agreement. In this context we discuss the role of higher twist corrections for small Q 2 values. Finally, we analyse the radial distribution of dipole intermediate states involved between the virtual photon and the ρ meson, and discuss the role of the saturation models on the specific example of the Golec-Biernat and Wüsthoff saturation model [52]. We also compare our radial distribution with the overlap of the γ * and ρ meson wave functions obtained in the approach of the dipole models [46,47,55], where the ρ meson wave function is modeled and the parameters are fitted to HERA data.

Helicity amplitudes of the hard ρ meson leptoproduction in the impact factor representation
In the impact factor representation at the Born order, the amplitude of the exclusive process as illustrated in figure 1. The γ * (λ γ ) → ρ(λ ρ ) impact factor Φ γ * (λγ )→ρ(λρ) is defined through the discontinuity of the S matrix element for γ * (λ γ ; q)g(k) → g(k − ∆)ρ(λ ρ ; p ρ ) as where κ = (k + q) 2 . In eqs. (2.1) and (2.2) the momenta q and p ρ are parameterized via Sudakov decompositions, in terms of two lightlike vectors p 1 and p 2 such that 2 p 1 .p 2 = s, as where Q 2 = −q 2 >> Λ 2 QCD is the virtuality of the photon being the large scale which justifies the use of perturbation theory, and m ρ is the mass of the ρ meson. Here −t min denotes the minimal value of −t . The nucleon impact factor Φ N →N in eq. (2.1) cannot be computed within perturbation theory, and is related 2 at Born order to the unintegrated gluon density F(x, k). In the forward limit ∆ ⊥ = 0, the helicity amplitudes read, The impact factors Φ γ * (λγ )→ρ(λρ) (k, Q, µ 2 F ) and the nucleon impact factor vanish at k → 0 or k → ∆, which guarantees the convergence of the integral in eq. (2.4) on the lower limit 3 The computation of the γ * → ρ impact factor is performed within collinear factorization of QCD. The dominant contribution corresponds to the γ * L → ρ L transition (twist 2), while the other transitions are power suppressed. The γ * L → ρ L and γ * T → ρ L impact factors were computed a long time ago [41], while a consistent treatment of the twist-3 γ * T → ρ T impact factor has been performed only recently in ref. [40]. It is based on the light-cone collinear factorization (LCCF) beyond the leading twist, applied to the amplitudes γ * (λ γ )g(k) → g(k − ∆)ρ(λ ρ ), symbolically illustrated in figure 2. Each of these scattering amplitudes is the sum of the convolutions of a hard part (denoted by H and H µ for two-and three-parton contributions, respectively) that corresponds to the transition of the virtual photon into the constituents of the ρ meson and their interactions with off-shell gluons of the t channel, and a soft part (denoted by Φ and Φ µ ) describing the hadronization of the constituent partons into the ρ meson and parameterized up to a given twist by DAs.
2 Normalization of the impact factors differs from [56] by a factor 2π, Φ [56] = 2πΦ Here . For clarity, we have restored in eq. (2.4) the colored indiced carried by the impact factor Φ γ * λγ →ρ λρ ab . 3 This property of the impact factors is universal in the case of the scattering of colorless objects and is related to gauge invariance [57,58]. Figure 2: Two-and three-parton correlators attached to a hard scattering amplitude in the specific case of the γ * → ρ impact factor, where vertical lines are hard t− channel gluons in the color singlet state.

Helicity amplitudes and polarized cross-sections
3.1 Impact factors γ * L,T → ρ L,T In the Sudakov basis, the longitudinal and transverse polarizations of the photon are 4 For t = t min the same parametrization will be used for the ρ-meson polarization with obvious substitutions Q 2 → −m 2 ρ and Q → m ρ . Let us introduce p and n, two light-cone vectors such that p ρ ≈ p at twist 3 and p · n = 1. The polarization of the out-going ρ meson is denoted by e * . For further use, we denote R * ⊥α = ε αλβδ e * λ ⊥ p β n δ . We use the following convention for the transverse euclidean polarization vectors e ± as in eq. (3.1) The DAs needed to parameterize the ρ meson productions involved in the results of [40] for the impact factors with ∆ = 0, according to [40] are: the twist 2 DA ϕ 1 (y; µ 2 F ) associated to the production of the longitudinal ρ meson, the two-parton twist 3 DAs ϕ T 1 (y; µ 2 F ), ϕ T A (y; µ 2 F ) and the three-parton twist 3 DAs B(y 1 , y 2 ; µ 2 F ) and D(y 1 , y 2 ; µ 2 F ) that parameterize the production of a transversely polarized ρ meson. Following [42,44], we also use the combinations, where ζ V 3 (µ 2 F ) and ζ A 3 (µ 2 F ) are the dimensionless coupling constants The scale µ F is the factorization scale involved in the production of the ρ meson, that we put equal to the renormalization scale of the evolution of the DAs. We recall in appendix A.1 4 In ref. [40] we took ǫ ± = ∓ i √ 2 (0, 1, ±i, 0), which we change here for consistency with the usual experimental conventions [59]. some basics features of the chiral even twist 2 and twist 3 DAs present in this approach; more details can be found on the DAs in ref. [40]. We recall also in appendix A.2 the dependence of the DAs on the collinear factorization scale µ F which is driven by the renormalization evolution equations given in ref. [50].
We have recently shown in [44] that these impact factors, expressed in the impact parameter space, read, and where we have respectively denoted in the case of two-parton exchange as y andȳ = 1 − y the longitudinal momentum fractions of the quark and the antiquark, and in the threeparton exchange y 1 ,ȳ 2 = 1 − y 2 and y g = y 2 − y 1 the longitudinal momentum fractions of the quark, the antiquark and the gluon. The transverse displacement vectors of the color dipole configurations are denoted as r for the interacting dipole (two-and three-parton Fock component) and r ′ for the spectator dipole (present only in the three-parton Fock component). Note that |r| in the case of the two-parton component is the transverse size of the quark anti-quark colorless pair. The functions denoted are the Fourier transforms in the transverse plane of the two-parton component hard parts H and the three-parton component hard parts H µ (illustrated in figure 2), projected on the appropriate Fierz structures Γ µ . The computations of the hard parts lead to the following generic expressions, where the function is the scattering amplitude of a quark-antiquark color dipole with the two t−channel gluons, putting apart the color factor T r(t a t b ) = δ ab /2, with a and b color indices and N c the number of colors. The functions ψ are respectively the amplitudes of production of a ρ meson from a quark-antiquark (quark-antiquark gluon) system produced far upstream the target in the fluctuation of the virtual photon. These functions are computed up to twist 3 in the collinear aproximation in ref. [44] and they contain information about the relevant color dipole system that interacts with the target. The functions ψ γ * L,T →ρ L,T (qq) can be expressed in terms of the virtual photon wave functions Ψ where h = ± 1 2 ,h = ± 1 2 denote respectively the helicities of the exchanged quark and anti-quark, and of the combinations of DAs of the ρ meson φ (3.14) In the two-parton approximation, these functions φ ρ L (h,h) (y; µ 2 F ) and φ ρ T , (λρ) (h,h) (y, r; µ 2 F ) parameterize the moments of the wave functions of the ρ meson, i.e. the first terms of the Taylor expansion of the wave functions at small r. Note that in the approach of ref. [38] φ ρ L (h,h) (y; µ 2 F ) and φ ρ T , (λρ) (h,h) (y, r; µ 2 F ) are replaced by the Taylor expansion for small r of the modeled wave functions of the vector mesons. Finally, the functions ψ γ * L →ρ L (qq) and ψ γ * T →ρ T (qq) with two-parton components read The function ψ with three-parton components reads where the function F γ * T describes the fluctuation of the transversely polarized photon into a quark-antiquark-gluon color singlet. The function F γ * T can be expressed in terms of the longitudinally polarized photon wave function

From impact factors to helicity amplitudes and polarized cross-sections
According to eq. (2.4), the helicity amplitudes read, Inserting the expressions for the impact factor Φ γ * λγ →ρ λρ of eqs. (3.8, 3.10), one gets where we have identified the dipole cross-section as defined in ref. [56] σ(x, r) = N 2 c − 1 4 Note that we can separate the T 11 as the WW contribution and the genuine contribution, The formulas (3.23, 3.24) allow us to combine various models of the scattering amplitude of a dipole on a nucleon with the results obtained by twist expansion of the γ * → ρ impact factor. At t = t min the contribution to the cross-sections σ L and σ T are respectively coming from the helicity amplitudes T 00 and T 11 , The t−dependency is expected to be governed by non-perturbative effects of the nucleon, which can be phenomenologically parameterized by an exponential dependence of the differential cross-section This results in the polarized cross-sections The b(Q 2 ) slope has been measured by ZEUS and H1. We will use here quadratic fits of the b(Q 2 ) slope data of ref. [4] to determine the cross-section. In the following section, we will briefly present the dipole models we shall use to compare our predictions with HERA data. Note that we could in principle use a dipole model taking into account skewness effects. This skewness dependence can be implemented along approaches of refs. [60,61], but this subleading physical effects will be neglected in the present study.

Dipole models
In the dipole picture, the DIS cross-section reads [23,24] The low-x saturation dynamics of the nucleon target was first introduced in refs. [52,62] by Golec-Biernat and Wüsthoff (GBW) model to describe the inclusive and diffractive structure functions of DIS, which inspired many phenomenological descriptions of DIS HERA data [47,55,[63][64][65][66][67]. In this model the dipole cross-section readŝ where and it involves three independent parameters {σ 0 , x 0 , λ p }. One can see from eq. (4.1) that since the wave functions are peaked at r ∼ 1 Q , the domain in which the saturation effects are significant is given by To make contact with photoproduction, it is customary [52] to make the following modification of the definition of the Bjorken variable x where m f is an effective quark mass which depends on the flavour f and of the model used to fit the data. This modification is adopted in the following.
In the GBW saturation model, light quark masses are taken to be 0.14 GeV in order to get a good fit of the photoproduction region. The inclusion of the charm contribution, with m c = 1.5 GeV has also been performed in ref. [52]. The normalisation σ 0 results from the integration over the impact parameter b, assuming that the b−dependence of the dipole amplitude factorises as is related to the density of gluon within the target, leading tô This model provides a good description of inclusive and diffractive structure functions for the values of the parameters presented in table 1.  The small x evolution of the dipole cross-section can be modeled as in the refs. [52], or computed numerically from the running coupling Balitsky-Kovchegov (rcBK) equation. Indeed, the x−dependence is driven at small x by perturbative non-linear equations, the Jalilian-Marian-Iancu-McLerran-Weigert-Leonidov-Kovner (JIMWLK) equation [68][69][70][71][72][73] and the Balitsky-Kovchegov (BK) equation [74,75]. The solutions of the two equations are not significantly different and the BK equation is simpler to solve numerically than the JIMWLK equation. The LO-BK equation cannot be used to obtain the low x evolution as it predicts a growth of the saturation scale much faster than the growth expected from the anlysis based on phenomenological models. It was shown in [53,54] that taking into account only the running coupling corrections of the evolution kernel allows to get the main higher order contributions, and to solve the discrepancy between the growths of the saturation scale. A numerical solution of the running coupling BK (rcBK) equation was obtained in refs. [76,77] based on initial conditions inspired by the GBW model [52] and the McLerran-Venugopalan (MV) model [78]. We will denote these numerical solutions for the dipole scattering amplitudes as the Albacete-Armesto-Milhano-Quiroga-Salgado (AAMQS)-model.
where β 0,n f = 11 − 2 3 n f , Λ n f is the QCD scale and C is one of the free parameters of the model. As usual, the scales Λ n f are determined by the matching condition α s,n f −1 (r 2 ⋆ ) = α s,n f (r 2 ⋆ ) at r 2 ⋆ = 4C 2 /m 2 f and an experimental value of α s . In eq. (4.9), the scale Λ QCD is identified with Λ 3 .
The parameters are fitted to the experimental data for the inclusive structure function of DIS and give a good description for the data of the longitudinal structure function Two sets of dipole cross-sections are available for each of the initial conditions. The first set of dipole cross-sections respectively denoted as (a) and (e) for the GBW and the MV initial conditions, are fitted to data by taking only into account the light quarks u, d, s . The second set, respectively denoted as (b) and (f) for the GBW and the MV initial conditions, are fitted taking into account the dipole cross-section of the light quarks σ 0 N l. (x, r) and of the heavy quarks σ h 0 N h. (x, r) . The normalization of the dipole cross-section σ h. 0 for the charm and the bottom quarks are assumed to be equal. Thus the cross-sections read (4.14) In the following parts we will focus mostly on the sets (a) and (b) using the GBW initial condition respectively denoted as AAMQSa and AAMQSb, as the effects due to different initial conditions do not involve significant changes in the numerical results of present study. We present in Tabs. 2 and 3 values of the parameters of the fits obtained in ref. [77].    Note that the models of dipole cross-section we use here involve only a quark anti-quark intermediate state between the initial and final wave functions used to describe DIS process while the impact factor γ * T → ρ T involves also qqg partonic state. The error induced by this approximation should be subleading compare to higher twist corrections or the choice of µ F , but still it is difficult to appreciate quantitatively how much it could impact on the predictions.
For completeness, we will also display predictions using the Gunion and Soper model (GS-model) [43]. This model was used in our first phenomenological study of ref. [42], and it assumes that the hadron impact factor takes the form where A and M are free parameters that correspond to the soft scales of the proton-proton impact factor. As it was discussed in [42], the ratios of helicity amplitudes are well described for M ≃ 1 GeV and the result is not very sensitive to this parameter. Note that this impact factor indeed vanishes when k → 0 or ∆ − k → 0 in a minimal way. Such a model was the basis of the dipole approach of high energy scattering [79] and used successfully for describing deep inelastic scattering (DIS) at small x [80]. The dipole-proton scattering amplitude in the two-gluon exchange approximation, computed in details in appendix A.3, resulting from the GS-model for the proton impact factor is:

Comparison with the HERA data
Our main results are the polarized cross-sections σ T and σ L of the process (1.1), that we compare with the data of the H1 collaboration [4].
In the following plots, experimental errors are taken to be the quadratic sum of statistical and systematical errors. The collinear factorization scale µ F only appears in the scattering amplitudes through the DAs and the coupling constants and unless specified we will assume that µ F depends on the virtuality Q 2 as The figures 3(a) and 3(b) show the full twist 3 predictions respectively for σ T and σ L , using the different dipole models. The results of the predictions for σ T (figure 3(a)) and σ L (figure 3(b)) are in agreement with the data for values of Q 2 larger than approximately 7 GeV 2 for σ T and 7 GeV 2 for σ L depending on the considered dipole model. The results based on the AAMQSa model are giving a slightly better description than with the AAMQSb although both sets give good quality fits of DIS data. The results obtained by GBW model are close to the one obtained by the AAMQSa model; both are in good agreement with the data for Q 2 above 7 GeV 2 . Let us insist on the fact that the agreement for Q 2 ≥ 7 GeV 2 of our predictions with data is non-trivial since all free parameters of the dipole cross-section (normalization, R 0 (x), etc...) are completely fixed by the fit of DIS data, while on the ρ side the normalizations are given by decay constants obtained from QCD sum rules.
Keeping in mind that the results are not too sensitive to the precise choice of the dipole model, below we will focus on the predictions of the AAMQSa model. Some results of the GBW, AAMQSb and GS models are presented in appendix A.4.
In figure 4, we show separately three different contributions: • the full twist 3 (Total) contribution, involving both the WW and the genuine solutions of the DAs.
• the WW contribution, only involving the WW solutions of the DAs.  • the asymptotic (AS) contribution (for µ F → ∞), involving the asymptotic solutions of the DAs. In this limit the genuine contribution vanishes and the WW DAs can be expressed as functions of the asymptotic DA ϕ 1 (y) = 6yȳ.
The results of the predictions for σ T (figure 4(a)) and σ L ( figure 4(b)) are in agreement with the data for values of Q 2 respectively larger than Q 2 min T ∼ 6.5 GeV 2 and Q 2 min L ∼ 5 GeV 2 which confirms that the amplitude factorizes into a universal colour dipole scattering amplitude and that the truncated twist expansion of the ρ meson soft part is justified.
The two scales Q 2 min  indication that the discrepancy between data and our predictions at low Q 2 are mainly due to the higher twist contributions to the impact factor from the meson structure rather than an effect of the saturation dynamics of the nucleon which should be well described at this scale by the saturation models. Let us mention, again that these saturation models are known to fit very well inclusive DIS as well as diffractive DIS data (for GBW) at these low Q 2 values.
The saturation scale, given by Q S = 1/R 0 (x) is of order 1 GeV in the kinematics of HERA. Since our predictions are only consistent with data in the region Q 2 > Q 2 min L,T > Q 2 S , this limitation do not allow us to access the domain Q 2 S Q 2 where saturation effects can be essential.
The predictions are dominated by the WW-contribution and are not very sensitive to the choice of the collinear factorization scale. This will be further discussed in section 6. An estimation of the error on the cross-sections caused by the error bars on the b−slope measurements is obtained by fitting the upper and lower bounds of the b−slope as shown in figure 5, and then by computing the predictions based on these fits as shown in figure 6. Note that we have assumed that the longitudinal b L and the transverse b T slopes are equal to the b−slope of the total cross-section. This assumption is supported by H1 data where the measurements of the difference b L −b T for Q 2 = 3.3 GeV 2 and Q 2 = 8.6 GeV 2 are much smaller than the b−slope value. Let us also emphasize that in this approach we compute the polarized differential cross-sections in the limit t = t min ≈ 0 where only the s-channel helicity conserving (SCHC) amplitudes T 00 and T 11 are non-zero. The contributions of other helicity amplitudes are encoded in the phenomenological t−dependence given in eq. (3.30) and it turns out that data for the total differential cross-section are dominated by a t-region of very small values, with a typical spread given by the scale 1 b ≈ 1 6 GeV 2 . We now compare our predictions with the data for the total cross-section σ, given by the sum σ = σ L + σ T according to ZEUS convention in ref. [3] or σ = εσ L + σ T following H1 notation [4], where ε is the photon polarization parameter 5 ε ≃ (1 − y)/(1 − y + y 2 /2) .  Finally, our analysis also provides predictions for the ratios R and for the spin density matrix element r 04   correction to the ratio r 04 00 due to the amplitude T 01 , for the t−range of H1, to be below 1%. The results are shown in figure 10 for the ratio R and in figure 11 for the spin density matrix element r 04 00 , using AAMQSa model.   represents the amplitude of probability for these states to dissociate into a quark anti-quark color dipole of size r, the quark having a longitudinal momentum fraction y. We define thus the probability amplitude W λρλγ as the corresponding parts of the impact factors appearing in eqs.   for the amplitude T 00 and for the Total, the WW and the genuine contributions to T 11 . The probability amplitudes W λρλγ permit in turn to define the radial distributions P λρλγ of the interacting dipole, as P λρλγ (r, Q 2 , µ 2 F ) = 1 N λρλγ r dy W λρλγ (y, r; µ 2 F , Q 2 ) , (6.5) where N λρλγ are normalization factors N λρλγ = ∞ 0 dr r dy W λρλγ (y 1 , y, r; µ 2 F , Q 2 ) . (6.6) Expressed in terms of these functions, the scattering amplitudes read The probability amplitude for a dipole of size r to scatter on the nucleon is then proportional to P λρλγ (r, Q 2 , µ 2 F )σ(x, r), which justifies the inclusion of the factor r in (6.5). Below we will also use the rescaled radial distributions P λρλγ (λ, µ 2 F ), which depend on r and Q only through the variable λ = r Q and we choose to put µ 2 F = µ 2 F (Q 2 ), see eq. (5.1). Note that the rescaled asymptotic distributions,  This change of variable leads to the formulas The average value of a function f (y) depending of the longitudinal fraction of momentum y carried by one of the partons, will be estimated by f (y) λρλγ = 1 N λρλγ dr dy f (y) r W λρλγ (y, r; µ 2 F , Q 2 ) .   (r, Q 2 ) vs the size r of the interacting dipole, for Q 2 = 1 GeV 2 (12(a)) and Q 2 = 10 GeV 2 (12(b)), and the dipole cross-sectionσ(x, r) normalized by the factor 4σ 0 for W = 50 GeV and W = 150 GeV.

The radial distribution of the
The distributions P 00 (r, Q 2 , µ 2 F ) and P (AS) 00 (r, Q 2 ) ≡ P 00 (r, Q 2 , ∞) are close to each other, as it is shown in figures 12(a) and 12(b), which indicates that the distribution P 00 (r, Q 2 , µ 2 F ) is not sensitive to µ 2 F . We can then restrict the study of P 00 (r, Q 2 , µ 2 F ) by considering only P (AS) 00 (r, Q 2 ), which is also simpler for the analytic treatment. At first glance we see that the distributions are peaked around r ∼ 1.3 Q and consequently the peak moves to the right and the distribution becomes wider as Q 2 decreases. Note that the dependency ofσ(x, r) with respect to Q 2 , which can be seen in figure 12, only occurs through the dependency of R 0 (x), according to eq. (4.4).  In figure 13 we show the Total and AS rescaled radial distributions P 00 (λ, µ 2 F (Q 2 )) and P (AS) 00 (λ, µ 2 F (Q 2 )). This last one reads P (AS) 00 The average value of λ estimated with P (AS) 00 (6.14) About half of the dipoles are contained in the region 1 < λ < λ (AS) 00 , the peak of the distribution being at λ peak ∼ 1.3. The typical transverse scale µ = yȳ Q 2 entering the wave functions overlap can be estimated using eq. (6.12), The choice of the factorization scale µ F (Q 2 ) given by eq. (5.1) is then a good approximation of the transverse dynamical scale µ (AS) 00 involved in the process. The dipole scattering amplitude plays the role of a filter that selects dipoles of λ > λ Sat. (Q 2 , W ) = 2R 0 (x) Q. Note that the critical saturation line eq. (4.5) is given by in accordance with eq. 4.6. In the kinematics of HERA, the energy in the center of mass W varies roughly from 50 GeV to 150 GeV, leading to the following bounds for the two values Q 2 = 1 GeV 2 and Q 2 = 10 GeV 2 , as one can see in figure 13 (plotted in logarithmic scale). The large difference between N λ>λ Sat. (1,90) and N λ>λ Sat. (10,90) indicates that the integrand of T 00 shown in figure 14, is very sensitive, when Q 2 varies between 1 GeV 2 and 10 GeV 2 , to the overlapping of the dipole cross-section bandwidth (λ > λ Sat. (Q 2 , W )) and the radial dipole distribution P (AS) 00 (λ); we are then probing with a high accuracy the quality of the shape of the dipole cross-section.
At high Q 2 the tails of the distributions plays a dominant role. The tail of the distribution corresponds to the region where the integrand of the radial distribution can be approximated by an exponential fall, where typically (6.23) In the case Q 2 = 1 GeV 2 , the bandwidth of the dipole cross-section mostly overlaps with the peak of the distributions as λ Sat. (1, 90) ∼ λ peak , while it only overlaps with the tail of the distribution when Q 2 = 10 GeV 2 , λ Sat. (10, 90) ∼ λ tail . Figure 14 shows the normalized integrand of T 00 . This summarize our discussion on the respective role of the radial distribution and of the dipole cross-section. This integrand is peaked near the saturation radius r ∼ 2R 0 (x). Comparing the cases Q 2 = 10 GeV 2 and Q 2 = 1 GeV 2 , we see that this peak is moving to the right as Q 2 decreases, going through the bandwidth of the dipole cross-section.
• the WW contributionP W W
Contrary to the γ * L → ρ L transition, we see that the dependence on µ 2 F is quite strong if one compares P 11 (λ, µ 2 F ) to P (AS) 11 (λ). We cannot then restrict ourselves to only study the asymptotic case.
It is interesting to estimate the average λ obtained with the different distributions P 11 , P  which are normalized separately.
The results in table 4 show that the γ * T → ρ T transition is more sensitive to saturation effects than the γ * L → ρ L transition as λ 11 is about twice larger than λ 00 . Indeed it means that more dipoles are produced in the bandwidth of the dipole cross-section by the radial distribution P 11 (λ, µ 2 F ) than P 00 (λ, µ 2 F ). As a consequence the polarized cross-section σ T should be a more sensitive observable to probe features of the saturation regime than σ L .
The transverse dipole scale associated to the WW contribution, using eq. (6.12), is which is not so far from the values of the function µ 2 that is used here. The tail of the distributionP W W 11 (λ, µ 2 F ) can be defined as, is approximated by √ yȳ AS 11 ≈ 0.37. The genuine contribution P gen 11 (λ, µ 2 F ), which vanishes in the limit µ 2 F → ∞, depends strongly on the factorization scale µ 2 F (Q 2 ), as one can see in figure 15. The distribution P gen 11 (λ, µ 2 F ) can be split into two contributions, whereP gen (qq) 11 (λ, µ 2 F ) involves the exchange of two-parton in the hard part and the genuine solutions of the two-parton DAs, while P gen (qqg) 11 (λ, µ 2 F ) involves the three-parton exchange hard part F γ * T (y 1 , y 2 , r; Q) defined in eq. (3.19). The transverse scale µ present inP gen (qq) 11 (µ 2 F ) ∼ 0.5 and is not very sensitive to µ 2 F . The choice of µ 2 F (Q 2 ) ∼ Q 2 is then a good choice for this contribution.
The contributionP gen (qqg) 11 (λ, µ 2 F ) have several transverse scales which are µ 1 , µ 2 , µq g , µ qg and µ qq defined by eqs. (3.20), each corresponding to a dipole configuration involving two of the three partons available in the process. In order to estimate these transverse scales, we first evaluate the average fraction of momentum carried by the quark y 1 gen (qqg) 11 (µ 2 F ), the antiquark ȳ 2 gen (qqg) 11 (µ 2 F ) and the gluon y g gen (qqg) 11 .  (µ 2 F ) 41.8% 41.8% 16.3% (µ 2 F ) of the fraction of longitudinal momentum of the ρ meson carried by each of the parton i.
An estimation of the transverse scales can be made, using where we assume that the transverse scale values are roughly approximated by using the average fractions of longitudinal momentum in eqs. (3.20). These values are evaluated at µ 2 F (1 GeV 2 ). Other values of µ 2 F have also been used, leading to approximately the same results. We note that the function µ F (Q 2 ) ≈ Q 2 is close to the µ 1 and µ qq values, and of the same order of magnitude than µ qg but in principle, one should adapt the choice of the factorization scale to the relevant transverse scales at stake for each part of the process. The WW and the genuine contributions to the radial distribution P 11 (λ, µ 2 F ) are of the same order of magnitude for Q 2 ∼ 10 GeV 2 , and the genuine contribution becomes even more important for Q 2 = 1 GeV 2 . At Q 2 = 10 GeV 2 , as λ gen 11 (∼ 3.5) < λ Sat. (10 GeV 2 )(∼ 7.7) < λ W W

11
(∼ 8.7) , most of the dipoles in the bandwidth of the dipole cross-section are provided by the WW contribution, which explains why the predictions for σ T are dominated by the WW predictions, as shown in figure 4(a), and consequently, why the results depend weakly on the factorization scale. However the fact that the genuine contribution is important even at large Q 2 , indicates that the three-parton exchange between the γ * T and ρ T states is important in such relations as the normalization eq. (6.36) of the ρ meson wave function [7,81], and the electronic decay width eq. (6.37) [7,82], where the exchange of only two-parton is assumed. Indeed, the r.h.s. of eq. (6.37), if one expands at large Q 2 the ρ meson wave function around r = 0 , is the WW approximated result, which therefore misses the genuine contributions arising from three-parton correlators, which can have a significant effect even for large Q 2 values, see figure 15. These relations are usually used to constrain the parameters of the wave function models of the ρ meson, assuming that the meson is solely constituted of a quark and an antiquark, which then consist in neglecting the higher Fock state contributions like the genuine contribution. At Q 2 = 1 GeV 2 , the peak of the genuine part of the radial distribution enters the dipole cross-section bandwidth (λ Sat. (1, 90) < λ gen 11 ) and gives an important contribution to the integrand of T 11 . When increasing Q 2 , as shown in figure 16 where we display the product P 11 (r, Q 2 , µ 2 F )σ(x, r) for Q 2 = 1 GeV 2 and Q 2 = 10 GeV 2 , we can note that the difference between the AS and the Total results is a consequence of the genuine contribution growing when Q 2 decreases. It is the convolution with the dipole cross-section which washes-out the effect of these genuine twist-3 contributions.
6.3 Comparison with the radial distributions obtained from models of the ρ meson wave function.
It is instructive to compare shapes of the radial distributions P 00 and P 11 used in our analysis with those used in two other approaches which involve the overlap of the virtual photon wave functions and of the ρ meson wave functions: • the "Boosted Gaussian" (BG) model [46], • the "Gaus-LC" model [55].
We use here the convention and the parameter values of ref. [47], which for completeness are shown in table 6. The scalar parts of the wave functions are given by,  Table 6: Parameter of the "Gaus-LC" and the "Boosted Gaussian" models taken from ref. [47], for M ρ = 0.776 GeV, f ρ = 0.156 GeV, m f = 0.14 GeV and with f L ρ = f ρ .
φ Gauss-LC The overlaps with the virtual photon wave function are,   and AS (red, dashed) radial distributions for the γ * L → ρ L transition (top) and for the γ * T → ρ T transition (bottom), vs r for Q 2 = 1 GeV 2 (left) and Q 2 = 10 GeV 2 (right), as well as the dipole cross-sectionσ(x, r) rescaled by the factor 5σ 0 for W = 90 GeV (black, dot-dashed).
at Q 2 = 1 GeV 2 (figures 17(a) and 17(c)) and Q 2 = 10 GeV 2 (figures 17(b) and 17(d)) we see that at large Q 2 , in the bandwidth of the dipole cross-section, our distributions are converging with the distributions obtained from the Gaus-LC and BG models.
In the γ * T → ρ T case in figures 17(c) and 17(d), we see that the BG and Gaus-LC models are closer to the distribution P 11 than to the asymptotic distribution. When compared to other distributions, the asymptotic distribution P AS 11 is shifted to the right, thus selecting larger dipole sizes.
In the γ * L → ρ L case, in figure 17(a) we see that the distributions from the Gaus-LC and BG models are not close to our predictions, indicating that the higher twist corrections are presumably more important at small Q 2 than in the γ * T → ρ T transition. These qualitative remarks remains the same when using the AAMQS dipole models, and we expect that they are model independent.

Conclusions
We performed phenomenological analysis of experimental data from HERA on ρ meson electroproduction within the approach based on the recently derived [44] impact parameter representation of γ * → ρ impact factor up to twist-3 accuracy. The important feature of this representation consists in the inclusion of contributions coming from both two-and three-partonic Fock states, maintaining a close connection with the dipole model picture. Consequently, it was possible to include in our framework the saturation effects. Our predictions show that we can get simultaneously good predictions for the polarized crosssections σ T and σ L . The abilility of the model to reproduce the data is the confirmation of the following points: • the factorization of the dipole cross-section in the helicity amplitudes of the electroproduction of the ρ meson works and, as an universal quantity, is the same for T 00 and T 11 , giving the good energy dependence and normalizations of the polarized cross-sections, • the collinear factorization procedure of the ρ meson is justified and works successfully beyond the leading twist.
As expected the model has some limits due to the truncation of the twist expansion. Thanks to HERA data, we have identified the virtuality Q 2 min ∼ 5 GeV 2 were the higher twist corrections become important, which is a motivation to compute impact factors beyond the twist 3 accuracy in order to probe the genuine saturation regime which starts at Q 2 S ∼ 1 GeV 2 .
Other helicity amplitudes could be computed keeping the same approach, they would be useful in the t = t min regime. The kinematics of the impact factor can be also extended to take into account the t−dependence of the impact factors, which would be a test for the dipole models which include the impact parameter dependence, providing a probe of the proton shape [7], in particular through local geometrical scaling [83,84].
Data also exist for φ leptoproduction. In this case quark-mass effects should be taken into account, in particular, because this allows the transversely polarized φ to couple through its chiral-odd twist-2 DA. Indeed, as it was pointed out in [42], the fact that the ratio T 11 /T 00 is not the same (after trivial mass rescaling) for ρ and φ mesons points to the importance of this effect. This is also beyond the scope of our present study, but may open an interesting way for accessing chiral-odd DAs.
The next-to-leading order effects -both on the evolution and on the impact factorshould be studied, since it is now known that both may have a important phenomenological effect [85][86][87][88][89].

A.1 Distribution amplitudes in the LCCF parametrization
The seven chiral-even 7 ρ-meson DAs up to twist 3 are defined through 8 the following matrix elements of nonlocal light-cone operators [40], for two-parton and three-parton correlators where we used the standard notation . DAs are linked by linear differential relations derived from equations of motion and n−independency condition [39,40]. The solutions for ϕ P (y) ≡ {ϕ 3 , ϕ A , ϕ T 1 , ϕ T A } are the sum of the solutions in the so-called WW approximation and of genuine solutions, ϕ P (y) = ϕ W W P (y) + ϕ gen P (y) . (A.7) 7 The chiral-odd twist-2 DA for the transversely polarized ρ meson does not contribute to the process considered at the accuracy discussed here. This is also true in the approach based on collinear factorization of generalized parton distributions [100,101]. 8 In the approximation where the mass of the quarks is neglected with respect to the mass of the ρ meson.
The WW approximation consists in neglecting the contribution from three-parton operators, thus taking B(y 1 , y 2 ; µ 2 F ) = D(y 1 , y 2 ; µ 2 F ) = 0. Then, ϕ W W P (y) become functions of ϕ 1 (y) only, and their explicit expressions are given by Genuine solutions only depend on {B(y 1 , y 2 ; µ 2 F ), D(y 1 , y 2 ; µ 2 F )} or equivalently on the combinations {S(y 1 , y 2 ; µ 2 F ), M (y 1 , y 2 ; µ 2 F )} defined by eq. (3.3), namely where A(u; µ 2 F ) has the compact form and it obeys the conditions The correspondence between our set of DAs and the one defined in ref. [50] is achieved through the following dictionary derived in ref. [40]. It reads, for the two-parton vector DAs, and for the axial DA, For the three-parton DAs, the identification is Explicit forms for ϕ 1 , B, and D are obtained with the help of the results of ref. [50] obtained within the QCD sum rules approach. The first terms of the expansion in the momentum fractions of the three independent DAs thus have the form The dependences on the renormalization scale µ F of the coupling constants a 2 , ω A {1,0} , ζ A 3 , and ζ V 3 are given in ref. [50]. In appendix A.2 we present both the evolution equations and the values of these constants at µ 2 F = 1 GeV 2 used in our analysis, as well as the dependence on µ F of the DAs.

A.2 Evolutions of DAs and coupling constants with the renormalization scale
The parameters entering the DAs at µ 2 0 = 1 GeV 2 are updated 9 in ref. [51] and their evolution equations are given in ref. [50], we recall in table 7 their values for the ρ meson.

A.4 Results using the GBW and AAMQSb models
We present some of the predictions obtained by using the GBW or the AAMQSb model for the dipole cross-section. As expected the results are not so far from the ones obtained with the AAMQSa model. In figure 18 and 19 are respectively shown for the AAMQSb and the GBW model, the polarized cross-sections σ T and σ T . The spin density matrix element r 04 00 predictions using these dipole models are shown in figure 20, for completeness we show also the prediction obtained with the GS model.