Prompt charmonia production and polarization at LHC in the NRQCD with kt-factorization. Part I: psi(2S) meson

In the framework of kt-factorization approach, the production and polarization of prompt psi(2S) mesons in pp collisions at the LHC energies is studied. Our consideration is based on the non-relativistic QCD formalism for bound states and off-shell amplitudes for hard partonic subprocesses. The transverse momentum dependent (unintegrated) gluon densities in a proton were derived from Ciafaloni-Catani-Fiorani-Marchesini evolution equation or, alternatively, were chosen in accordance with Kimber-Martin-Ryskin prescription. The non-perturbative color-octet matrix elements were first deduced from the fits to the latest CMS data on psi(2S) transverse momentum distributions and then applied to describe the ATLAS and LHCb data on psi(2S) production and polarization at 7 TeV. We perform the estimation of polarization parameters $\lambda_\theta$, $\lambda_\phi$ and $\lambda_{\theta \phi}$ which determine psi(2S) spin density matrix and demonstrate that taking into account the off-shellness of initial gluons in the color-octet contributions leads to unpolarized psi(2S) production at high transverse momenta, that is in qualitative agreement with the LHC data.


Introduction
The production of quarkonium states in high energy hadronic collisions is under intense theoretical and experimental study [1][2][3] since two decades ago, when the measurements of prompt J/ψ and Υ production cross sections at the Tevatron revealed a more than one order-of-magnitude discrepancy with theoretical expectations obtained in the framework of the color singlet model [4]. This fact had induced extensive theoretical activity mainly connected with modeling the formation of quarkonium states from unbound heavy quark pairs produced in hard interation. There exist two competing theoretical approaches known in the literature under the names of color-singlet (CS) and color-octet (CO) models. In general, a quark-antiquark pair is produced in a state 2S+1 L (a) J with spin S, orbital angular momentum L, total angular momentum J and color a, which can be either identical to the respective quantum numbers of the resulting quarkonium (as accepted in the CS model) or different from those. In the latter case, the heavy quark pair transforms into physical quarkonium state by means of soft (nonperturbative) gluon radiation, as considered in the formalism of non-relativistic Quantum Chromodynamics (NRQCD) [5,6]. The quarkonium formation probability is then determined by the respective nonperturbative matrix elements (NMEs), which are assumed to be universal (process-independent) and not depending on the quarkonium momentum. Though not strictly calculable within the theory, the NMEs are assumed to obey certain hierarchy in powers of the relative quark velocity v. To the leading order in v, an S-wave vector meson such as ψ(2S) can be formed from a quark pair produced as color singlet 3 S with J = 0, 1 or 2. We know already that the CS model with leading order (LO) hard scattering matrix elements fails to describe the experimental data on J/ψ and Υ production at the Tevatron and LHC. Including the next-to-leading order (NLO) [7] and dominant tree-level next-tonext-to-leading order (NNLO * ) [8] corrections to the CS mechanism significantly improves the description of the collider data [9]. In the NRQCD formalism, a reasonably good agreement with the measured quarkonia cross sections can be achieved by adjusting the NME values which play the role of free fitting parameters [10][11][12][13][14][15]. This was already demonstrated by comparing the calculations with the ATLAS [16], CMS [17] and LHCb [18] experimental data taken at √ s = 7 TeV. However, the values of the extracted NMEs dramatically depend on the minimal quarkonium transverse momentum p T used in the fits [19] and are incompatible with each other (both in size and even in sign!) when obtained from fitting the different sets of data. Furthermore, none of the fits is able to accommodate the polarization measurements [20,21]. The fits involving low p T measurements lead to the conclusion that the production of S-wave quarkonia at high p T must be dominated by CO contributions with transverse polarization (namely, by the 3 S channel). The latter contradicts to the unpolarized production seen by the CDF [22,23] Collaboration at the Tevatron 1 and CMS [20,24] and LHCb [21] Collaborations at the LHC. To obtain an unpolarized state it is necessary to either assume that quarkonium production is dominated by the scalar 1 S (8) 0 channel [12] or restrict the NRQCD fits to very high p T region [19]. This problem is known as "quarkonium polarization puzzle" and still far from understanding.
A new solution to the polarization puzzle has been guessed in the k T -factorization approach of QCD [25,26], where studies of quarkonia production and polarization have their own long history (see, for example, [27][28][29][30][31][32][33][34] and references therein). The k T -factorization approach is based on the Balitsky-Fadin-Kuraev-Lipatov (BFKL) [35] or Ciafaloni-Catani-Fiorani-Marchesini (CCFM) [36] gluon evolution equations and provides theoretical grounds for including the effects of initial gluon radiation and intrinsic gluon transverse momentum 2 . A combination of the usual CS scheme and the k T -factorization formalism results in reasonably good description of the data on J/ψ, χ c and Υ production at HERA [27,29,30,33], Tevatron [28,32] and LHC [34]. These results are also in good agreement with more complicated explicit NNLO * calculations [7][8][9] performed in the collinear approximation of QCD. The longitudinal polarization of directly produced quarkonia in the k T -factorization is an immediate consequence of initial gluon off-shellness [27]. Adding the feed-down from P -wave states (χ c and χ b ) leads to essentially unpolarized prompt J/ψ and Υ mesons [32,34]. In the k T -factorization approach at LO, the P -wave states are produced in 2 → 1 gluon-gluon fusion and are expected to dominate at high transverse momenta. However, the latest LHC data show that the feed-down contributions χ c → J/ψ+γ and χ b → Υ+γ do not constitute more than only 20 − 30% of the visible cross section at large p T values.
So, the production of χ c and χ b mesons requires a dedicated study which will be the subject of our forthcoming papers, while in the present analysis we concentrate on the direct mechanism and only restrict to ψ(2S) mesons having no contamination from higher states. Here we present a systematic analysis of ATLAS [16], CMS [17] and LHCb [18] data collected at √ s = 7 TeV regarding the transverse momentum distributions and polarization parameters λ θ , λ φ and λ θφ which describe the spin density matrix of the produced ψ(2S) mesons. The outline of our paper is the following. In Section 2 we briefly recall the NRQCD formalism and the k T -factorization approach. In Section 3 we perform a numerical fit to the latest CMS data and extract the color-octet NMEs using three different sets of transverse momentum dependent (TMD) gluon distributions. Later in this section we check the compatibility of the extracted papameters with ATLAS and LHCb data on the production and polarization of ψ(2S) mesons. The comparison is followed by a discussion. Our conclusions are collected in Section 4.

