Twist-3 contributions to γγ→π+π-,K+K-\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma \gamma \rightarrow \pi ^+\pi ^-,K^+K^-$$\end{document} processes in perturbative QCD approach

As one of the simplest hadronic processes, γγ→M+M-\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma \gamma \rightarrow M^{+}M^{-}$$\end{document} (M=π,K\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$M=\pi ,K$$\end{document}) could be a good testing ground for our understanding of the perturbative and nonperturbative structure of QCD, and will be studied with high precision at BELLE-II in the near future. In this paper, we revisit these processes with twist-3 corrections in the perturbative QCD approach based on the kT\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$k_{T}$$\end{document} factorization theorem, in which transverse degrees of freedom as well as resummation effects are taken into account. The influence of the distribution amplitudes on the cross sections are discussed in detail. Our work shows that not only the transverse momentum effects but also the twist-3 corrections play a significant role in the processes γγ→M+M-\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\gamma \gamma \rightarrow M^{+}M^{-}$$\end{document} in the intermediate energy region. Especially in the few GeV region, the twist-3 contributions become dominant in the cross sections. And it is noteworthy that both the twist-3 result of the π+π-\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\pi ^{+}\pi ^{-}$$\end{document} cross section and that of the K+K-\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$K^{+}K^{-}$$\end{document} cross section agree well with the BELLE and ALEPH measurements. For the pion and kaon angular distributions, there still exist discrepancies between our results and the experimental measurements. Possible reasons for these discrepancies are discussed briefly.


Introduction
A meaningful and historic subject of the perturbative QCD, which has received considerable attention in the past few decades, is the study of the exclusive processes at large momentum transfer. The pioneer works in this area are performed by Efremov and Radyushkin [1] and, independently, Brodsky and Lepage [2]. They pointed out that the exclusive processes at large momentum transfer can be factora e-mail: wangj@mails.ccnu.edu.cn b e-mail: hejk@mails.ccnu.edu.cn c e-mail: zhoumz@swu.edu.cn ized into perturbatively calculable kernels and hadronic wave functions. Although there is general agreement that perturbative QCD is able to make reliable predictions for the exclusive processes in the large energy region [3], the applicability of perturbative QCD to these processes in the intermediate energy region has been developed in controversy. In the field of exclusive processes, two-photon collisions γ γ → π + π − , K + K − [4][5][6][7][8] are the specially ones with the initial states simple and controllable and the strong interactions only in the final states. These characteristics make the processes γ γ → π + π − , K + K − be a good testing ground for our understanding of the perturbative and nonperturbative structure of QCD.
In the perturbative QCD approach based on the collinear factorization, the first investigation of the processes γ γ → π + π − , K + K − was carried out by Brodsky and Lepage [9]. At leading-twist leading-order level, they gave detailed calculations of the helicity amplitudes of these processes and presented the known relation between the differential cross section and the time-like electromagnetic form factor: Several years later, Nižić [10] performed the leadingtwist next-to-leading-order perturbative QCD calculations for these processes. In Refs. [10,11], it was found that both the leading-order and next-to-leading-order perturbative QCD predictions are almost an order of magnitude smaller than the experiment data [12][13][14][15][16] in the intermediate energy region. Motivated by the calculation of the pion electromagnetic form factor [17] in which the nonleading twist contribution becomes the dominant one in a few GeV region, Gorsky [18] investigated the two-photon process γ γ → π + π − with the higher-twist corrections. It was pointed out that the twoparton twist-3 contributions are logarithmically larger than the leading-twist contributions in the intermediate energy region due to the chiral enhancement effects. And this conclusion had also been confirmed in Ref. [19] with the Brodsky-Huang-Lepage (BHL) prescription [20]. But in a few GeV region, the predicted cross sections with higher-twist corrections [18,19] are several times larger than their experimental measurements [12][13][14][15][16].
As it is known that the higher-power effects from the intrinsic transverse momentum play a crucial role in the pion electromagnetic form factor at the scale of few GeV [21,22], one may expect that the same situation can also be found in the two-photon processes γ γ → π + π − , K + K − . Based on the k T factorization theorem [21,23] in which the transverse momentum dependence was retained, the leadingtwist perturbative QCD predictions for the cross sections of these two-photon processes were obtained in Refs. [24][25][26][27][28][29]. From the analysis of the differential cross sections obtained in both perturbative QCD and QCD sum rule, the conclusion that the transition from nonperturbative to perturbative QCD in γ γ → π + π − , K + K − occurs at Q ≈ 2 GeV was drawn in Refs. [25,26,29]. However, the numerical results in Refs. [27,28] show that the twist-2 cross sections of the γ γ → π + π − , K + K − in the k T factorization are still much smaller than the experiment data in the intermediate energy region.
To consider the higher-twist effects for the processes γ γ → π + π − , K + K − in k T factorization would be a natural next step, and seems necessary because the highertwist effects and the higher-power effects from the transverse momentum become important simultaneously in intermediate energy region. Compared with the previous calculations in collinear factorization [18,19], the higher-twist calculations in k T factorization have several major improvements: Firstly, both the transverse momentum dependence of hadronic wave functions and that of hard kernels are included. The transverse momentum dependence of the hard kernel regularizes the end-point singularities of the internal propagators. While the transverse momentum dependence of the wave function makes a sizable suppression to the perturbative QCD contributions, and cannot be ignored as pointed out in Ref. [22], especially in the few GeV region. Secondly, the corrections from the Sudakov resummation [21] substantially suppress the nonperturbative contributions from both the soft end-point regions and the large-b regions. Moreover, the corrections from the threshold resummation [30] suppress the nonperturbative contributions from the end-point regions further and eliminate the end-point singularity without a artificial cut-off [18] or the BHL prescription [20] in the twist-3 calculations. Thirdly, the scales of the coupling constant α s are chosen to be momentum-fraction dependent in order to avoid large logarithms from higher-order perturbative QCD corrections. All of the aforementioned improvements make the perturbative calculation become more selfconsistent, even for momentum transfer as low as a few GeV.
In this paper, we present a detailed twist-3 calculations for the two-photon processes γ γ → π + π − , K + K − in the perturbative QCD approach based on the k T factorization theorem. At the twist-2 level, our results are consistent with the predictions given in Refs. [27,28]. While, with the twist-3 corrections, it is found that both the π + π − and the K + K − cross sections agree well with the experiment data [13][14][15]. This paper is organized as follows: In Sect. 2, the related conventions to these processes are introduced at first, next the twist-2 and twist-3 light-cone wave functions are described briefly, then we perform the calculation of the hard kernels with twist-3 corrections. The numerical analysis and discussion are given in Sect. 3. The last section is our summary and conclusion. The expression of the Sudakov function can be found in "Appendix".

