Joint resummation for pion wave function and pion transition form factor

We construct an evolution equation for the pion wave function in the $k_T$ factorization theorem, whose solution sums the mixed logarithm $\ln x\ln k_T$ to all orders, with $x$ ($k_T$) being a parton momentum fraction (transverse momentum). This joint resummation induces strong suppression of the pion wave function in the small $x$ and large $b$ regions, $b$ being the impact parameter conjugate to $k_T$, and improves the applicability of perturbative QCD to hard exclusive processes. The above effect is similar to those from the conventional threshold resummation for the double logarithm $\ln^2 x$ and the conventional $k_T$ resummation for $\ln^2 k_T$. Combining the evolution equation for the hard kernel, we are able to organize all large logarithms in the $\gamma^{\ast} \pi^{0} \to \gamma$ scattering, and to establish a scheme-independent $k_T$ factorization formula. It will be shown that the significance of next-to-leading-order contributions and saturation behaviors of this process at high energy differ from those under the conventional resummations. It implies that QCD logarithmic corrections to a process must be handled appropriately, before its data are used to extract a hadron wave function. Our predictions for the involved pion transition form factor, derived under the joint resummation and the input of a non-asymptotic pion wave function with the second Gegenbauer moment $a_2=0.05$, match reasonably well the CLEO, BaBar, and Belle data.


INTRODUCTION
Great efforts have been devoted to the extension of the k T factorization theorem for exclusive processes [1,2,3,4,5,6] to subleading levels recently. The next-to-leading-order (NLO) corrections to the pion transition (electromagnetic) form factor associated with the πγ * → γ(π) scattering have been calculated at leading power [7,8]. Those to the B → π transition form factors involved in B meson semileptonic decays were derived in [9]. Up to subleading power, the three-parton contributions to the pion electromagnetic form factor, to the B → γ transition form factor, and to the B → π transition form factors have been studied in [10], [11], and [12], respectively. A k T -dependent hard kernel is defined as the difference between QCD diagrams and effective diagrams for transverse-momentumdependent (TMD) hadron wave functions. Therefore, to obtain a NLO hard kernel, both QCD diagrams and effective diagrams need to be evaluated at the same level. The NLO analysis of the B meson and pion wave functions have revealed various important logarithms, which stimulate corresponding resummation formalisms for their organization to all orders in the strong coupling constant.
A TMD hadron wave function contains the light-cone singularity from the region with a loop momentum parallel to a Wilson line on the light cone [13]. To regularize the light-cone singularity, one may rotate the Wilson line away from the light cone to an arbitrary direction u with u 2 = 0 [13,14]. The higher-order wave function then depends on u 2 through the scale ζ 2 P ≡ 4(P · u) 2 /u 2 , where P denotes the hadron momentum. The variation of u, namely, of ζ 2 P introduces a factorization-scheme dependence into the hadron wave function. The evaluation of the NLO effective diagrams for the B meson wave function indicates the existence of the logarithms ln 2 (ζ 2 P /m 2 B ) and ln x ln(ζ 2 P /m 2 B ) [9], m B being the B meson mass and x being the momentum fraction of the spectator. The NLO diagrams for the pion wave function produce the mixed logarithm ln x ln(ζ 2 P /k 2 T ) [8], k T being the parton transverse momentum. All the above logarithms become large as u 2 → 0, and as x and k T are small, which is the dominant kinematic region in the k T factorization theorem for exclusive processes. The logarithms in the B meson wave function have been organized under the rapidity resummation [15], whose effect was shown to diminish the B meson wave function at the end point x → 0.
The above observation hints that the resummation of the mixed logarithm ln x ln(ζ 2 P /k 2 T ) for the pion wave function would modify both the x and k T dependencies. It then calls for the joint resummation [16,17,18,19], which was proposed to unify the conventional threshold resummation for ln 2 x [20,21,22] and the conventional k T resummation for ln 2 k T [14,23]. For a recent review on this subject, see [24]. In this paper we will construct an evolution equation in the scale ζ 2 P following the idea in [17], whose solution resums the mixed logarithm in the Mellin (N , conjugate to x) and impact-parameter (b, conjugate to k T ) spaces. The inverse Mellin transformation is then applied to get the x dependence of the pion wave function. It will be demonstrated that the joint resummation induces suppression which is stronger at small x than at moderate x, and intensifies with increase of b. This effect, similar to those of the threshold and k T resummations, improves the applicability of perturbative QCD (PQCD) to hard exclusive processes. Combining the evolution equation for the hard kernel of the γ * π 0 → γ scattering, we organize all the relevant large logarithms, and remove the factorization-scheme dependence on ζ 2 P mentioned before. This is the first time that the k T factorization for a simple exclusive process can be made scheme independent in the presence of the light-cone singularity.
It has been known that γ * π 0 → γ serves as an ideal process for the determination of the pion wave function, and the involved pion transition form factor F (Q 2 ), Q 2 being the momentum transfer squared, has been investigated thoroughly. In particular, it was claimed that the quantity Q 2 F (Q 2 ) (including those for the η and η ′ meson transition form factors) begins to saturate at relatively low Q 2 as calculated in the hard scattering approach [25], in QCD sum rules (QCDSR) [26], in light-cone sum rules (LCSR) [27], and in the light-front holographic QCD [28]. We will analyze the leading-order (LO) and NLO contributions to the pion transition form factor with inputs of different model wave functions, including the asymptotic model, the flat model, and the model with the second Gegenbauer moment a 2 . The results are compared to those from the PQCD approach [29], that incorporates the conventional threshold and k T resummations. It will be observed that the significance of the NLO correction to and the saturation behavior of Q 2 F (Q 2 ) differ under the joint resummation and the conventional resummations. It implies that QCD logarithmic corrections to a process must be handled appropriately, before its data are used to extract a hadron wave function. Our predictions for Q 2 F (Q 2 ) from a non-asymptotic pion wave function with a 2 = 0.05 match reasonably well the CLEO, BaBar, and Belle data, which seem to indicate scaling violation at currently accessible Q 2 .
In Sec. 2 we construct the evolution equation for the resummation of the mixed logarithm in the pion wave function, and then solve it in the Mellin and impact-parameter spaces. The inverse Mellin transformation of the solution is performed in Sec. 3, with different initial conditions of the evolution. Note that the running of the strong coupling constant down to the low energy region has to be modified in order to avoid the Landau pole. The joint resummation effect on the x and b dependencies of the pion wave function is then examined. In Sec. 4 the pion transition form factor is evaluated for a given model wave function at the LO and NLO levels under the joint resummation and the conventional resummations. The different outcomes for the NLO contributions and for the saturation behaviors at high energy are compared. We summarize our findings, and discuss potential extension of our formalism to more complicated processes in Sec. 5. The explicit expressions for the solutions of the evolution equation are collected in Appendix A.