Theoretical framework
We start with briefly recalling the essential calculation steps. Our consideration is based on the following leading-order off-shell partonic subprocesses [5,6]: where J = 0, 1 or 2, and the four-momenta of all particles are indicated in parentheses. The subprocesses (1) and (2) represent the CS and CO contributions, respectively. The corresponding production amplitudes can be obtained from the one for an unspecified cc state by applying the appropriate projection operators, which guarantee the proper quantum numbers of the cc state under consideration. These operators for the different spin and orbital angular momentum states read [4]: where m = 2m c is the mass of the considered cc state, and p c and pc are the four-momenta of the charmed quark and anti-quark. States with various projections of the spin momentum onto the z axis are represented by the polarization four-vector ǫ µ (S z ). The probability for the two quarks to form a meson depends on the bound state CS and fictituous CO wave functions Ψ (a) (q), where the relative four-momentum q of the quarks in the bound state is treated as a small quantity in the non-relativistic approximation. So, we represent the quark momenta as Then, we multiply the hard subprocess amplitude A (depending on q) by the meson wave function Ψ (a) (q) and perform integration with respect to q. The integration is done after expanding the amplitude A around q = 0: Since the expressions for A| q=0 and ∂A/∂q α | q=0 are no longer dependent on q, they may be factored outside the integral sign. A term-by-term integration of this series employs the identities [4]: where R (a) (x) are the radial wave functions in the coordinate representation.The first term in (7) contributes only to S-waves, but vanishes for P -waves. On the contrary, the second term contributes only to P -waves, but vanishes for S-waves. States with various projections of the orbital angular momentum onto the z axis are represented by the polarization fourvector ǫ µ (L z ). The NMEs of S-wave states are directly related to the CS and fictituous CO wave functions: where N c = 3 and J = 1. A similar relation holds for R ′ (a) if P -wave states are involved. The CS wave function at the origin of coordinate space is known from the measured ψ(2S) leptonic decay width. The color-octet NMEs are extracted from experimental data and obey the relation coming from heavy quark spin symmetry at LO in v. The polarization vectors ǫ µ (S z ) and ǫ µ (L z ) are defined as explicit four-vectors. In the frame where the z axis is oriented along the quarkonium momentum vector p µ = (E, 0, 0, |p|), these polarization vectors read The states with definite S z and L z are translated into states with definite total momentum J and its projection J z using the Clebsch-Gordan coefficients: Further evaluation of partonic amplitudes is straightforward and is done using the algebraic manipulation system form [38]. Our results for perturbative production amplitudes squared and summed over polarization states agree with ones [39]. Here we only mention several technical points. First, according to the k T -factorization prescription [25,26], the summation over the incoming off-shell gluon polarizations is done using the gluon spin density matrix in the form ǫ µ ǫ * ν = k µ T k ν T /k 2 T , where k T is the gluon transverse momentum orthogonal to the beam axis. In the collinear QCD limit, when |k T | → 0, this expression converges to the ordinary ǫ µ ǫ * ν = −g µν /2 after averaging over the azimuthal angle. Second, the ψ(2S) spin density matrix is expressed in terms of the momenta l 1 and l 2 of the decay leptons and is taken as This expression is equivalent to the standard one ǫ µ ǫ * ν = −g µν + p µ p ν /m 2 but is better suited for studying the polarization observables because it gives access to the kinematic variables describing the orientation of the decay plane. In all other respects the evaluation follows the standard QCD Feynman rules. An important point in the NRQCD formalism is connected with the emission of soft gluons taking place after the hard interaction is over. It is usually assumed that the emitted soft gluons bring away the unwanted color and change other quantum numbers of the produced CO system, but do not carry any energy, thus keeping the kinematics intact. However, such an emission contradicts to the basic QCD property: soft gluons can never be radiated as they are confined. In order that the quantum numbers get changed, one needs to radiate a real gluon with some energy E ∼ Λ QCD , giving us the confidence that we do not enter into the confinement or perturbative domains. When considering the gluon radiation from 3 P states, we rely upon the dominance of electric dipole E1 transitions, which is supported by the E835 experimental data [40]. The corresponding invariant amplitudes can be written as follows [41]: where p (CO) µν are the four-momenta and polarization fourvectors (tensor) of corresponding particles and e µναβ is the fully antisymmetric Levi-Civita tensor. The gluon radiation from other CO states is generated according to the phase space.
The cross section of ψ(2S) production at high energies in the k T -factorization approach is calculated as a convolution of the off-shell partonic cross sections and the TMD gluon densities in a proton. The contribution from the CS production mechanism can be presented in the following form: where f g (x, k 2 T , µ 2 ) is the TMD gluon density, p T and y are the transverse momentum and rapidity of produced ψ(2S) meson, y g is the rapidity of outgoing gluon and √ s is the pp center-of-mass energy. The initial off-shell gluons have a fraction x 1 and x 2 of the parent protons longitudinal momenta, non-zero transverse momenta k 1T and k 2T (k 2 1T = −k 2 1T = 0, k 2 2T = −k 2 2T = 0) and azimuthal angles φ 1 and φ 2 . For the CO production we have: According to the general definition [42], the off-shell gluon flux factor in (18) is defined 3 The squares of the corresponding off-shell partonic amplitudes, as being too lengthy, are not presented there but the full C++ code is available on request 4 . The multidimensional integration have been performed by means of the Monte Carlo technique, using the routine vegas [43].
Numerically, we have tested several different sets of the TMD gluon densities. Two of them (namely, JH and A0 sets) have been obtained [44,45] from the CCFM equation where all input parameters have been fitted to describe the proton structure function F 2 (x, Q 2 ). Besides the CCFM-evolved gluon densities, we applied the one obtained from the Kimber-Martin-Ryskin (KMR) prescription [46]. The KMR approach is a formalism to construct the TMD quark and gluon distributions from well-known conventional ones. For the input, we have used leading-order Martin-Stirling-Thorn-Watt (MSTW'2008) set [47].

Numerical results
We now are in a position to present our numerical results. First we describe our input and the kinematic conditions. Having the TMD gluon distributions chosen, the cross sections (17) and (18) depend on the renormalization and factorization scales µ R and µ F . We set µ 2 R = m 2 +p 2 T and µ 2 F =ŝ+Q 2 T , where Q T is the transverse momentum of the initial off-shell gluon pair. The choice of µ R is the standard one for studying the charmonia production whereas the special choice of µ F is connected with the CCFM evolution [44,45]. Following [48], we set ψ(2S) mass m = 3.686 GeV, branching fraction B(ψ ′ → µ + µ − ) = 0.0077 and use the LO   [17]. The NMEs obtained in the NLO NRQCD fits [11,15] are shown for comparison.
formula for the coupling constant α s (µ 2 ) with n f = 4 quark flavours at Λ QCD = 200 MeV, such that α s (M 2 Z ) = 0.1232. In Table 1 we list our results for the NMEs fits obtained for three different TMD gluon densities. We have fitted the transverse momentum distributions for prompt ψ(2S) mesons measured recently by the CMS Collaboration at √ s = 7 TeV [17]. These measurements were done at moderate and high transverse momenta 10 < p T < 100 GeV, where the NRQCD formalism is believed to be most reliable. In contrast with [11,15], we performed the fitting procedure under requirement that the NME values be positive only. The color singlet NMEs were not fitted, but just taken from the known ψ(2S) → µ + µ − partial decay width [48]. For comparison, we also present in Table 1 two sets of NMEs, obtained within the NLO NRQCD in [11,15]. The main difference betwen [11] and [15] is in that these fits were based on differently selected sets of data points.
In the k T -factorization approach, the fitted NME values strongly depend on the choice of TMD gluon density. We find that the 1 S is very close to the one obtained in the NLO NRQCD analysis [11]. Both the NLO NRQCD fits [11,15] significantly (by one order of magnitude) exceeds the values of O ψ 3 S , obtained with all of the TMD gluon densities. It is almost consistent with estimates performed by other authors [49,50]. In contrast with the results [39], our fit leads to non-zero O ψ 3 P values. Summing up, we can conclude that the NME values obtained by the different authors on the basis of different data sets or by the same authors using different gluon densities are widely spread, that spoils the belief in the universality of the matrix elements. Now we turn to comparing our predictions with the data collected by the ATLAS [16], CMS [17] and LHCb [18] Collaborations. The ATLAS Collaboration has measured prompt ψ(2S) transverse momentum distribution at 10 < p T < 100 GeV at central rapidities |y| < 2 [16]. The CMS Collaboration probes the transverse momentum in the range 10 < p T < 100 GeV at |y| < 1.2 [17], and the LHCb Collaboration works in the kinematic range p T < 16 GeV and 2 < y < 4.5 [18]. In all cases the data were obtained at √ s = 7 TeV. The results of our calculations are shown in Figs. 1 -3. With our sets of NMEs, we achieve reasonably good description of the data with any of the considered TMD gluon density. We observe dominance of the CO contributions in the whole p T range. In particular, the 3 S 1 contribution dominates at the large p T > 25 GeV, whereas the 1 S 0 channel is mostly important at low p T values. Taken solely, the CS contributions (even incorporated with the k T -factorization) are unable to describe the data. They are important at relatively low p T only and are comparable there with the 3 P (8) J contributions. At moderate and high transverse momenta the CS contributions are below the data by about one order-of-magnitude. The predictions obtained with the chosen TMD gluon densities are very close to each other at p T > 6 GeV, while the difference becomes only sizable at low p T < 6 GeV (see Fig. 3). Therefore, similar to the collinear QCD factorization, including the low p T data to the fit procedure can change the relative weight of different NMEs in the k T -factorization approach, that is out of our present study. Now we turn to the the ψ(2S) polarization issue, which is the most interesting part of our study. In general, the spin density matrix of a vector particle decaying into a lepton pair depends on three angular parameters λ θ , λ φ and λ θφ which can be measured experimentally. The double differential angular distribution of the decay leptons can be written as [51]: where θ * and φ * are the polar and azimuthal angles of the decay lepton measured in the charmonium rest frame. The case of (λ θ , λ φ , λ θφ ) = (0, 0, 0) corresponds to unpolarized state, while (λ θ , λ φ , λ θφ ) = (1, 0, 0) and (λ θ , λ φ , λ θφ ) = (−1, 0, 0) refer to fully transverse and longitudinal polarizations. The CMS [20] and LHCb [21] Collaborations have measured all these parameters as functions of ψ(2S) transverse momentum in two complementary frames: the Collins-Soper and helicity ones. In addition, the CMS Collaboration provided measurements in the perpendicular helicity frame. In the Collins-Soper frame the polarization axis z bisects the two beam directions whereas the polarization axis in the helicity frame coincides with the direction of ψ(2S) momentum in the laboratory frame. In the perpendicular helicity frame the z axis is orthogonal to that in the Collins-Soper frame and lies in the plane spanned by the two beam momenta. Additionally, the frame-independent polarization parameter [51,52] λ * = (λ θ + 3λ φ )/(1 − λ φ ) was investigated. Below we estimate the polarization parameters λ θ , λ φ , λ θφ and λ * for the CMS and LHCb conditions. Our calculation generally follows the experimental procedure. We collect the simulated events in the kinematical region defined by the CMS and LHCb experiments, generate the decay lepton angular distributions according to the production and decay matrix elements, and then apply a three-parametric fit based on (19). In Figs. 4 -8 we confront our predictions for polarization parameters λ θ , λ φ , λ θφ and λ * with the latest CMS [20] and LHCb [21] data. We find slight transverse polarization (λ θ ∼ 0.2) in the Collins-Soper frame and slight longitudinal polarization (λ θ ∼ −0.2) in the helicity frame at low transverse momenta covered by the LHCb experiment. These results are practically independent on the ψ(2S) rapidity. At higher p T the polarizations of ψ(2S) mesons, calculated in the Collins-Soper frame, go from slight transverse (λ θ ∼ 0.15) to almost zero values (λ θ ∼ 0.05) as the transverse momentum increases from p T ∼ 10 GeV to 50 GeV. In the helicity and perpendicular helicity frames ψ(2S) polarization changes from longitudinal (λ θ ∼ −0.2) to slight longitudinal (λ θ ∼ −0.1). Here we arrive at the key point of our paper. Figs. 4 -8 show that ψ(2S) production, calculated in the k T -factorization approach, tends to be unpolarized at high p T , in agreement with the CMS data [20]. Indeed, as a strict consequence of the initial gluon off-shellness, a large fraction of ψ(2S) mesons with zero helicity is produced in the partonic subprocesses, including both CS and CO contributions. Moreover, the fraction of such events increases when p T grows up (see, for example, [27][28][29][30][31][32][33][34]). The only exception refers to 1 S (8) 0 channel, which is produced unpolarized due to its spinless nature. It is a remarkable property of the k T -factorization scheme that the gluon fragmentation to 3 S 1 states produces nearly unpolarized ψ(2S) mesons (in contrast with conventional collinear NRQCD where the mesons carry strong transverse polarization). Thus, we can conclude that the problem of ψ(2S) spin alignment can be solved if the initial gluon off-shellness is taken into account.
A comparison of our predictions with the LHC data [20,21] shows that the latter seem to support the trend observed in the k T -factorization formalism. However, while our predictions for λ φ and λ θφ parameters agree with the data, the description of λ θ and λ * is still rather qualitative than quantitative, due to the huge experimental uncertainties.
Finally, we would like note that there are significant theoretical uncertainties connected with the choice of the renormalization and/or factorization scales, the inclusion of NLO subprocesses and exact definition of NMEs. The detailed study of these uncertainties is out of our present paper.

Conclusions
We have considered prompt ψ(2S) production and polarization in pp collisions at the LHC energy √ s = 7 TeV in the framework of k T -factorization approach. We have used the LO nonrelativistic QCD formalism including both color-singlet and color-octet contributions. Using the TMD gluon densities in a proton derived from the CCFM equation and from the Kimber-Martin-Ryskin prescription, we extracted the color-octet NMEs O ψ 1 S for ψ(2S) mesons from fits to transverse momentum distributions provided by the latest CMS measurements. Using the fitted NMEs, we have successfully described the data presented by the ATLAS, CMS and LHCb Collaborations. We estimated the polarization parameters λ θ , λ φ and λ θφ and demonstrated that taking into account the offshellness of the initial gluons in the color-octet contributions leads to unpolarized ψ(2S) production at high transverse momenta, that is in qualitative agreement with the LHC data.