Kinematics and conventions
We consider the simplest hadronic processes, where the incoming photons are real with momenta p 1 , p 2 and polarization vectors ε λ 1 1 , ε λ 2 2 with the photon helicities λ 1 , λ 2 = ±1, M denotes the out-going pseudoscalar mesons with momenta p 3 , p 4 . The contributions to the γ γ → (qq) + (qq) amplitude arise from 20 tree-level Feynman diagrams, which can be grouped together into four basic types (a, b, c and d) based on the virtualities of the internal gluon. We depict the four basic Feynman diagrams in Fig. 1. The remaining diagrams arise from permutations of the gluon and quark lines.
For convenience, the amplitude is calculated in the centerof-momentum frame with the momenta of the final mesons along z-axis. θ is the scattering angle and Q denotes the twophoton center-of-mass energy. In the light-cone coordinate, a momentum A μ has the form Under these conventions, the momenta of photons and mesons read and the corresponding polarization vectors of the initial photons can be written as The total number of tree-level diagrams in groups a, b, c and d are 6, 6, 4, and 4, respectively. The grouping is based on the virtualities of the internal gluon with ω = Q 2 , s = sin θ 2 , c = cos θ 2 and p = (2ωsc, 0) for abbreviation. The momenta of the quark and antiquark k 1 , k 2 in the M + meson and k 3 , k 4 in the M − meson as labeled in Fig. 1 can be written as . Here x (y) andx (ȳ) are the longitudinal momentum fractions of quark and antiquark, respectively, and the transverse momenta k ⊥1 , k ⊥2 are relate to the momenta of the two out-going mesons.

