Global Analysis of hadronic two-body B decays in the perturbative QCD approach

Based on the flavor structure of four-quark effective operators, we develop an automatic computation program to calculate hadronic two-body $B$ meson decay amplitudes, and apply it to their global analysis in the perturbative QCD (PQCD) approach. Fitting the PQCD factorization formulas for $B\to PP,VP$ decays at leading order in the strong coupling $\alpha_s$ to measured branching ratios and direct CP asymmetries, we determine the Gegenbauer moments in light meson light-cone distribution amplitudes (LCDAs). It is found that most of the fitted Gegenbauer moments of the twist-2 and twist-3 LCDAs for the pseudoscalar meson $P$ ($P=\pi$, $K$) and vector meson $V$ ($V=\rho$, $K^*$) agree with those derived in QCD sum rules. The shape parameter for the $B_s$ meson distribution amplitude and the weak phase $\phi_3(\gamma)=(75.1\pm2.9)^\circ$ consistent with the value in Particle Data Group are also obtained. It is straightforward to extend our analysis to higher orders and higher powers in the PQCD approach, and to the global determination of LCDAs for other hadrons.


I. INTRODUCTION
The study of heavy quark physics is firmly in the precision era nowadays. On the experimental side, the B factories, i.e., the BaBar, Belle, and LHCb have collected abundant data of exclusive B meson decays, which can be employed not only to explore involved rich QCD dynamics, but also to probe the origin of CP violation and potential new physics signals [1]. Vastly more data will be still accumulated by the upgraded LHCb and Belle-II Collbarotions [2][3][4]. On the theoretical side, tremendous progress on the development of QCD treatments of exclusive B meson decays with controllable uncertainties has been achieved. Strict confrontation between data and theoretical expectations has led to some mild tensions between experimental observations and the Standard Model [1], which may be vaguely attributed to new physics beyond the Standard Model. This undoubtedly motivates the attempt to gain deeper understanding of QCD dynamics in exclusive B meson decays and better control of hadronic uncertainties.
The b quark mass m b is much larger than the QCD hadronic scale Λ QCD , which renders QCD analyses of exclusive B meson decays possible. Nonperturbative dynamics in heavy meson decays is reflected by infrared divergences in radiative corrections. When a factorization theorem holds, infrared divergences are absorbed into hadron light-cone distribution amplitudes (LCDAs), so that the remnant, being infrared finite, is calculable at the parton level in perturbation theory. A physical quantity, such as a heavy-to-light transition form factor, is then factorized into a convolution of a b quark decay hard kernel with hadron LCDAs in parton momentum fractions. The corresponding factorization theorem should be proved to all orders in the strong coupling α s and to certain power in Λ QCD /m b . LCDAs, despite of being nonperturbative, are universal, i.e., process-independent. With this universality, LCDAs, determined by nonperturbative methods like QCD sum rules [5][6][7] and lattice QCD [8][9][10], or extracted from experimental data, can be employed to make predictions for other modes involving the same hadrons.
The theoretical approaches based on factorization theorems in the heavy quark limit include light-cone QCD sum rules (LCSR) [11][12][13], the QCD-improved factorization (QCDF) [14], the perturbative QCD (PQCD) factorization [15][16][17][18][19][20], and the soft-collinear effective theory (SCET) [21,22]. The collinear factorization applies to relevant correlators in LCSR, where some hadronic states are expanded into parton Fock states characterized by different twists. The QCDF approach is an extension of the naive factorization assumption for hadronic two-body B meson decays in the collinear factorization theorem. The SCET for kinematic regions with energetic final state hadrons is equivalent to the collinear factorization theorem, but formulated in terms of effective operators. The k T factorization theorem is the basis of the PQCD approach, which is more appropriate in the endpoint region of parton momentum fractions. Many efforts have been devoted to systematic investigation of hadronic two-body B meson decays at various orders in α s and powers in Λ QCD /m b [23][24][25]. In all the above formalisms nonperturbative hadron LCDAs provide a major source of theoretical uncertainties.
A hadron LCDA can be expanded into a series of Gegenbauer polynomials with the coefficients, namely, the Gegenbauer moments being determined by other methods as aforementioned. Though some attempts have been made to calculate Gegenbauer moments using lattice QCD [8][9][10], not all LCDAs are constrained in this way so far. Here we will perform a global fit of the Gegenbauer moments in light meson LCDAs to measured branching ratios and direct CP asymmetries in hadronic two-body B meson decays in the PQCD approach. The decay amplitudes at leading order (LO) of the strong coupling α s will be constructed automatically with a computation program by making use of flavor SU(3) properties. We establish a Gegenbauer-moment-independent database, by means of which each decay amplitude is expressed as a combination of the relevant Gegenbauer moments and Cabibbo-Kobayashi-Maskawa (CKM) matrix elements. The Gegenbauer moments in the leading-twist (twist-2) and next-to-leading-twist (twist-3) LCDAs for the pseudoscalar meson P (P = π, K) and vector meson V (V = ρ, K * ) are then fixed in the global fit, most of which are found to agree with those from QCD sum rules [5][6][7]. It should be noticed that the precision of extracted LCDAs depends on the accuracy of the involved hard kernels. As a by-product, the shape parameter for the B s meson distribution amplitude (DA) and the weak phase φ 3 (γ) = (75.2 ± 2.9) • consistent with the value in Particle Data Group [26] are also obtained. Though we have focused on the B → P P, V P decays at LO in PQCD, our work provides the first setup for a global analysis of exclusive B meson decays, and can be generalized straightforwardly to include other modes, and higher-order and/or higher-power corrections.
The rest of this paper is organized as follows. We give a brief overview of the theoretical framework for hadronic two-body B meson decays in Sec. II. The automatic derivation of the decay amplitudes in the PQCD approach is formulated in Sec. III, where the Gegenbauer-moment-independent database for the considered modes is established. We perform a global fit of meson LCDAs and the CKM angle φ 3 (γ) to a limited number of physical observables in the B → P P, V P decays, and present the numerical results in Sec. IV. We also compare our predictions for some other modes excluded in the fit with experimental data. A few remarks and future improvements on our analysis are outlined at the end of this section. Section V contains a summary of the present work. The explicit factorization formulas and their ingredients are collected in the Appendix.