EVOLUTION EQUATION
The TMD pion wave function Φ(x, k T ) is defined by the non-local hadron-to-vacuum matrix element 1 where µ f is the factorization scale, the coordinate y = (y + , 0, y T ) is off the light cone generally, and xP − and k T are the longitudinal and transverse momenta carried by the anti-quarkq, respectively. A TMD hadron wave function describes the distributions of a light parton in both light-ray and transverse directions. To maintain the gauge invariance of the definition in Eq. (2.1), the gauge-link operator W y (u) has been introduced, where g is the QCD coupling constant, and P denotes the pathordered exponential. The non-light-like vector u, different from the usual Wilson line direction n + = (1, 0, 0 T ), plays a role of the regulator for the light-cone divergences [13]. The transverse gauge link I u;y,0 , unraveling the cusp obstruction in the contour of the Wilson lines at infinity, does not contribute in the covariant gauge [31]. As determining a NLO hard kernel in the k T factorization formula for a pion-induced process, we perform the infrared substraction defined as the convolution of the NLO pion wave function with the LO hard kernel. The QCD correction to the pion wave function gives rise to the mixed logarithm ln x ln(ζ 2 P −2 /k 2 T ) [7,8,9], with the dimensionless rapidity parameter

3)
n − = (0, 1, 0 T ) being a light-like vector along the pion momentum P . The double rapidity logarithm ln 2 ζ 2 in the B meson case is absent here because of the color-transparency mechanism for an energetic pion, which suppresses soft gluon contributions. The goal of this section is to construct an evolution equation, whose solution sums the mixed logarithm in the pion wave function. Following [17], we trade the derivative with respect to the rapidity parameter ζ 2 for the variation of the Wilson link direction u, It is obvious that this chain rule simplifies the analysis dramatically as the u dependence appears only through the Wilson line interactions. Applying Eq. (2.4) to the Feynman rule associated with the Wilson link, we have with the special vertexû We will derive the rapidity evolution equation where ⊗ represents convolutions in the momentum fraction x and the transverse momentum k T , and the evolution kernel Γ involves the diagrams with the special vertex.

