Charmonia production from $b$-hadron decays at LHC with $k_T$-factorization: $J/\psi$, $\psi(2S)$ and $J/\psi + Z$

We consider the production of $J/\psi$ and $\psi(2S)$ mesons originating from the decays of $b$-flavored hadrons at the LHC using the $k_T$-factorization approach. Our analysis covers both inclusive charmonia production and production of $J/\psi$ mesons in association with $Z$ bosons. We apply the transverse momentum dependent (or unintegrated) gluon density in a proton derived from Catani-Ciafaloni-Fiorani-Marchesini (CCFM) evolution equation and adopt fragmentation functions based on the non-relativistic QCD factorization to describe the inclusive $b$-hadron decays into the different charmonium states. Our predictions agree well with latest experimental data taken by the CMS, ATLAS and LHCb Collaborations at $\sqrt s = 7$, $8$ and $13$~TeV. The contributions from double parton scattering to the associated non-prompt $J/\psi + Z$ production are estimated and found to be small.


Introduction
A few years ago, the CMS, ATLAS and LHCb Collaborations reported the measurements of the inclusive production of J/ψ and ψ(2S) mesons originating from b-hadron decays (so called "non-prompt" production) at the LHC energies √ s = 7 and 8 TeV [1][2][3][4]. More recently, the LHCb Collaboration has measured [5] the total and differential cross sections of the inclusive non-prompt J/ψ production at √ s = 13 TeV. Additionally, first experimental data on the associated non-prompt J/ψ meson and Z boson production have been presented by the ATLAS Collaboration at √ s = 8 TeV [6]. These measurements provide a new opportunity for studying b-flavor production in pp collisions at the LHC. A special point of interest because it involves both weak and strong interactions and therefore serves as a complex test of electroweak theory, perturbative Quantum Chromodynamics (pQCD) and our knowledge of parton distributions in a proton. The fixed-order next-to-leading-log (FONLL) approach [7] and the general-mass variable-flavor-number scheme (GM-VFNS) [8] are able to describe the inclusive production of non-prompt J/ψ and ψ(2S) mesons over the whole charmonia transverse momentum range. The role of the next-to-leading (NLO) pQCD corrections to the associated non-prompt J/ψ + Z production and contributions from the double parton scattering (DPS) mechanism were investigated very recently [9] using the MadGraph tool [10]. The main goal of our present study is to obtain first predictions for non-prompt J/ψ and ψ(2S) production, both inclusively and in association with a Z boson, employing the k T -factorization approach [11,12] and to give a systematic self-consistent analysis of recent CMS, ATLAS and LHCb data [1][2][3][4][5][6]. The k T -factorization approach has certain technical advantages in the ease of including higher-order radiative corrections that can be taken into account in the form of transverse momentum dependent (TMD) parton distributions 1 .
The calculations involve two main ingredients: the cross sections of b-hadron production and the partial widths of their subsequent decays into the different charmonia states. In our previous papers [14,15], we paid attention to the first of them, describing the inclusive production of b-hadrons and production of Z bosons in association with one or two b-jets (or rather b-hadrons) within the k T -factorization formalism. Our consideration below continues this line and extends to the second ingredient. In addition to the previously tested components of the theory, here we further inspect the evolution details of b-quarks fragmenting and decaying into the final state charmonia.
The outline of our paper is following. In Section 2 we briefly describe our theoretical input and calculation steps. In Section 3 we present a numerical results and a discussion. Section 4 contains our conclusions.