II. THEORETICAL FORMALISM
In exclusive processes, such as heavy-to-light transition form factors, the range of a parton momentum fraction x, contrary to that in an inclusive case, is not experimentally controllable, and runs from 0 to 1. Hence, the endpoint region with x → 0 is unavoidable. If no endpoint singularity is developed, implying that the endpoint region is likely power suppressed, the collinear factorization will work. If such a singularity occurs, the collinear factorization will break down, and the k T factorization should be adopted. In fact, the observation QF 2 (Q 2 )/F 1 (Q 2 ) ∼ const. [27,28], F 1 and F 2 being the proton Dirac and Pauli form factors, respectively, and Q being a momentum transferred, indicates that the k T factorization is an appropriate tool for studying exclusive processes [29]. It has been shown that infrared divergences appearing in loop corrections to exclusive processes can be absorbed into hadron LCDAs in the k T factorization without breaking the gauge invariance [30]. Since the k T factorization theorem was proposed [31,32], there had been broad applications to various processes [33].
The application of the collinear factorization theorem to exclusive B meson decays, for instance, the B → π transition form factors, suffers the endpoint singularities mentioned above [34][35][36]: the twist-2 and twist-3 contributions are logarithmically and linearly divergent, respectively. The inclusion of parton transverse momenta k T , regulating the endpoint singularities, induces soft logarithms in higher-order corrections. Their overlap with the existent collinear logarithms generates the double logarithms α s ln 2 k T , which must be organized in order not to spoil perturbative expansion. The basic idea for the k T resummation of the double logarithms into a Sudakov factor has been elaborated in [15-17, 31, 37], where the explicit expressions of the Sudakov exponents can be found. The resultant Sudakov suppression on the low k T contribution in the endpoint region renders the magnitude of k 2 T roughly O(m b Λ QCD ). The coupling constant α s ( m b Λ QCD )/π ∼ 0.13 is then small enough to justify the perturbative evaluation of heavy-tolight transition form factors at large recoil [18,38,39].
On the other hand, the double logarithms α s ln 2 x from radiative corrections were observed in the semileptonic decay B → πlν [40] and in the radiative decay B → γlν [41]. It has been argued that when the endpoint region is important, these double logarithms should be organized into a quark jet function systematically in order to improve perturbative expansion. The procedure is referred to as the threshold resummation [42]. The resultant jet function has been shown to vanish quickly as x → 0. It turns out that in a self-consistent perturbative evaluation of the heavy-to-light transition form factors, where the original factorization formulas are further convoluted with the jet function, the endpoint singularities do not exist [42]. The threshold resummation for the jet function has been pushed to the next-to-leading-logarithm accuracy recently [43]. Note that either the threshold or k T resummation smears the endpoint singularities. To suppress the soft contribution sufficiently, both resummations are required, such that reliable results for the heavy-to-light transition form factors can be attained.
We emphasize that the power counting for a parton transverse momentum k T is nontrivial, compared to the power counting for the fixed scales like m b and Λ QCD . The k T factorization is suitable for a multi-scale process, like a heavy-to-light transition form factor, to which the region of a small momentum fraction x dominates. The small x introduces an additional intermediate scale A parton k T , being an integration variable in a k T factorization formula, can take values of orders of the above scales. The k T factorization should apply, as a hard kernel depends on the large scale m 2 b and the intermediate scale m b Λ QCD , but not on the small scale Λ 2 QCD , and the factorization of hadron wave functions hold for a parton k T at both the intermediate and small scales. Once these criteria are satisfied, the k T dependence in a hard kernel is not negligible [30], and a convolution between the hard kernel and the wave functions in k T is demanded. If a hard kernel involves only the large scale, the k T dependence of the hard kernel can be neglected. It is then integrated out in the wave functions, and one is led to the collinear factorization.
Since a wave function contains the contributions characterized by both the intermediate and small scales, it is legitimate to further factorize the former out of the wave function, as the intermediate scale is regarded as being perturbative. This gives the aforementioned k T resummation, which is justified perturbatively for the scale k 2 T ∼ m b Λ QCD . After this organization, the remaining piece, ie., the initial condition for the Sudakov resummation, involves only the small scale Λ 2 QCD , and corresponds to a hadron DA. Note that a more sophisticated formalism, called the joint resummation, which organizes the mixed logarithms formed by the above two different scales, has been developed in [44]. Similarly, it is also legitimate to further factorize the contribution characterized by an intermediate scale out of a hard kernel in the k T factorization. This re-factorization yields the jet function, through which the logarithms of xm 2 b are resummed to all orders. The effective Hamiltonian for hadronic two-body B meson decays is given by with the Fermi constant G F , the CKM matrix elements V qb(D) , D = d, s, the local four-quark operators O i , and the Wilson coefficients C i . All the factorizable, nonfactorizable and power-suppressed annihilation contributions resulting from the above four-quark operators are calculable in the PQCD approach without the endpoint singularities. The arbitrary cutoffs introduced in QCDF [45,46] are not necessary, and PQCD factorization formulas involve only universal and controllable inputs. The B → M 2 M 3 decay amplitude is generically factorized into the convolution of the Wilson coefficient C, a six-quark hard kernel H, the jet function J t , and the Sudakov factor S with meson LCDAs φ [47][48][49][50], all of which are well defined and gauge invariant. The partition of nonperturbative and perturbative contributions depends on factorization schemes. However, a decay amplitude, as a convolution of the above factors, is independent of factorization schemes in principle.