Evolution kernel
It is easy to see that the structure of the special vertex suppresses a collinear gluon contribution to Γ [15]. The evolution kernel is then dominated by soft and hard gluon exchanges, usually denoted as the functions K and G, respectively. The soft and hard gluon radiations off the active quark, as shown in Fig. 1, lead to , (2.8) for the function K, and , for the function G with the variablex ≡ 1 − x.  Adopting the dimensional regularization for the ultraviolet divergence and regularizing the infrared divergence with the gluon mass λ, we obtain Since the soft divergences cancel between K 1 and K 2 , and between G 1 and G 2 , the gluon mass λ will approach to zero eventually. For the evaluation of K 2 , we apply the Mellin and Fourier transformations b being the impact parameter. Equation (2.9) then givesK 2Φ (N, b, ζ 2 , µ f ) with the soft kernel (2.13) in which the terms suppressed by powers of 1/ζ 2 have been dropped, and K 0 is the zerothorder modified Bessel function of the second kind. Hence, the bare functionK (b) is written asK where the large-N expansion of Eq. (2.13) has been made, and the superscript (b) labels the bare function explicitly. The bare hard function G (b) can be calculated following the same line, and reads (2.15) We will adopt the approximationx ≈ 1 in the small x region, where the mixed logarithm plays a significant role. It is found that both the soft and hard functions depend on the factorization scale µ, and such a dependence cancels in their sum. This fact is attributed to the µ independence of the mixed logarithm that we are going to resum. Applying the modified minimal substraction (MS) scheme to the ultraviolet renormalization yieldsK where the additive counterterms δK (δG) of the functionK (b) (G (b) ) can be read from Eq. (2.14) (Eq. (2.15)). The renormalization-group (RG) equations for the soft and hard functions are then given, in terms of the anomalous dimensions λK and λ G , by which lead to the RG improved evolution kernel We choose the scales to diminish the initial conditionsK (r) (µ 0 ) and G (r) (µ 1 ). The evolution kernel Γ also contains the diagrams with gluon radiations from the spectator quark in principle. However, these diagrams contribute at the next-to-leading logarithmic level, because the NLO effective diagrams with gluon radiations off the spectator quark do not generate the mixed logarithm ln x ln(ζ 2 P −2 /k 2 T ) as indicated by Eqs. (36) and (37) in [9]. The corresponding soft and hard functions are expressed as which generate (2.21) The logarithm ln x in the soft function K ′ 2 ⊗ Φ can be extracted by implementing the approximation [16] Φ under which we obtain The cancelation of the soft divergences between K ′ 1 (= K 1 ) in Eq. (2.11) and K ′ 2 in Eq. (2.23) is evident, whose sum gives (2.24) Applying the MS scheme to the bare soft and hard functions gives the renormalized ones (2.26) Obviously, K ′ (r) and G ′ (r) are characterized by the same scale xζP − , implying that a RG treatment is not necessary here, and that the sum K ′ (r) + G ′ (r) produces only a nextto-leading logarithm as stated above. Hence, this contribution can be absorbed into the solution of the evolution equation by tuning the initial rapidity parameter ζ, whose variation within the order-unity range causes a next-to-leading logarithmic effect. We will take advantage of the freedom in choosing the bounds of ζ to achieve the matching between the resummation formula and the NLO results of the pion transition form factor. That is, the summation of the above next-to-leading logarithms can be taken care of by the matching procedure, and the kernel K ′ (r) + G ′ (r) will be neglected below.

Equation (2.7) under the Mellin and Fourier transformations becomes
with the evolution kernel (2.28) Solving the differential equation (2.27), we get which constitutes one of the main technical results of this paper. The initial rapidity parameter ζ 0 will be specified later, and the step function in the exponent will become effective as we perform the inverse Mellin transformation. Apart from the mixed logarithm, the NLO pion wave function contains the single logarithm ln(µ f /Q), which can be summed via the standard RG equation with the anomalous dimension [8] γ Combining the joint resummation and the solution to Eq. (2.30) leads tõ where µ i is the initial scale of the RG evolution. Note that the physical form factor is independent of the factorization scheme and the factorization scale µ f , whereH represents the hard kernel in the Mellin and impact-parameter spaces. Therefore, we have the evolution equation for the joint resummation, and the RG equation The solution of the above two differential equations gives the resummation improved hard kernelH with the final rapidity parameter ζ 1 and the characteristic hard scale t. We make use of the freedom of choosing the bounds ζ 2 0 and ζ 2 1 for the joint resummation, such that the NLO logarithmic enhancements inΦ(N, b, ζ 2 0 , µ i ) andH(N, b, ζ 2 1 , Q 2 , t), shown in Eqs. (39) and (40) of [7], respectively, are eliminated. This requires with the constants which recapitulates the joint-resummation improved k T factorization formula. The exponential factor in Eq. (2.39) describes the evolution from the initial conditionΦ(N, b, ζ 2 0 , µ i ) to the resummation improved wave functionΦ(N, b, ζ 2 1 , t). We have confirmed that the expansion of the exponential factor up to O(α s ) reproduces the mixed logarithm and the single logarithm ln(1/N ) in the NLO pion transition form factor [7]. Note that our resummation formalism was established in the conjugate space, while the calculation in [7] was performed in the momentum space. Hence, the correspondence between ln(1/N ) in the former and ln x in the latter is not precise, and the matching condition confirmed above in fact suffers order-unity uncertainty at the next-to-leading-logarithmic level.
At last, we point out that the double logarithm ln 2 (Q 2 /k 2 T ) in the TMD wave function and ln 2 x in the hard kernel were resumed in the conventional PQCD approach [7,32]. Besides, the rapidity parameter ζ 2 was fixed to a specific value for convenience. Compared to the joint resummation, the double logarithm ln 2 x in the TMD wave function has been ignored (see Eq. (37) in [7] for its existence), and the PQCD approach is not factorizationscheme independent, strictly speaking. That is, the formalism presented in this work represents a complete treatment of the logarithmic enhancement in the pion transition form factor, and the first scheme-independent k T factorization formula.