The model
In the present note we strictly follow the approach described earlier [14,15]. Here we only briefly recall the main points of the theoretical scheme.
First, to calculate the cross sections of inclusive b-hadron production and associated Z + b production in pp collisions the k T -factorization approach was applied, mainly based on the off-shell gluon-gluon fusion subprocesses: where the Z boson further decays into a lepton pair. The corresponding gauge-invariant off-shell (dependent on the transverse momenta of initial gluons) production amplitudes were calculated earlier (see [14,15] and references therein) and implemented into the Monte-Carlo event generator cascade [16]. Then, b-flavor production cross sections were calculated as a convolution of the off-shell partonic cross sections and TMD gluon distributions in a proton.
In addition to the off-shell gluon-gluon fusion, several subprocesses involving quarks in the initial state are taken into account. These are the flavor excitation, quark-antiquark annihilation and quark-gluon scattering subprocesses: with subsequent decays of Z bosons into a a lepton pair. These processes become important at large transverse momenta p T (or, respectively, at large parton longitudinal momentum fraction x, which is needed to produce high p T events) where the quarks are less suppressed or can even dominate over the gluon density. We find it reasonable to rely upon the collinear Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) factorization scheme [17], which provides better theoretical grounds in the large-x region. We consider a combination of two techniques where each of them being used at the kinematic conditions where it is best suitable (off-shell gluon-gluon fusion subprocesses at small x and quark-induced subprocesses at large x values). An essential point of our consideration is using a numerical solution of the Ciafaloni-Catani-Fiorani-Marchesini (CCFM) evolution equation [18] to derive the TMD gluon density in a proton. The CCFM equation provides a suitable tool since it smoothly interpolates between the small-x Balitsky-Fadin-Kuraev-Lipatov (BFKL) [19] gluon dynamics and high-x DGLAP dynamics. Following [15], we adopt the latest JH'2013 parametrization [20], taking JH'2013 set 2 as the default choice. This TMD gluon density was fitted to high-precision DIS data on the proton structure functions F 2 (x, Q 2 ) and F c 2 (x, Q 2 ). For the conventional quark and gluon densities we use the MSTW'2008 (LO) set [21]. The fragmentation of the produced b quarks into b-hadrons is described with the Peterson fragmentation function with b = 0.0126 [22]. Using these parameters, we successfully reproduce the latest CMS experimental data on inclusive B + meson production 2 at √ s = 13 TeV [24].
The approach above provides the necessary starting point for the theoretical description of J/ψ or ψ(2S) production from b-hadrons. At the next step, the obtained b-hadron cross sections have to be convoluted with b → J/ψ and b → ψ(2S) fragmentation functions. These fragmentation functions are the longitudinal momentum distributions of the J/ψ and ψ(2S) mesons from b-hadron decay, appropriately boosted along the b-hadron flight direction. The latter have been calculated [25] in the framework of nonrelativistic QCD (NRQCD) factorization [26,27] using the approach [28]. In this approach, the decays of b-hadrons into the final charmonia are represented by a sum of products. Each of them consists of a perturbative coefficient for the production of a cc pair in a specific angular-momentum and color state and nonperturbative NRQCD matrix element (NME). It parametrizes the subsequent transition of this intermediate cc state into the final physical charmonium with the help of nonperturbative soft gluon radiation and it has to be extracted from the experimental data. The resulting predictions [25] for the b → J/ψ+X and b → ψ(2S)+X three-momentum distributions reasonably agree with the measurements performed by the CLEO [29] and BABAR [30] Collaborations. Then, the developed formalism [25] is implemented into our calculations without any changes 3 . To be precise, we employ the asymptotic expression [25] for the b-hadron decay distribution differential in the longitudinal momentum fraction z carried by the produced charmonium state, obtained in the limit |p b | m b , where p b and m b are the momentum and the mass of decaying b-hadron. This approximation is valid within 11% and 5% accuracy for |p b | = 10 GeV and 20 GeV, respectively, that is suitable for our phenomenological study.
We confront our predictions with the available data in Figs. 1 -4. The solid histograms represent our central predictions by fixing both the renormalization µ R and factorization µ F scales at their default values (see [14,15] for the detailed description of our input), while the shaded regions correspond to the scale uncertainties of our predictions. To estimate the latter we used the JH'2013 set 2+ and JH'2013 set 2− sets instead of default one JH'2013 set 2. These two sets represent a variation of the renormalization scale used in the off-shell production amplitude [20]. We achieve good agreement with the LHC data on the transverse momentum distributions for both J/ψ and ψ(2S) mesons within the experimental and theoretical uncertainties. The slight disagreement is only  [3][4][5]. The cross section reported by the CMS Collaboration has been multiplied by the J/ψ → µ + µ − branching fraction.
observed at low transverse momenta, p T ≤ 10 GeV. It could be attributed to the finite |p b | terms in the longitudinal momentum distributions of the charmonia (important at |p b | ∼ m b ), which are not taken into account in our calculations. The rapidity distributions of J/ψ mesons measured by the LHCb Collaboration at different energies are well reproduced. Our results are consistently close to the FONLL [7] predictions presented in [1][2][3][4][5], which is explained by the fact that the main part of collinear QCD higher-order corrections (namely, NLO + NNLO + N 3 LO + ... contributions which correspond to the log 1/x enhanced terms in the perturbative series) are effectively taken into account as a part of the CCFM gluon evolution. Here we demonstrate again the main advantage of the k T -factorization approach, which gives us the possibility to estimate the size of these higher-order corrections and reproduce in a straighforward manner the main features of rather cumbersome higher-order pQCD calculations.
The estimated non-prompt J/ψ production fiducial cross sections at √ s = 7, 8 and 13 TeV are listed in Table 1 in comparison with the available LHC data. One can see that our predictions agree well with these data within the theoretical and experimental uncertainties.
Now we turn to associated non-prompt J/ψ + Z production. First experimental data were obtained recently by the ATLAS Collaboration at √ s = 8 TeV [6]. The leptons originating from the Z boson decay are required to have pseudorapidities |η l | < 2.5, transverse momenta p l T > 25 GeV (leading lepton) and p l T > 15 GeV (sub-leading lepton) and invariant mass of the lepton pair M ll should lie within the interval |M ll − m Z | < 10 GeV. The produced J/ψ meson is required to have the transverse momentum 8.5 < p T < 100 GeV and rapidity |y| < 2.1. Fiducial cross sections were measured as ratios to the inclusive Z boson production rate in the same fiducial volume. We implemented this experimental setup into our numerical program. The results of our calculations are shown in Fig. 5. One can see that, similar to inclusive J/ψ production case, our predictions for associated J/ψ + Z production reasonably agree with the ATLAS data [6]. We find that both the transverse momentum and rapidity distributions of the produced J/ψ mesons are described within the theoretical and experimental uncertainties. However, we note that the overall description could be even improved if additional contributions from the DPS mechanism are taken into account. To estimate the latter we apply a simple factorization formula (for details see the reviews [32][33][34] and references therein): where σ eff is a normalization constant which incorporates all "DPS unknowns" into a single phenomenological parameter. A numerical value of σ eff 15 mb was obtained from fits to pp and pp data (see, for example, [35]). We find that the DPS mechanism gives ∼ 10% contribution to the production cross section in the considered kinematical region and populates mainly at low transverse momenta (see Fig. 5). Being added to the predictions of usual single parton scattering (SPS) mechanism, the DPS contributions slightly improve the description of the J/ψ rapidity distribution in a whole y range. Some reasonable variations in σ eff 15 ± 5 mb would affect DPS predictions, though without changing our basic conclusion 4 . The presented DPS estimation coincides with the earlier one performed by the ATLAS Collaboration [6].