A. Lightcone Distribution Amplitudes
The momenta p B , p 2 and p 3 of the B meson, emitted meson M 2 , and recoiling meson M 3 , respectively, and their associated parton momenta are chosen, in the light-cone coordinates, as which are labelled in Fig. 1 with m B being the B meson mass and x i being the momentum fractions. The light meson LCDAs are defined through the matrix elements where N c = 3 is the number of colors, m 0 is the chiral enhancement scale for the pseudoscalar meson P , the dimensionless vector v = √ 2p/M B lies in the direction of the meson momentum p, the dimensionless vector n lies in the direction of the quark coordinate z with n · v = 1, and m V (ǫ L ) is the mass (longitudinal polarization vector) of the vector meson V . The light meson LCDAs are expanded as in terms of the orthogonal Gegenbauer polynomials where the decay constants f P , f V and f T V can be extracted from leptonic decay widths such as Γ(π → µν) and Γ(τ → ρν), and the superscripts f of the Gegenbauer moments label the species of mesons.
The B meson DA is defined via the matrix element with the light spectator momentum k, and the dimensionless vectors n + = (1, 0, 0 T ) and n − = (0, 1, 0 T ). In this work we adopt the model for the B (s) meson DA, , the shape parameter ω B (s) will be determined in the next section, and b is the impact parameter conjugate to the transverse momentum k T . It has been argued [38] that the contribution fromφ B (k) is power suppressed, so it will be neglected in the numerical analysis below.