RESUMMATION IMPROVED WAVE FUNCTIONS
In this section we explore the detailed properties of the resummation improved wave func-tionΦ(N, b, ζ 2 1 , t). The factorization theorem for hard exclusive processes is usually formulated in the momentum-fraction space (see, however [33]). The inverse Mellin transformation forΦ(N, b, ζ 2 1 , t) gives where the parameter c is an arbitrary real number larger than the real part of the rightmost singularity ofΦ(N, b, ζ 2 1 , t) in the complex N plane, and the Cauchy theorem can be applied to deform the integration contour whenever necessary. We will not implement the inverse Fourier transformation, so that the joint-resummation effect can be compared with the Sudakov-resummation effect directly, which is usually studied in the impact-parameter space.
The parametrization of the TMD pion wave function has been extensively discussed in the literature (for a recent discussion, see [34]). For simplicity, factorization of the initial pion wave function in the longitudinal and transverse momentum spaces will be postulated. Keep in mind that the major task of this section is to illustrate the joint-resummation effect. For definiteness, the transverse momentum distribution is taken as where the prefactor is introduced to obey the normalization 4) and the shape parameter β is related to the root mean square of the transverse momentum via According to [35,36], the input k 2 T 1/2 = 350 MeV that fulfills various constraints (including the π → γ γ decay rate) leads to β = 2.0 GeV −1 .
The longitudinal momentum distribution φ(x, ζ 2 0 , µ i ) is assumed to be the same as the light-cone distribution amplitude (LCDA) ϕ(x, µ i ). The one-loop evolution equation indicates that the pion LCDA can be expanded in terms of the Gegenbauer polynomials C where the odd Gegenbauer moments a 2n+1 vanish due to symmetry prosperities. The dependence of a 2n on the scale µ i is governed by the well-known Efremov-Radyushkin-Brodsky-Lepage equation [37,38]. Along this line, we consider the following three models for the longitudinal momentum distribution with the Gegenbauer polynomial C 3/2 2 (x) = (3/2)(5x 2 − 1). The first model φ I corresponds to the pion LCDA in the asymptotic limit. The flat distribution φ II was proposed in [39,40], where a nonperturbative correction beyond the operator product expansion was also introduced to explain the scaling violation indicated by the BaBar data. As there is overwhelming evidence that the pion LCDA at energy scales accessible in current experiments is broader than the asymptotic model, we keep the sub-leading Gegenbauer term in φ III . The contribution from a higher Gegenbauer term to the pion transition form factor depends on the momentum transfer squared Q 2 and the shape parameter β [41]. Fitting to the Babar data, it has been realized that one can at best determine the second Gegenbauer moment and the shape parameter simultaneously in the framework of the k T factorization [34]. For this reason, also expecting the quick convergence of the conformal spin expansion of the pion wave function (see, however, [42]), we will confine the analysis to the second Gegenbauer moment. Our formalism can be extended to include higher Gegenbauer terms straightforwardly.