Conclusions
We investigated the production of J/ψ and ψ(2S) mesons, originating from the decays of b-flavored hadrons, both inclusively and in association with Z bosons, at the LHC conditions using the k T -factorization approach. Our consideration was mainly based on the dominant off-shell gluon-gluon fusion subprocesses where the transverse momenta of initial gluons are taken into account. In the case of associated non-prompt J/ψ + Z production, a number of subleading quark-induced subprocesses have been considered in the conventional collinear scheme. To describe the inclusive b-hadron decays into the different charmonium states, we applied fragmentation functions, based on the NRQCD factorization. The feed-down contributions from χ cJ (with J = 0, 1, 2) and ψ(2S) decays to J/ψ meson production were taken into account.
Using the TMD gluon densities derived from the CCFM gluon evolution equation, we achieved a good agreement between our predictions and latest CMS, ATLAS and LHCb data collected at √ s = 7, 8 and 13 TeV. We provided first theoretical expectations for the associated non-prompt J/ψ + Z production and estimated the contribution of the double parton scattering mechanism to the production cross sections. The latter is found to be of order of ∼ 10%.      The differential cross sections of associated non-prompt J/ψ + Z production at √ s = 8 TeV as functions of J/ψ transverse momentum and rapidity. Notation of curves is the same as in Fig. 1. The estimated DPS contributions are shown separately and not summed with SPS ones. The experimental data are from ATLAS [6].