B. SU(3) Flavor Structure
To calculate hadronic two-body B meson decay amplitudes systematically, we introduce the following SU(3) matrix elements for various species of mesons, which reflect the internal structure of the flavor SU(3) group. The isosinglet mesons like η q , η s , ω, and φ will not be considered in the global analysis below, but their properties are listed here for completeness. We will take into account these hadrons, as extending the database in the future. The matrices relevant for the heavy-to-light transitions are given by The factorization formula for a B → P P decay amplitude in Eq. (2)  We decompose the total B → M 2 M 3 decay amplitude into the combination where A u (B → M 2 M 3 ) denotes the tree contribution with the product V ub V * uf of the CKM matrix elements, and A t (B → M 2 M 3 ) denotes the penguin contribution with the product V tb V * tf . These amplitudes are written, in terms of the matrices in Eqs. (9) and (10), as The B → V P decay amplitudes can be simply inferred from the B → P P amplitudes through the replacements of light meson LCDAs and of a chiral enhancement scale by a vector meson mass. For the B → V 2 P 3 emission and B → P 2 V 3 annihilation, we apply the rule φ  , and further flip the signs of the LR and SP amplitudes, where the subscript (2) means twist 2, and (3) means twist 3. For the B → V 2 P 3 annihilation, we apply φ 3 , and further flip the signs of the LR and SP amplitudes. The above rule holds for both the pseudoscalar P and vector V mesons, and for the factorizable and nonfactorizable diagrams. For the B → P 2 V 3 emission, we apply φ 3 , and further flip the signs of the nonfactorizable LR amplitudes and the factorizable SP amplitudes.
As shown in Eq. (5), there are 9 Gegenbauer moments a f 's in total for the twist-2 pseudoscalar LCDA φ P (x) and twist-3 LCDAs φ P P (x) and φ T P (x). Note that the Gegenbauer moment a π 1 vanishes due to the isospin symmetry, and a f P 4 are not included in the fit, because they cannot be constrained effectively under the current limited experimental accuracy. Thus, a B → P P decay amplitude contains 9 × 9 combinations of the Gegenbauer moments, where the product of the Gegenbauer moments a f n a f m has been factored out explicitly. We compute the factorization formula M nm , which involves only the Gegenbauer polynomials associated with a f n and a f m , to establish a 9 × 9 database. Each database has 20 sets of values, corresponding to the Wilson coefficients a 1 · · · a 10 and C 1 · · · C 10 in Eq. (12). To analyze the B → V P decays, we construct a 9 × 4 × 2 database for M nm in a similar manner, where ×2 is attributed to the two possible final states P 2 V 3 and V 2 P 3 . The inputs for the Fermi constant, the meson decay constants, the meson masses, and the chiral enhancement scale are the same as in [51], and the magnitudes of the CKM matrix elements are referred to [26] in the above computations.