Resummation with fixed α s
To make our discussion more transparent, we start from the inverse Mellin transformation with a frozen coupling constant, and then generalize it to the case with a running coupling constant. For a frozen coupling α s , the joint-resummation improved wave functioñ The exponential contains a branch cut on the negative real N axis and a singularity at N = 0. The analytical property of the wave functionΦ(N, b, ζ 2 1 , t) also depends on the initial conditionΦ(N, b, ζ 2 0 , µ i ). The Mellin and Fourier transformations defined in Eq. (2.12) lead the three models tõ As inserting the asymptotic modelΦ I (N, b, ζ 2 0 , µ i ) into Eq. (3.8),Φ(N, b, ζ 2 1 , t) develops two additional poles at N = −1 and −2. According to the Cauchy theorem, the contour for the inverse Mellin transformation is deformed as displayed in Fig. 2, which (i) runs from minus infinity towards N = −2 below the branching cut, (ii) slides into an infinitesimal semicircle around N = −2, (iii) continues toward N = −1 below the branching cut, (iv) slides into another infinitesimal semicircle around N = −1, (v) runs to N = −r with 0 < r < 1, (vi) revolves around the origin along a finite circle of radius r, and (vii) runs back to minus infinity in a way that reverses the steps (i)-(v) above the branching cut. It is trivial to derive the joint-resummation improved pion wave function The first term in the above expression comes from the contributions of the N = −1 and −2 poles, the second term corresponds to the integration along the circle at N = 0 with radius r, and the last term arises from the discontinuity of the integrand along the branching cut. It can be verified that Φ I (x, b, ζ 2 1 , t) is independent of the radius 0 < r < 1 as it should. The inverse Mellin transformation is performed along the same line in the case of the flat model. The initial condition introduces only a single pole at N = 0, which is also a branching point of the exponential in Eq. (3.8). The contour is shown in Fig. 2, with the two infinitesimal circles around N = −1 and = 2 being removed. This immediately yields (3.14) The resummation improved wave function with the initial conditionΦ III (N, b, ζ 2 0 , µ i ) is calculated similarly, albeit with more involved analytical structures of the integrand; we need to modify the contour in Fig. 2, so that two additional poles at N = −3 and −4 are circumvented. It implies that including higher Gegenbauer terms in the initial condition of the wave function generates a longer sequence of poles to be avoided in the contour integration. We derive where the coefficients κ n and the function f (x 1 , x 2 ) are given by Apparently, Eq. (3.15) reduces to Eq. (3.12), when the second Gegenbauer moment a 2 is set to zero. A more complicated model for the pion wave function with 0 < α < 1, has been advocated in [43]. Specifically, the model with α = 1/2 was derived in the light-front holographic QCD approach [44]. The corresponding initial conditionΦ develops infinitely many poles at N + α = 0, −1, −2, ... in the complex N plane. The study of the joint-resummation effect on this model is similar, and will not be performed here.

Resummation with running α s
We are now in a position of computing the joint-resummation improved pion wave function Φ(x, b, ζ 2 1 , t) with a running coupling α s . To avoid the Landau singularity in the inverse Mellin transformation, the parametrization [45] α s (µ) = 4π β 0 is adopted at one loop, with the QCD scale Λ QCD and the one-loop QCD β function β 0 = (11N c − 2N f )/3, N c and N f being the numbers of colors and flavors, respectively. In the above expression the first term preserves the ultraviolet behavior of the standard QCD coupling, and the second term cancels the ghost pole at µ = Λ QCD . The substitution of Eq. (3.19) into Eq. (2.32) produces where the functions A 1 and C 1 are written as The exponential in Eq. (3.20) still contains a branching cut along the negative real N axis, so the contour in Eq. (3.1) is deformed in the way exactly the same as in the case with a frozen coupling. For the asymptotic pion wave function, we get where the explicit expressions of the functions F i (λ 1 , λ 2 , λ 3 , λ 4 , η) are collected in Appendix A, and the discontinuity of the polylogarithm function has been inserted. It has been also verified that the r dependence of Φ I (x, b, ζ 2 1 , t) cancels between the last two terms for arbitrary r in the range 0 < r < 1.
The same procedure leads to the joint-resummation improved pion wave function for the flat model, and for the non-asymptotic model. The joint-resmmation effect on the pion wave function with a frozen coupling, for the example set parameters α s = 0.3, b = 2ã P − /a and 4ã P − /a, and Q 2 = 5 GeV 2 is displayed in Fig. 3, where the RG evolution effect is not included. For all the three considered models, the modification appears as the impact parameter b is greater than b min =ã/(a P − ), which is easily understood from the exponentiation of the mixed logarithm − ln (ã P − b/a) (ln N + 2) in Eq. (3.8). If the rapidity and factorization-scale evolutions are switched off, it is confirmed that the pion wave function obeys its normalization. A crucial consequence of the joint resummation, as read from Fig. 3, is that the small x region receives stronger suppression compared to the moderate x region as expected, while the large x region almost remains intact. Moreover, the suppression strengthens with the transverse separation b between the valence quarks at a given longitudinal momentum fraction. Therefore, the joint-resummation effect does not allow a significant contribution from soft gluon exchanges. This well known Sudakov mechanism, first formulated in QED, improves the applicability of PQCD to hard exclusive processes.
Turning to the case with a running coupling, we adopt the parameter Λ QCD = 250 MeV and the flavor number N f = 6. As observed from Fig. 3, the resummation improved pion wave function takes on a behavior rather similar to that for a frozen coupling. A minor difference is that the small x region is even more suppressed in the former case, which further boosts our confidence on the applicability of PQCD to exclusive processes at moderate momentum transfer.
Before closing this section, we highlight the distinction between the pion wave functions including the Sudakov resummation and including the joint resummation. For simplicity, we confine ourselves to the asymptotic model, because the other two models exhibit a similar b dependence. As indicated in Fig. 4, both resummation formalisms lead to suppression on the wave function in the large b region, which becomes more significant as the momentum transfer Q increases. This observation fulfills the concept that an energetic pion is a compact hadronic bound state. A striking feature is that the joint-resummation improved wave function concentrates on the small b region more than the Sudakov-resummation improved one, which falls off smoothly b region is attributed to the exponentiation of the mixed logarithm − ln (ã P − b/a) ln N , where the suppression with b is magnified by the large coefficient ln N . However, the Sudakov-resummation improved wave function vanishes quickly as b → 1/Λ QCD , since the running coupling α s (1/b) hits the Landau pole. In the present derivation the Landau pole has been avoided as shown in Eq. (3.19), such that the joint resummation does not diminish the wave function as b → 1/Λ QCD . We emphasize that this distinction, due to the different treatments of the Landau-pole contribution, is not physically crucial. It is found that the large b region is more suppressed with the growing of x in Figs. 4(b), but not in Fig. 4(d): the small-x approximation has been adopted in the joint resummation, so its effect is insensitive to the variation of x. The phenomenological consequences on the pion transition form factor from the two resummations will be elaborated in the next section.