Light-cone wave functions
The Fierz identity is generally employed to factorize the fermion flow: Here q andq are the quark and antiquark field. The Dirac structure γ μ γ 5 in Eq. (5) contributes at the twist-2 level, while iγ 5 and σ μν γ 5 contribute at the twist-3 level.
With the light-cone expansion to twist-3 accuracy, the twoparton light-cone wave functions of π − can be expressed via the matrix element of the bilocal operator [31][32][33][34] where f π refers to the decay constant of pion and the mass parameter μ π is defined as where m u , m d are the current quark masses and m π denotes the pion meson mass. π π is the twist-2 wave function and p π , σ π are two-parton twist-3 wave functions. Due to the SU (2) isotopic symmetry, the wave functions of π + are identical with the wave functions of π − . Generally, the twist-3 contributions are power-suppressed. However, compared with the twist-2 contributions, the suppression factor μ π /Q in the two-parton twist-3 contributions is usually not small enough because of the chiral enhancement effects, i.e., through the relations between light quark masses in the chiral perturbation theory, μ π (2 GeV) = 2.43 GeV [35]. Especially, in the intermediate energy region, such as the BELLE experimental range 2.4 GeV < Q < 4.1 GeV [15], one can obtain the suppression factor μ π /Q ∼ O (1). That is the reason why, in some processes [18,19,[36][37][38][39][40], the two-parton twist-3 contributions could become even more important than the twist-2 contributions in a few GeV region.
Transforming the matrix element of the bilocal operator in Eq. (6) from the coordinate space into the momentum space, one can obtain the light-cone projection operator [34,41,42] . From the first principles in QCD, the exact transverse momentum dependence of the above wave functions is still unknown. For the convenience of calculations, we assume that π π , p π and σ π have the same transverse momentum dependence. Following the ansatz made in Refs. [22,43], the wave functions of the pion can be separated into two parts One part is the distribution amplitudes φ π π , φ p π and φ σ π . Using the conformal symmetry, the distribution amplitudes can be expanded in the Gegenbauer polynomials [33] with the normalization conditions 1 0 dxφ π π (x, μ) = 1, Here a π 2,4 are the Gegenbauer moments, the twist-3 parameters η 3π , ρ 2 π are defined as (12) and the parameters f 3π , ω 3π are defined in terms of matrix elements of the local three-parton twist-3 operators [33] 0|ūσ μα γ 5 The other part is the k ⊥ -dependence function π (x, k ⊥ ).
With the constraint from the π − → μ −ν μ decay [20], one can obtain i.e., the function π (x, k ⊥ ) satisfy the normalization condition An appropriate candidate for the function π (x, k ⊥ ) is usually assumed to be a simple Gaussian form [22,44] π (x, k ⊥ ) = where the oscillator parameter β π is determined by the mean square transverse momentum The detailed discussions about the choice of the oscillator parameter β π are shown in Refs. [22,44]. In this work, we will take the value k 2 ⊥ 1/2 π = 0.35 GeV, which is compatible with the π 0 → γ γ constraint [20].
For the K − case, the matrix element of the bilocal operator can be expressed as [45] where f K refers to the decay constant of kaon, and the chiral enhancement parameter is defined as μ K = m 2 K /(m u + m s ), m u and m s are the current quark masses and m K denotes the kaon meson mass. We should note that the wave functions of K + are given by i [45]. Similar to the pion case, the wave functions of the kaon meson can also be written as with Here we have adopted the approximation with the normalization conditions Here the SU (3) symmetry breaking effects have already been taken into account. The parameters a K 1,2 are the first two nonvanishing Gegenbauer moments, and the twist-3 parameters where the parameters f 3K , ω 3K , λ 3K are defined in terms of matrix elements of the local three-parton twist-3 operators [33] 0|ūσ Using the leading-order renormalization group equations, the scale dependence of the various parameters is given by for the pion and for the kaon with the evolution factor L = α s (μ 2 )/α s (μ 2 0 ) and β 0 = (33 − 2n f )/12. To the leading logarithmic accuracy, the QCD running coupling constant can be expressed as and the anomalous dimensions are given by , with C F = 4/3 and C A = 3.
In perturbative QCD approach, the convolution of wave functions and hard kernels needs to be performed in the transverse configuration b-space. With the definition one can obtain that the full soft wave functions of the pion and kaon in the transverse configuration b-space have the formŝ where the Fourier transformed functionsˆ M (x, b) read

The helicity amplitudes with twist-3 contributions
In the perturbative QCD approach, the helicity amplitudes can be expressed as the convolution between the soft hadron wave functions and the hard kernels with respect to both the longitudinal momentum fractions and the transverse separation of quark and antiquark. To the π + π − production in two-photon collision, the helicity amplitudes have the form [26,29] where T λ 1 λ 2 i j with the transverse momentum k ⊥1 , k ⊥2 kept can be calculated perturbatively by Feynman rules. In the helicity amplitudes, the scripts i, j = π correspond to the twist-2×twist-2 term, and the scripts i, j = p, σ, σ correspond to the two-parton (qq) twist-3×twist-3 terms. It is noteworthy that the two-parton twist-2×twist-3 terms vanish because of the spin structure of the corresponding hard kernels. Compared with the twist-2×twist-2 contributions, the contributions from the two-parton twist-3×twist-3 ones are generally regarded as power-suppressed. However in a few GeV region, the suppression factor is μ 2 π /Q 2 ∼ O(1) as a consequence of the chiral enhancement effects. So the two-parton twist-3×twist-3 contributions become important in this energy region. While, compared with the twist-2×twist-2 contributions, the contributions from the threeparton (qqg) twist-3×twist-3 terms and that from the twoparton twist-2×twist-4 terms are strongly suppressed by the factors f 3π /( f π Q) and m 2 π /Q 2 respectively. According to the factors f 3π /( f π Q) and m 2 π /Q 2 being at the order of 10 −2 in the intermediate energy region, these two contributions would be much smaller than the twist-2×twist-2 and twoparton twist-3×twist-3 contributions in this energy region. A similar analysis have also been given in the studies of the time-like electromagnetic form factors [38,46]. Recently, the explicit calculations of the three-parton twist-3×twist-3 [40,47] and two-parton twist-2×twist-4 contributions [40] for the pion electromagnetic form factors have been performed, and the numerical results have shown that these two contributions are at least an order of magnitude smaller than the twist-2×twist-2 and the two-parton twist-3×twist-3 contributions, which is consistent with the previous analysis [46]. So in this work, we will neglect the contributions from the three-parton twist-3×twist-3 terms and the two-parton twist-2×twist-4 terms.
With the transverse momentum k ⊥ retained, an important property of the hard kernels needed to make clear is the gauge invariance. Here we present a brief analysis. As pointed out in the Refs. [48,49], one may encounter the problem of gauge dependence caused by the off-shellness of partons in a k ⊥ -dependent hard kernel. For exclusive processes, it have been proven that the hard kernels derived in the k T factorization are gauge invariance to all orders in α s at twist-2 level [50,51]. While at high twist level, the situation becomes more complex. As the results shown in the π 0 [42] and ρ 0 [52] photoproduction processes, the two-parton twist-3 contributions must be combined with the three-parton twist-3 contributions in order to guarantee gauge invariance. In fact, the gauge dependence of the two-parton twist-3 hard kernel is due to the terms proportional to the transverse momentum in the hard kernel, and it have been shown in Ref. [47] that the gauge dependence proportional to the transverse momentum in the numerator of the two-parton hard kernel and the gauge dependence associated with the three-parton Fock state exactly cancel each other. In other words, if the terms proportional to the transverse momentum in the numerators of the two-parton hard kernels (power-suppressed by 1/Q) and the three-parton Fock state (strongly suppressed by the nonperturbative parameter f 3π ) are precluded simultaneously, the two sources of gauge dependence would disappear, and the k ⊥ -dependent two-parton hard kernels would become gauge invariance in the same way as the collinear case [53]. In this work, we employ the k T factorization theorem proposed in Refs. [30,54], in which only the transverse momentum dependence in the denominators of the two-parton hard kernels is retained. As analyzed above, the k ⊥ -dependent two-parton hard kernels in our calculations are gauge invariance. In the following calculations of the hard kernels, the covariant gauge is adopted.
Next we will show our procedure to the Fourier transformation of the hard kernels by means of a specific Feynman diagram in the Group a in Fig. 1 as an example. The hard kernels T λ 1 λ 2 a 1 i j with the subscript "a 1 " specifically refer to the results of one of the six diagrams in Group a. Its leading-twist contributions are illustrated as with the factor κ 1 = αα s (μ R )π 2 f 2 π . Here the terms contained the transverse momentum in the numerators, which are power-suppressed, have been dropped. And the denominators of the corresponding quark and gluon propagators arẽ with The convolution formula in Eq. (32) involves two independent scales of this process, the factorization scale μ F and the renormalization scale μ R . In analogy to the cases of the pion electromagnetic form factor [21,37,55] and transition form factor [56][57][58], the transverse separation b 1 and b 2 of the pion wave functions provide the factorization scales μ Fi = 1/b i (i = 1, 2), below which the QCD dynamics is regarded as being nonperturbative and can be absorbed into the the soft wave function. While the renormalization scale μ R = t is set to the largest mass scale associated with the momentum of the internal hard gluon propagators in order to minimize the higher-order QCD corrections.
To simplify the calculation, the hierarchy Q 2 x Q 2 ∼ y Q 2 x y Q 2 , k 2 1⊥ , k 2 2⊥ is postulated in the small-x, y region, as elaborated in Refs. [55,59]. Under this hierarchy, one can ignore the transverse momentum dependence in the quark propagators and just keep it in the gluon propagator to regulate the endpoint singularity, then the hard kernels in Eq. (35) are simplified as By using the Fourier transformation of the propagators regularized with the i prescription one can obtain the Fourier transformed hard kernels (b = |b 1 |) Here H (1) 0 and K 0 denote Hankel and modified Bessel function. As a consequence of the above simplification, the hard kernels in Eq. (40) only depend on a single b parameter. The underlying physics picture is that, the virtual quark lines involved in the hard kernels are thought of as being far from mass shell, and shrunk to a point [21,26].
There are two types of resummation of the higher-order effects: the region with large transverse separation b for Sudakov resummation and the small longitudinal momentum fraction x region for threshold resummation. The Sudakov resummation of the double logarithms α s ln 2 [x/(Q 2 b 2 )] produced by overlapping collinear and soft divergencies for massless quarks as well as the renormalization group equation transformation from the factorization scale 1/b to the renormalization scale t, are incorporated into the Sudakov factor exp[−S]. In next-to-leading logarithm approximation, the Sudakov exponent S reads [21,26,60] S(x, y, b, where the function s(x, b, Q) is originally derived by Botts and Sterman [23] and later on slightly improved, for instance, by Li [21] and Kroll [60] et al. The explicit expression of s(x, b, Q) is given in "Appendix". One can find that the Sudakov factor exhibits a strong fall off at large b. And this property makes the nonperturbative contributions from large b, no matter what x is, less important.
Since there exists the end-point enhancement in the twoparton twist-3 contributions (logarithmical enhancement) for the scattering process γ γ → π + π − [18,19], the Sudakov factor exp[−S] is still not effective enough to suppress the nonperturbative contributions from small x region, which would spoil the perturbative calculation. It had been argued in Ref. [30] that as the end-point region is important, the corresponding large double logarithms α s ln 2 x that arise from the higher-order corrections also need to be organized to all orders, and into a jet function S t (x) as a consequence of threshold resummation. In next-to-leading logarithm accuracy, the jet function S t (x) can be parameterized into a universal form [30,61] In Ref. [57], Li and Mishima have proposed a parabolic parametrization for the parameter c, and frozen c at c = 1, when it exceeds unity: Since the factor S t (x) drops rapidly as x → 0, 1, the endpoint singularities are eliminated and the nonperturbative contributions from the dangerous end-point regions are suppressed effectively. After including the above two types of resummation, both the twist-2 and twist-3 contributions are well-behaved, and that results in the perturbative calculation more selfconsistent. Performing the integration over b 2 , the helicity amplitudes in Eq. (32) can be then rewritten as a much clearer and simpler form in a single-b form Making the substitution d 2 b → bdbdϕ and performing the integral over the polar angle ϕ, we can obtain the helicity amplitudes for γ γ → π + π − process in the conjugate b space by summing up all of the 20 Feynman diagrams with the twist-2 hard kernelŝ F(s a , b), and the nonzero twist-3 hard kernelŝ (s c , b), (s d , b), (s c , b), (s d , b), (s a , b), (s d , b), (s a , b), (s c , b), (s c , b), (s d , b), (s c , b), (s d , b), (s c , b), (s d , b), y  F(s b , b), (s c , b), F(s a , b), Here the factor κ 2 is defined as Compared with the twist-2 hard kernels, the twist-3 ones are easily found to be suppressed by the factor κ 2 /κ 1 = μ 2 π /ω 2 . In the above expressions for theT λ 1 λ 2 ni j , J 0 represents the Bessel function with the variable pb = 2ωscb and the func- The renormalization scales t n (n = a, b, c, d), enteringT H as the arguments of strong coupling constant α s and depending on the kinematics of the specific group, are chosen as The factors s n (n = a, b, c, d) come from the gluon propagators in each group In addition, it is worth noticing the hard kernels with different helicity have the relationŝ For the γ γ → K + K − process, the hard kernels can be obtained directly with the replacements of e d → e s , f π → f K , μ π → μ K . And the corresponding helicity amplitudes can be expressed as

Numerical analysis
The differential cross sections of the processes γ γ → M + M − (M = π, K ) can be expressed as In order to make comparison with the experimental measurements, we integrate over the scattering angle in the region | cos θ | < 0.6 to get the cross sections σ 0 for both π + π − and K + K − processes. In this work, we focus on the intermediate energy region 1.0 GeV < Q < 7.0 GeV to check the applicability of the perturbative QCD in γ γ → π + π − , K + K − processes.
In one-loop accuracy, the QCD running coupling is given by α s (μ 2 ) = π/[β 0 ln(μ 2 / 2 QC D )] with the QCD scale QC D = 0.2 GeV and the interaction scale μ = t n  Table 2 Gegenbauer moments of five sample models of the pion and kaon distribution amplitudes at the scale μ 0 = 1 GeV  a, b, c, d) corresponding to each group in Fig. 1. The most important parameters used in our numerical analysis are listed in Table 1. In addition, the decay constants are taken as f π = 130.4 MeV and f K = 159.8 MeV. For the meson masses and the quark mass entering the pion and kaon distribution amplitudes, we take m π = 139.6 MeV, m K = 493.7 MeV and m s (2 GeV) = 95 MeV in MS scheme. All of the values are quoted from PDG [64]. The Gegenbauer moments a π 2,4 and a K 1,2 in the pion and kaon distribution amplitudes come from QCD sum rule [35,62]. For the chiral enhancement parameters μ π and μ K , we employ the well known chiral perturbative theory relations [65] R = 2m s m u + m d = 24.4 ± 1.5, and obtain The remaining parameters of the twist-3 distribution amplitudes in our calculations are taken from Ref. [45] and shown in Table 1 for brevity. Fig. 2 Twist-2 results of the cross sections σ 0 (π + π − ) and σ 0 (K + K − ) with the five sample models of the distribution amplitudes listed in Table. 2.
The points with errors are the experiment data from Refs. [12,13,15] Fig. 3 Twist-3 results of the cross sections σ 0 (π + π − ) and σ 0 (K + K − ) with the five sample models of the distribution amplitudes listed in Table. 2.
The points with errors are the experiment data from Refs. [12,13,15] The theoretical uncertainty in this work mainly comes from the nonperturbative inputs of meson distribution amplitudes: one is the Gegenbauer moments a M n and the other is the chiral enhancement parameters μ M (n = 2, 4 for pion and n = 1, 2 for kaon, M = π, K ). In Table 2, we present five different models of the distribution amplitudes. And the corresponding Gegenbauer moments are in the region which can cover the ones given in the recent lattice determinations [66].
With the five models of the pion and kaon distribution amplitudes, the cross sections σ 0 (M + M − ) at twist-2 and twist-3 levels are illustrated in Figs. 2 and 3 respectively. For the π + π − cross sections shown in Fig. 2, one can find that the lines of model I, IV and V with a π 2 + a π 4 close to the same value (∼ 0.23) almost coincide with each other. And the line of model II (III) with the smallest (largest) a π 2 and a π 4 have the largest (smallest) slope. † It indicates that the twist-2 result of σ 0 (π + π − ) depends on the Gegenbauer momentums almost in the combination of a π 2 + a π 4 . Moreover, by fixing a π 2 (a π 4 ), i.e., comparing the result of model II with that of model IV (V), the slope of the pion lines is found to become large as a π 4 (a π 2 ) decrease. The same conclusion can also be obtained by comparing the result of model III with that of model V (IV). While, for the K + K − cross sections shown in Fig. 2, one can find that the lines of model II and V (III and IV) almost coincide with each other, which implies the twist-2 result of σ 0 (K + K − ) is insensitive to the variations of a K 1 . Furthermore, the lines of model II and V (III and IV) with the smallest (largest) a K 2 have the largest (smallest) slope, i.e., the slope of the kaon lines become large as a K 2 decrease. The points with errors are the experiment data from TPC [12], BELLE [15] and ALEPH [13] However, with all the five distribution amplitudes, we find that the twist-2 result of σ 0 (K + K − ) is several times smaller than its experiment data [12,13,15], and the twist-2 result of σ 0 (π + π − ) is almost an order of magnitude smaller than the data [12,13,15]. Similar conclusion was already drawn in Refs. [27,28]. At twist-3 level, to assess the uncertainties related to the Gegenbauer moments, we draw the figures at specific chiral enhancement parameters. From Figs. 2 and 3, one can see that the relative uncertainties come from the Gegenbauer moments in the twist-3 results are much smaller than that in the twist-2 results. With the Gegenbauer moments chosen at their central values and the chiral enhancement parameters being in their reasonable ranges, the twist-3 cross sections are reshown in Fig. 4. Comparing Fig. 3 with Fig. 4, it is found that the main uncertainties in the twist-3 contributions come from the chiral enhancement parameters. For both the π + π − and K + K − processes, the twist-3 cross sections are found to be in good agreement with the BELLE measurements [15] in the energy region 2.4 GeV < Q < 4.1 GeV and the ALEPH measurements [13] in the energy region 2.0 GeV < Q < 6.0 GeV. Especially, when the chiral enhancement parameters are taken the upper value μ π (2 GeV) = 2.80 GeV and μ K (2 GeV) = 2.75 GeV, the agreements are remarkable. But, in the relatively low energy region 1.0 GeV < Q < 2.0 GeV, there still exists considerable discrepancies between our results and the experiment data [12]. For the π + π − process, the discrepancy in 1.0 GeV < Q < 2.0 GeV may be ascribed to the interference of the continuum with the resonances, such as the f 2 (1270). For the K + K − process, the neglection of the twist-4 contributions, which are suppressed by the factor m 2 K /Q 2 , may lead to the discrepancy in 1.0 GeV < Q < 2.0 GeV.
In Fig. 5, we make a comparison between the twist-2 and the twist-3 results of the cross sections σ 0 (π + π − ) and σ 0 (K + K − ) with all the nonperturbative input parameters chosen at the central values. Here we should note that the twist-2 results in the collinear factorization could be enhanced by about 2-3 times if modern distribution amplitudes are used, which would be similar with the situation in the k T factorization (shown in Fig. 2). Compared with the previous twist-2 results [19] obtained in the collinear factorization, the twist-2 results in the k T factorization are enhanced several times in the energy region 2.0 GeV < Q < 7.0 GeV, but suppressed significantly in the relatively low energy region 1.0 GeV < Q < 2.0 GeV. The enhancement in 2.0 GeV < Q < 7.0 GeV may be due to the higher-power corrections from the transverse momentum k ⊥ . While the suppression in 1.0 GeV < Q < 2.0 GeV may be ascribed to the corrections from the sudakov resummation, which have been organized into the suppression factor exp[−S]. As pointed out in Ref. [24], the asymptotic scattering amplitudes for the processes γ γ → π + π − , K + K − are insensitive to the Sudakov corrections. It indicates that the Sudakov factor exp[−S] affects the cross sections slightly in the relatively high energy region, but may lead to a strong suppression at the relatively low energy scale. Comparing the twist-2 results with the twist-3 ones in the k T factorization, one can conclude that the twist-3 corrections and their interference with the k ⊥ effects are very significant for both the π + π − and K + K − processes in the intermediate energy region 1.0 GeV < Q < 7.0 GeV. Especially in the relatively low energy region 1.0 GeV < Q < 2.0 GeV, the twist-3 results of σ 0 (π + π − ) and σ 0 (K + K − ) are enhanced by almost an order of magnitude as compared with the twist-2 results. Similar situations can also be found in the calculations of the pion-photon transition form factor [59] and the pion electromagnetic form factor [37,39,55]. A more detailed investigation has shown that the twist-3 corrections fall off rapidly with increasing Q and, beyond 10 GeV region for both pion and kaon, fall below the twist-2 contributions. As is known that, the twist-2 contributions are expected to dominate in the cross section at the asymptotically large momentum transfers. However, as evident from Fig. 5, the twist-3 contributions should be taken into account in the intermediate energy region 1.0 GeV < Q < 7.0 GeV. At this moment, one may worry about possible contributions from the three-parton Fock state, which can also contribute to the cross sections at twist-3 level. As pointed out in Refs. [38,40,46,47], the three-parton twist-3 contributions are strongly suppressed by the factor f 3M /( f M μ M ) (M = π, K ) and at least an order of magnitude smaller than the two-parton twist-3 contributions for the time-like pion and kaon electromagnetic form factors in the intermediate energy region. Due to the similarity of the calculations † between the helicity amplitudes of the two-photon process and the time-like electromagnetic form factors, one can neglect the three-parton contributions safely in γ γ → π + π − , K + , K − processes. † At twist-2 level, the similarity of the calculations between the helicity amplitudes of the two-photon process (γ γ → π + π − ) and the time-like pion electromagnetic form factors (γ * → π + π − ) leads to the relation in Eq. (1), i.e., the famous BL approximation [9].  [12], ALEPH [13] and BELLE [15] The cross section ratio σ 0 (K + K − )/σ 0 (π + π − ) is presented in Fig. 6. We find that the twist-2 results both in the collinear factorization and in the k T factorization are almost equal to ( f K / f π ) 4 ≈ 2.26. Similar theoretical results have also been given in Refs. [9,11]. While the fitted result from BELLE measurements [15] is 0.89 ± 0.04 ± 0.15 for 3.0 GeV < Q < 4.1 GeV energy region, which is clearly smaller than the twist-2 results obtained in our work and others [9,11]. Including the twist-3 corrections, our result of the ratio σ 0 (K + K − )/σ 0 (π + π − ) is found to be in line with the ALEPH measurements [13] in 2.5 GeV < Q < 4.0 GeV and the BELLE measurements [15] in 3.5 GeV < Q < 4.1 GeV. In the energy region 1.0 GeV < Q < 2.5 GeV, there still exists discrepancy between our prediction and the experimental measurements. The reason may be due to the contributions from the resonances and the m 2 M /Q 2 terms, which have been neglected in our calculations, may play a significant role in this energy region. Another interesting issue we should note here, is related to the cross section ratio between the neutral and charged channels, such as the ratios σ 0 (π 0 π 0 )/σ 0 (π + π − ) and σ 0 (K 0 , which exhibit the difference of these two types of cross sections on the energy dependencies. In literature, these ratios have been predicted by different methods, such as the perturbative QCD [9] and the handbag model [67,68]. However, the predictions from both two methods are inconsistent with the experimental measurements [69,70]. More discussions about this issue can also be found in Refs. [6,7,71]. In this work, one can see that the twist-3 hard kernels shown in Eqs. (48)-(57) bring new energy dependence into the cross sections. Moreover, the Feynman diagrams in groups a (b) make contributions to the hard kernels with e 2 u (e 2 s ) for the K + K − channel, and with e 2 d (e 2 s ) for the K 0 S K 0 S channel. While the Feynman dia-  [15] grams in groups c and d make contributions to the hard kernels with e u e s = −2/9 for the K + K − channel, and with e d e s = 1/9 for the K 0 S K 0 S channel. It means that c and d diagrams enter the charged channel with a additional minus sign as compared to the neutral channel. And these two points may provide possible explanation about the different energy dependencies between these two types of channels observed by the experimental measurements [69]. Similar analysis can also apply to the ππ channels [70]. In addition, as pointed out by Brodsky [4,9], the ratio σ 0 (π 0 π 0 )/σ 0 (π + π − ) is sensitive to the shape of pion distribution amplitude, which can be treated as a possible constraint about the pion distribution amplitude through comparing the theoretical prediction of the ratio with its experimental measurement.
At last but not least, we exhibit the pion and kaon angular distributions, i.e., the differential cross sections σ −1 0 dσ (γ γ → M + M − )/d| cos θ |, which have been normalized to the cross section σ 0 (M + M − ). Owing to the normalization, the theoretical uncertainties of the pion and kaon angular distributions from the Gegenbauer moments and the chiral enhancement parameters cancel to a large extent. Figure 7 shows our results for Q = 3.2−3.3 GeV which is in the middle range of the BELLE measurements. One can find the twist-2 results in the collinear factorization and those in the k T factorization give nearly the same angular distributions. While the angular distributions with the twist-3 corrections in the k T factorization become more steady over the region | cos θ | < 0.6. However, the agreement between the angular distributions obtained in this work and those of the BELLE measurements [15] is not perfect. The reason may be related to the simplification made in Eq. (38), which is illegal in some angular coverage. In Fig. 7, we also show the known BL result [9], corresponding to the blue dot-dashed line 1.227/ sin 4 θ , which is obtained by using the approximate relation in Eq. (1) with F M (Q 2 ) ≈ 0.4 GeV 2 /Q 2 . However, the CLEO measurements [72] gave a larger form factor F M (Q 2 ) ≈ 1.01 GeV 2 /Q 2 , which would put the known BL result [9] in question. As a matter of fact, the approximate relation in Eq. (1) was obtained at the twist-2 level, while the twist-3 corrections play a key role in both the two-photon processes and the time-like electromagnetic form factor in the intermediate energy region. It implies that the BL relation shown in Eq. (1) may be not a good approximation in the intermediate energy region.

Summary and conclusion
Our work represents the first investigation of the twist-3 corrections to the two-photon processes γ γ → π + π − , K + K − in the perturbative QCD approach based on the k T factorization theorem. The transverse momentum dependence and the resummation effects are included in the perturbative QCD approach. We use the twist-2 and twist-3 light-cone wave functions incorporating transverse degrees of freedom as the nonperturbative dynamical inputs. The nonperturbative contributions from the end-point regions can be effectively suppressed by the Sudakov factor exp[−S] and the threshold resummation factor S t (x). Consequently, the perturbative QCD calculation becomes more self-consistent and applicable, especially in the few GeV region.
Within the uncertainties from the distribution amplitudes, the twist-2 results of the cross sections σ 0 (π + π − ) and σ 0 (K + K − ) are much smaller than the experiment data [12,13,15]. As shown in our numerical analysis, it is found that both the transverse momentum effects and the twist-3 corrections play a significant role in the two-photon processes γ γ → π + π − , K + K − in the intermediate energy region. In the relatively high energy regions, the twist-3 results of the cross sections σ 0 (π + π − ), σ 0 (K + K − ) and their ratio σ 0 (K + K − )/σ 0 (π + π − ) are in good agreement with their corresponding experiment data [13,15]. But in the relatively low energy region, the predicted cross sections and their ratio are still in disagreement with the experimental measurement [12]. From the above analysis, we can conclude that the cross sections for γ γ → π + π − , K + K − with angular region | cos θ | < 0.6 may be dominated by the perturbative QCD contributions when the collision energy Q 2 GeV. And it is noteworthy that, by analyzing the differential cross sections, the authors in Refs. [25,26] have concluded that the transition from nonperturbative to perturbative QCD in the processes γ γ → π + π − , K + K − is at Q ≈ 2 GeV and θ ≈ 40 • . Obviously, these two conclusions are consistent with each other.
However, our results of the pion and kaon angular distributions are different from the BELLE measurements [15]. The reason may be partly due to the simplification made in Eq. (38), which is illegal in some angular coverage. In addition, in the relatively low energy region 1.0 GeV < Q < 2.0 GeV, two-photon processes γ γ → π + π − , K + K − have been investigated by other methods, such as the QCD sum rule [25,73] and the chiral perturbation theory [74,75], which can be seen as the complementation of the perturbative QCD approach for a full understanding of these processes.