A. Least-Square Fit and Bayesian Analysis
We determine the Gegenbauer moments and the weak phase φ 3 (γ) by fitting the branching ratios and direct CP asymmetries formulated from the decay amplitudes in Eq. (13) to experimental data using the nonlinear least-χ 2 (lsq) method [52]. The lsq method minimizes the summed residual S, where r i is the residual at the i-th point x i , y i (ŷ i ) represents the experimental data (response value), and n is the number of data points. One defines a model functionŷ i = f (x i , β), where the vector β contains the m adjustable parameters considered in the fit. The minimum of Eq. (14) is obtained by equating the gradient to zero, For a linear model, f (x i , β) can be decomposed into a sum of multiple linear functions G j (x i ) with the corresponding coefficients β j , f (x i , β) = j β j G j (x i ). One then regards the function G j (x i ) as a matrix, and solves for β j from Eq. (15) with the data y i .
A non-linear model is more subtle, to which there is no closed-form solution in general. A possible approach is to select some initial values of the parameters, and refine the parameters by iteration. At each iteration, the model function f (x i , β) is linearized through the first-order Taylor series expansion at β k j , with k being an iteration number, and J being a Jacobian function. The minimum of the residual at this iteration is achieved by equating the gradient to zero, where ∆β j can be solved by inverting the Jacobian matrix. The parameters then take the values β k+1 j = β k j + ∆β j for the next iteration.
Because of the approximations in the Taylor expansion and in the matrix inversion (if the Jacobian matrix is not a square one), no algorithm works for all nonlinear models, and fit results may be sensitive to initial conditions. To stabilize a complicated nonlinear fit, one can perform a Bayesian analysis (conditionally biased fit). Instead of the single χ 2 term corresponding to Eq. (14), with δy i being the errors of experimental data, we employ a modified version The second term χ 2 prior is a stabilizing function, whereβ j are artificially introduced default values for the fitted parameters with some errorsσ j . When the summed residual in Eq. (19) are minimized, the χ 2 prior term favors β j in the range (β j −σ j ,β j +σ j ). The values ofβ j andσ j should be chosen reasonably according to prior physical knowledge on the fitted parameters. In the present study we select the Gegenbauer moments derived in QCD sum rules as our Bayesian dataβ, and five times of the errors from QCD sum rules for ourσ, which reduce the weight of χ 2 prior in the fit.

B. Bs Meson Distribution Amplitude
We point out that the B (s) meson DA appears in all the B (s) → M 2 M 3 decay amplitudes, so it is difficult to extract the shape parameters ω B (s) in a global fit. The investigation in [38] shows that ω B = 0.4 GeV for the B meson DA leads to reasonable results for the B → π transition form factors, which agree with those from light-cone sum rules and lattice QCD. Therefore, we choose this value as the input, and perform the global fit to determine the Gegenbauer moments of the light meson LCDAs. As to the shape parameter ω Bs in the B s meson DA, it is observed that the Gegenbauer moments fitted from the B s → P P data are not sensitive to its variation: we scan the range of ω Bs from 0.4 GeV to 0.6 GeV, and make sure that the fitted Gegenbauer moments are relatively stable in the range 0.45 GeV < ω Bs < 0.55 GeV. We then select the data of several precisely measured channels, B s → K + K − , K + π − , K 0K 0 , and π + π − , and compare them with the reconstructed data from the PQCD factorization formulas. The comparison displayed in Fig. 4, where the red bands are from the experimental data and the blue bands are from the reconstructed data, indicates that the choice ω Bs = 0.48 GeV is preferred: with this value of ω Bs , the PQCD results from the fit accommodate the selected four piece of data simultaneously.