PION TRANSITION FORM FACTOR
The pion transition form factor F (Q 2 ) involved in the γ * π 0 → γ process is defined by the following matrix element γ(P ′ , ǫ * )|j em µ (q)|π 0 (P ) = i g 2 em ε µναβ ǫ * ν P α P ′β F (Q 2 ) , where j em µ (q) is an electromagnetic current, q = P ′ − P is the momentum transfer, and ǫ denotes the polarization of the outgoing photon. The form factor F (Q 2 ) (Q 2 = −q 2 ) was written, in the collinear factorization, as [38] with the NLO hard kernel [47,48,49] It is seen that F (Q 2 ) scales as 1/Q 2 from the power counting of the hard kernel, and is determined by the inverse moment of the pion LCDA at LO.

k T factorization formula
To suppress the end-point contribution (soft gluon exchanges) from the small x region in the collinear factorization, the k T factorization has been developed for hard exclusive processes, and continually refined by including the resummations of important logarithms and power corrections as stated in the Introduction. This more sophisticated factorization theorem can be derived diagrammatically [50] by applying the eikonal approximation to collinear particles and the Ward identity to the diagram summation in the leading infrared regions. For the rapidity parameter ζ 2 = 2, the k T factorization formula at leading power of 1/Q 2 under the conventional resummations was given by [7,29] The Sudakov factor S(x, b, Q, t) sums the double logarithm ln 2 (k 2 T /Q 2 ) and the single logarithm ln(t 2 /Q 2 ) through the RG equation. The threshold factor from the resummation of ln 2 x has been parameterized as for convenience, in which the power c(Q 2 ) was determined to be c(Q 2 ) = 0.04Q 2 − 0.51Q + 1.87 , (4.6) by fitting to the exact threshold resummation formula in the Mellin space. It was then observed that the nontrivial Q 2 dependence of c(Q 2 ) is important for accommodating both low and high Q 2 data from BaBar. Note that the self interactions of the Wilson links have been included into the NLO hard kernel in Eq. (4.4), such that the coefficient of the first term in the brackets has been changed from "1" to "3", compared to Eq. (40) in [7]. As argued in [51], the additional contribution from these self interactions can be canceled by the soft subtraction in an alternative definition for the TMD pion wave function, as the involved gauge parameter is tuned appropriately. In this work we have adopted the definition of the TMD pion wave function with off light-cone Wilson links.
To minimize the factorization-scheme dependence of the pion transition form factor, the resummation of the mixed logarithms in both the pion wave function and the hard kernel has been performed in the previous section. The large logarithms in the initial conditions of the pion wave function and the hard kernel were eliminated by choosing the bounds ζ 2 0 and ζ 2 1 in Eqs. (2.32) and (2.36), leading to (4.7) We then arrive at the joint-resummation improved factorization formula for the pion transition form factor