C. Global Fit
We fit the PQCD factorization formulas with the database constructed in the previous section to the measured branching ratios and direct CP asymmetries A CP in the B (s) → P P, V P decays, which are collected in the left columns of Table III. The A CP data marked in red, which have larger errors, do not provide a strong constraint in the fit. Note that the LHCb Collaboration has updated their measurement of A CP in the B 0 →K 0 π 0 mode, which reads −13.8 ± 2.5 [53]. The data of those modes, which are greatly affected by subleading contributions according to the existent PQCD calculations [51,54,55], namely, suffer significant theoretical uncertainties, are excluded in our fit. The Gegenbauer moments of both the twist-2 and twist-3 LCDAs from a joint fit, corresponding to the shape parameters w B = 0.4 GeV and w Bs = 0.48 GeV, are listed in Tables I with χ 2 /d.o.f. = 0.77. The errors in our fit mainly arise from the experimental uncertainty. The χ 2 prior term in the Bayesian analysis introduces little error to the fit results. Some higher-order moments or moments of higher-twist LCDAs cannot be constrained effectively due to the current limited experimental accuracy. This is the reason why the values of the moments a K P 2 , a K T 2 and a K * 1 are not presented in Table I. It is seen that some fitted Gegenbauer moments, like a ρ 2 , a K 2 and a K * 2 , agree well with those from QCD sum rules [5,7] within 1σ error, which are listed in Table II for comparison. We stress that our fit is based on the LO PQCD factorization formulas, and that next-to-leading-order (NLO) corrections change the heavy-to-light transition form factors by about 30% [56,57]. It is difficult to estimate how much systematic error is caused by NLO effects for the fitting at LO, because NLO corrections to the non-factorizable amplitudes in hadronic two-body B meson decays have not yet been completed in the PQCD approach. Higher-power contributions (for example, the power-suppressed contribution from another B meson DAφ B in Eq. (8) was shown to be of the same order as the NLO one in [58]) have not been taken into account either. Therefore, it is likely that some fitted Gegenbauer moments, such as a π 2 , differ more significantly from those in QCD sum rules. The outcome of a π P 2 , slightly larger than unity, can be reduced by including the higher moment a π P 4 into the fit. As explained before, a π P 4 is not considered here, because it cannot be constrained effectively under the current experimental accuracy. It is worth mentioning that the weak phase φ 3 (γ) is found to be (75.2 ± 2.9) • , consistent with the value (72.1 +4.1 −4.5 ) • in Particle Data Group [26], and (69.8 ± 2.1 ± 0.9) • from the factorization-assisted topological diagram approach [60]. The agreements of our results with the Gegenbauer moments from sum rules and with φ 3 (γ) extracted in other methods support the PQCD factorization for hadronic two-body B meson decays. With the fitted Gegenbauer moments in Table I, we calculate the branching ratios and A CP in the LO PQCD approach, and present the results in the right columns of Table III. It is observed that all the considered data, except the B − → π 0 K * − branching ratio, are well reproduced. The observables removed from the fit, i.e., those suffering significant theoretical uncertainties, are also predicted in the LO PQCD formalism, and compared with the data in  Table IV. The predicted branching ratios are very close to the values obtained in the previous PQCD calculations, so the deviation from the data remains. In particular, A CP in the B − → π − ρ 0 mode has been predicted to be large and negative in most QCD approaches [61,62], but its data are as small as 0.009 ± 0.019 [26]. The inclusion of higher-order and higher-power contributions to hadronic two-body Bmeson decays may improve the consistency. The theoretical errors given in Tables III and IV arise only from those of the fitted Gegenbauer moments and φ 3 (γ). More precise measurements are urged, and subleading contributions should be included into the PQCD framework to strengthen the constraint on the Gegenbauer moments and to sharpen the confrontation between theoretical predictions and experimental data.

D. Remarks and Future Developments
A few remarks are given as follows.
• We have focused only on the branching ratios and direct CP asymmetries in the B → P P, V P decays, and neglected those modes involving isosinglet mesons in the above analysis. Other observables, such as mixing- Bs → π − K * + 2.9 ± 1.1 − 9.53 ± 0.24 −25.5 ± 1.0 7.6 [51] induced CP asymmetries and polarizations in B → V V decays, can be included straightforwardly. Though more parameters will be introduced through LCDAs for transversely polarized vector mesons, sufficient precise measurements on polarization observables can be achieved at LHCb and Belle-II.
• LCDAs also appear in the factorization formulas for heavy-to-light transition form factors that govern semileptonic B meson decays. One can take into account experimental constraints from these decays in the future, in particular those from their dependence on the lepton-pair invariant mass squared q 2 .
• As we have pointed out, decay widths of a few modes are suppressed at LO in PQCD, and may be well described with the inclusion of higher-order contributions [54,63]. Some sources of power corrections have been explored in Refs. [64,65]. A new database will be established in a similar way by using the flavor structure for these radiative and power corrections, via which the precision of a global analysis can be enhanced.
• We did not consider all systematic and parametric uncertainties in the current analysis, such as the ones originating from the variations of factorization scales and nonperturbative QCD parameters.
• If a high-precision global study reveals notable tensions between theoretical results and experimental data in the future, it may hint that new physics effects are inevitable. One is then motivated to include new physics contributions, which can also be analyzed according to the flavor structure of new physics operators.

V. SUMMARY
As stated in the Introduction, nonperturbative hadron LCDAs provide a major source of theoretical uncertainties in all the factorization-based approaches to hadronic two-body B meson decays. In this paper we have performed a global analysis of the light meson LCDAs by fitting the LO PQCD factorization formulas for B (s) → P P, V P decays to available data of branching ratios and direct CP asymmetries. A computation code was developed based on the flavor structure of the four-quark effective operators to establish the database, which contains the part of decay amplitudes without the Gegenbauer moments. This database facilitates the global fit, from which the Gegenbauer moments of the twist-2 and twist-3 LCDAs for the pseudoscalar meson P (P = π, K) and vector meson V (V = ρ, K * ) were determined. Most of our fit results agree with the moments derived in QCD sum rules, and those with discrepancies deserve more thorough investigation that takes into account higher-order and higher-power corrections in the PQCD approach. The weak phase φ 3 (γ) = (75.2 ± 2.9) • in consistency with the value in Particle Data Group was also extracted. Predictions for the modes, which were excluded in the fit due to large theoretical uncertainties, are close to the existent PQCD results, and still deviate from the data. To improve the consistency, subleading contributions to hadronic two-body B meson decays need to be included, when their evaluation is completed in the future.
Since the B (s) meson DA appears in all the B (s) → M 2 M 3 decay amplitudes, it is difficult to constrain the shape parameters ω B (s) in this DA in a global fit. The shape parameter ω B = 0.4 GeV is an input, and ω Bs = 0.48 GeV is subject to a discretionary choice in the present study. The difficulty is expected to be overcome, when data for exclusive processes other than hadronic two-body B meson decays are considered in the global fit. This is a straightforward extension of the framework proposed here, through which the global determination of LCDAs for other hadrons is also feasible. The explicit LO PQCD factorization formulas for the B → P P decay amplitudes from the various current operators and topologies are presented in this appendix, with C F = 4/3, the Wilson coefficients a i , and r i = m 0i /m B , where m 0i is the chiral enhancement scale: The hard scales involved in the above decay amplitudes are defined by The hard kernels h in the decay amplitudes consist of two parts, the jet function J t (x i ) derived in the threshold resummation and the Fourier transformation of the virtual particle propagators: with the Bessel function H (1) 0 (z) = J 0 (z) + i Y 0 (z). The following approximate parametrization for the jet function has been proposed for convenience [38], with the parameter c ≈ 0.3. The prefactor in the above expression is chosen to obey the normalization 1 0 J t (x)dx = 1. The jet function J t (x) gives a very small numerical effect to the nonfactorizable amplitude [42], so it is dropped from h n and h na .
The evolution factors E in which the Sudakov exponents are given by with the quark anomalous dimension γ q = −α s /π. Replacing the kinematic variables of M 2 by those of M 3 in S 2 , we get the expression for S 3 . The function s(Q, b) is expressed as where n f is the number of the quark flavors and γ E is the Euler constant. We adopt the one-loop running coupling constant, so only the first four terms of Eq. (A24) are picked up in the numerical analysis.