Numerical analysis
The first issue in the numerical analysis concerns the choice of the hard characteristic scale t. One choice would be t 2 = √ xQ/b that removes the remaining logarithm in Eq. (4.8). Another choice characterizing the typical quantum fluctuation of hard scattering is t = max( √ xQ, 1/b) as widely adopted in the PQCD approach [5]. In both scenarios the hard scale runs into the nonperturbative region at small x and large b, but take different values. We have confirmed that the two scenarios do not generate practical difference in our formalism for the pion transition form factor. It implies that the joint resummation has suppressed the contribution from the nonperturbative region effectively. Below we will adopt the second scenario as the default choice. Another important issue is the determination of the Gegenbauer moment a 2 . QCDSR calculations of moments of the pion wave function can be traced back to 1980s, pioneered by Chernyak and Zhitnitsky [52], where a rather high value a 2 (µ) ∼ 0.58 at the scale µ 2 ∈ [1, 2] GeV 2 was derived. This estimate was improved gradually by including NLO QCD corrections and refining the "internal" parameters of the QCDSR approach. The most recent update gave a 2 (1 GeV) = 0.15 ± 0.03 [53]. Following the strategy of LCSR, we will not use the Gegenbauer coefficient a 2 computed from QCDSR directly. Instead, we employ the value [46] a 2 (1 GeV) = 0.17 ± 0.08 , (4.9) extracted from matching the LCSR evaluation of the pion electromagnetic form factor, which includes NLO twist-2 corrections and higher power terms up to twist 6, to the experimental data from the Jefferson Lab Collaboration.
We start our numerical analysis with the input of the asymptotic pion wave function, choosing the initial scale µ i = 1/b for the RG evolution in Eq. (2.32) [5]. As observed from Fig. 5(a), the predicted Q 2 F (Q 2 ) with the conventional resummations at both LO and NLO levels saturates rapidly as Q 2 > 5 GeV 2 , and the NLO QCD correction enhances the form factor by (6 − 14)%. It is clear that the asymptotic pion distribution generally accommodates the Belle data except the first two bins. However, it cannot describe the CLEO and BaBar data in both small and large Q 2 region. Note that the incompatibility between the BaBar and Belle data on the pion transition form factor has been elaborated quantitatively in [57]. Besides, the impact at low energy of the BaBar and Belle high-energy data was analyzed by means of the Padé approximation inspired parametrization in [58]. The joint-resummation effect decreases the LO and NLO predictions in the conventional approach by (11 − 16)% and (8 − 27)%, respectively. Such decrease can be understood via the stronger reduction at small x from the joint resummation as shown in Figs. 4(b) and 4(d), which is the dominant region owing to the hard kernel K 0 ( √ xQ b). It is found that the saturation behavior of Q 2 F (Q 2 ) changes slightly at NLO in the joint-resummation improved k T factorization: the NLO correction brings about 6% suppression (15% enhancement) to the LO result in the small (large) Q 2 region. The above decrease of the NLO form factor at small Q 2 is explained as follows. The contribution to the pion transition form factor under the joint resummation mainly comes from the small b region as indicated in Fig. 4. We then have t 2 b ∼ 1/b > √ xQ at small Q, for which the logarithm of the NLO hard kernel in Eq. (4.8) flips sign. The failure of describing the experimental data suggests that the pion wave function might be broader. This observation is in agreement with the particular feature of the pion as a Nambu-Goldstone boson of dynamical chiral symmetry breaking and with the recent lattice simulations [43].
The computed pion transition form factor Q 2 F (Q 2 ) with the flat pion wave function is exhibited in Fig. 5(b). It is seen that the form factor grows steadily with Q 2 at LO and NLO in both resummation formalisms. This is easily realized from the scaling Q 2 F (Q 2 ) ∼ ln(Q 2 /k 2 T ) implied by the tree-level k T factorization formula with the flat pion wave function [29]. The LO curve from the conventional approach reasonably describes the scaling violation at large Q 2 observed by BaBar and the low Q 2 data from CLEO and Belle. The NLO correction increases the form factor by approximately (15 − 18)%, such that the agreement with the data deteriorates a bit. Compared to the conventional approach, the predictions from the joint resummation brings about 17% enhancement and (8 − 16)% suppression at the LO and NLO levels, respectively. The enhancement at LO, opposite to what was observed in Fig. 5(a) for the asymptotic model, arises from the weaker suppression than the parameterized threshold factor in Eq. (4.5) at small x. The NLO correction becomes destructive under the joint resummation in the whole range of Q 2 , decreasing the LO result by (10 − 15)%. This behavior, different from that in the case of the asymptotic model, is also traced back to the logarithmic term in Eq. (4.8): ln(t 2 b/( √ xQ)) remains positive in the small x region, which is probed more by the flat pion wave function. It is interesting to notice that the NLO curves from the conventional resummations and from the joint resummation turn out to be similar. The input of the third model of the pion wave function with a nonvanishing second Gegenbauer moment a 2 leads to Q 2 F (Q 2 ) displayed in Fig. 5(c). The form factor under the conventional resummations behaves in a way similar to that in Fig. 5(a): Q 2 F (Q 2 ) saturates as Q 2 > 10 GeV 2 , and the magnitude is larger; namely, it goes between the BaBar and Belle data. These features are attributed to the broader pion wave function, which enhances the small-x contribution. Compared to the form factor in the conventional approach, it shows (4 − 18)% enhancement at LO and (13 − 32)% correction at NLO under the joint resummation. With the strong suppression at the end point, the major contribution does not come from small x, but from moderate x, say, 0.1 < x < 0.2. In this range the pion wave function takes values Φ I < Φ II < Φ III as revealed in Fig. 3. It explains why the curves for Q 2 F (Q 2 ) under the joint resummation ascend fastest with Q 2 in Fig. 5(c), a result not expected from [25,26,27,28].
The above reasoning implies that a pion wave function with two humps, such as the Chernyak-Zhitnitsky model [52] or the Bakulev-Mikhailov-Stefanis model [59], which further lift the values in 0.1 < x < 0.2, will overshoot the BaBar data in our formalism. To improve the description of the data, instead, a smoother pion wave function, lying between the asymptotic and flat models in the range 0.1 < x < 0.2, serves the purpose as hinted by Figs. 5(a) and 5(b). The examples include a non-asymptotic model with a smaller Gegenbauer moment a 2 and a model in Eq. (3.17) with a fractional power α < 1. We present the LO and NLO results for Q 2 F (Q 2 ) under the joint resummation and the input of the pion wave function with a 2 = 0.05 2 in Fig. 6. Fairly speaking, this model wave function describes reasonably well the CLEO, BaBar, and Belle data in the whole range of Q 2 . In summary, the significance of the NLO contribution and the saturation behavior of the pion transition form factor are quite different under the joint resummation and the conventional resummations. Therefore, it is important to have appropriate treatment of QCD logarithmic corrections to a process, before its data are used to extract a hadron wave function. It is also crucial to clarify the high Q 2 data of the pion transition form factor on the experimental side, in order to acquire better understanding of the hadron structure and stringent scrutinization of perturbation theory.
To illustrate theoretical uncertainties, we vary the default choice of the hard scale into t = max(2 √ xQ, 1/b). As shown in Fig. 7, the scale variation increases the pion transition form factor in the large Q 2 region by approximately 8% and 1% at LO and NLO, respectively, for the asymptotic pion wave function. On the other hand, tuning the hard scale magnifies the QCD correction, as large as 7%, to the pion transition form factor in the flat model. This is certainly not unexpected, taking into account the highlighted role of the single logarithm ln(t 2 b/( √ xQ)) in the hard kernel. Similar observation also holds for the non-asymptotic model with a finite a 2 , albeit with the NLO correction being enhanced to 13% at NLO. The range of the above numerical results are basically consistent what was obtained in [57].
2 Such a value of the second Gegenbauer moment still lies in a very conservative bound 0 ≤ a2(1 GeV) ≤ 0.3 proposed in Eq. (9) of [60]. Moreover, a2 in the TMD pion wave function needs not to to be the same as in the pion LCDA.

CONCLUSION AND DISCUSSION
Applying the resummation technique with off-light-cone Wilson lines, we have constructed an evolution equation to resum the mixed logarithm ln x ln(ζ 2 P /k 2 T ) in the TMD pion wave function. The joint-resummation improved pion wave function modifies both the longitudinal and transverse momentum distributions. As a consequence, the moderate x and small b regions are more highlighted compared to the case with the conventional threshold and k T resummations. We stress that the joint resummation, organizing all the important logarithms in the pion wave function and in the hard kernel, is a treatment more appropriate and complete than the conventional resummations. In particular, Eq. (4.8) derived in this work represents the first scheme-independent k T factorization formula for the pion transition form factor in the presence of the light-cone singularity.
We have examined the significance of the NLO contribution and the saturation behavior of the pion transition form factor at high energy under the joint resummation. Differences from those under the conventional resummations were noticed, indicating that QCD logarithmic corrections to a process must be handled appropriately, before its data are used to extract a hadron wave function. Our predictions for the pion transition form factor have been confronted with the measurements from CLEO, BaBar and Belle by testing three models for the pion wave function. The comparison shows that a smooth pion wave function is favored over a pion wave function with two humps in our formalism. It turns out that a non-asymptotic pion wave function with a small second Gegenbauer moment a 2 = 0.05 describes reasonably well the CLEO, BaBar, and Belle data in the whole range of Q 2 . Resolving the discrepancy between the BaBar and Belle measurements will definitely improve our understanding towards the hadronic structure of a pion.
Our scheme-independent formalism can be extended to the k T factorization of more complicated exclusive processes. We will demonstrate this extension taking the pion electromagnetic form factor as an example. The first step is to verify that the choice of ζ 2 0 , with N and b being replaced by N 1(2) and b 1(2) for the incoming (outgoing) pion, defined in Eq. (2.37) diminishes the large logarithms in the NLO TMD pion wave function in Eq. (33) of [8]. It is indeed the case under the power counting k 2 1T ∼ k 2 2T ∼ x 1 x 2 Q 2 , confirming the universality of a TMD hadron wave function. To eliminate the large logarithms in the hard kernel given by Eq. (35) of [8], we may set , which arise from the simultaneous solution to the evolution equations for two TMD pion wave functions and one hard kernel. We will present the details of the joint-resummation improved factorization for the pion electromagnetic form factor elsewhere.
YLS acknowledges the support of National Science Foundation of China under Grant No. 11005100. YMW is supported by the DFG Sonderforschungsbereich /Transregio 9 "Computergestützte Theoretische Teilchenphysik".