Study of heavy meson production in p-Pb collisions at $\sqrt{S}$=5.02 TeV in the general-mass variable-flavour-number scheme

We study inclusive charm and bottom production, for both D and B mesons, in p-Pb collisions at the LHC. Numerical results for p_T-differential production cross sections are obtained at next-to-leading-order in the general-mass variable-flavor-number scheme. We compare our results with recent data from ALICE, LHCb and CMS at a center-of-mass energy of 5 TeV and find good agreement. A comparison with p-p cross sections does not reveal the presence of nuclear initial-state interaction effects that could be expected to become visible as deviations of the ratio of p-Pb and p-p cross sections from one.


Introduction
The study of heavy-quark (charm or bottom) production in p-p collisions at LHC energies is a useful test of perturbative Quantum Chromodynamics (QCD) since the heavy quark mass provides a hard scale that allows calculations within perturbation theory. The QCD calculations are based on the factorization approach, in which cross sections are calculated as a convolution of three terms: the parton distribution functions (PDF) of the incoming protons, the partonic hard scattering cross sections computed as a perturbative series in the strong interaction coupling constant, and the fragmentation functions (FF), describing the relative production yield and momentum distribution for a given heavy hadron (D or B meson) in a parton. Corresponding recent calculations at the perturbative level at next-to-leading order (NLO) with next-to-leading-log resummation (FONLL) [1,2] or in the framework of the general-mass-variable-flavour-number scheme (GM-VFNS) [3,4] have provided good descriptions for bottom meson production inp-p collisions at √ S = 1.96 TeV at the FNAL Tevatron Collider [5][6][7] and in p-p collisions at √ S = 7 TeV at the CERN Large Hadron Collider (LHC) by the CMS, ATLAS and the LHCb collaborations [8][9][10][11][12][13]. The production cross section of charmed hadrons (D mesons) at the Tevatron [14] or of the ATLAS collaboration at the LHC [15] is also reasonably well described within theoretical and experimental uncertainties [16,17].
The GM-VFNS is essentially the conventional NLO parton-model approach, supplemented with finite-mass effects, intended to improve the description at small transverse momentum p T . The original GM-VFNS prescription [3,4,17] is, however, not suitable for calculations of the cross section dσ/dp T for heavy-quark hadron production at very small transverse momentum p T . This is due to the specific choice of scale parameters for initial-state (µ I ) and final-state (µ F ) factorization. The original prescription was to set µ I = µ F = m 2 Q + p 2 T , where m Q is the mass of the heavy quark, charm or bottom. At p T = 0, the scale parameters approach µ I = µ F = m Q , and at this point the heavy quark PDFs are put to zero by construction in almost all available PDF parametrizations. Therefore the transition to the fixed-flavour-number-scheme (FFNS), which is the appropriate scheme for calculating dσ/dp T at rather small p T , is not reached for p T > 0, since the heavy quark PDF in the proton decouples at p T = 0, and not for finite p T > 0.
A smooth transition to the FFNS at finite p T can be achieved by exploiting the freedom to choose the factorization scale. In Refs. [18,19] we have studied the prescription to fix the initial-state factorization scale at µ = 0.5 m 2 Q + p 2 T instead of µ = m 2 Q + p 2 T . For simplicity we have chosen the scales for initial and final state factorization equal to each other, µ I = µ F . With this scale choice we could achieve a reasonably good description of the data for B meson production down to p T = 0 for the CDF data [6] inp-p collisions at the Tevatron and of the LHCb data [13] for p-p collisions at the LHC in the forward rapidity region at √ S = 7 TeV. A comparison of data for all D meson states D 0 , D + , D * + and D + s measured by the LHCb collaboration at √ S = 5, 7 and 13 TeV with predictions from the GM-VFNS scheme with the original scale choice for p T > 3 GeV can be found in [22][23][24].
The LHC Collaborations have also measured cross sections for heavy-quark production in p-Pb and Pb-Pb collisions. The ALICE collaboration, e.g., have performed detailed studies of the p T -differential and rapidity-differential cross sections dσ/dp T and dσ/dy for D-meson production in p-Pb collisions at √ S = 5.02 TeV [25,26], also for small p T , as well as in Pb-Pb collisions at √ S = 2.76 TeV [27]. Collisions with two heavy nuclei are of particular interest for studies of the Quark-Gluon Plasma (QGP), a high-density colourdeconfined medium. On the other hand, data from p-Pb collisions can be used to determine the nuclear modification factor R pPb , i.e., the ratio of p-Pb cross sections relative to the corresponding p-p cross sections scaled by the mass number of the Pb nucleus (A = 208). Data are in particular interesting at small p T where one expects the largest deviation from R pPb = 1. The value of R pPb is of interest for several reasons. First large deviations from one, in particular for larger p T , would indicate the presence of initial-state interaction effects which are needed to obtain a reliable interpretation of corresponding Pb-Pb collision data. Second, the value of R pPb is of interest by itself and could help to obtain information on the nuclear PDFs, which are modified compared to the proton PDFs in bound nucleons depending on the parton fractional momentum x and the atomic mass number A.
Ideally, measurements of the cross sections to determine the nuclear modification factor R pPb should be done at the same center-of-mass energy √ S. Unfortunately, this is not the case; data for p-p and p-Pb collisions at the same √ S are not available. Instead, the reference p-p cross section at √ S = 5.02 TeV was obtained from data at √ S = 7 TeV [28] by scaling the energy based on predictions from perturbative QCD. The scaling factor was determined for each D-meson species separately from the FONLL calculations [29]. In case of B meson production in p-Pb collisions at √ S = 5.02 TeV, measured by the CMS collaboration [34], the reference cross section dσ/dp T for p-p collisions was directly taken from the FONNL calculations at √ S = 5.02 TeV [29] without any extrapolation from their data at larger √ S.
Due to the interest in the nuclear modification factor R pPb for heavy quark hadron production, in particular as we expect to obtain important information about initial-state interaction effects in Pb-Pb collisions, it is desirable to study R pPb also within other factorization schemes. This is the purpose of the present work in which we provide results from calculations of p-p cross sections dσ/dp T for D and B meson production at √ S = 5.02 TeV in the framework of the GM-VFNS. We compare our results with data for the production of various D meson species at √ S = 5.02 TeV measured by the ALICE [25,26] and LHCb collaborations [30] and for B meson production at √ S = 5.02 TeV measured by the CMS collaboration [34]. Using our results for the p T -differential cross sections, we also study the nuclear modification factor R pPb .
The outline of our work is as follows. In the next section, Sect. 2, we give the details of the calculations for D mesons with the kinematic constraints of the ALICE and LHCb experiments. Section 3 contains our results for B meson production at √ S = 5.02 TeV and a comparison with the CMS data. Section 4 is reserved for a discussion of the results.

D meson production in p-p and p-Pb collisions
The theoretical background and explicit analytic results of the GM-VFNS approach were previously presented in detail, see Refs. [3,4] and the references cited therein. Here we only describe the input needed for the present numerical analysis.
Throughout this paper, we use the PDF set CTEQ14 [35] as implemented in the program library LHAPDF [36]. The fragmentation functions determined in Ref. [37] for D 0 , D + and D * + mesons and in Ref. [38] for the D + s meson were used. These FFs always refer to the average of charge-conjugated states. The data from ALICE and CMS are understood as averaged cross sections as well, (σ(D) + σ(D))/2 and (σ(B) + σ(B))/2, while the LHCb collaboration decided to present their data as the sum of charge-conjugated states.
Originally, the default value for the scale parameters for renormalization and factorization were set by the transverse mass m T = m 2 Q + p 2 T . By convention, variations around a default value by factors of two up and down were considered to obtain an estimate of unknown higher-order perturbative contributions and, thereby, assign a theoretical uncertainty to numerical results. We introduce the dimensionless parameters ξ i (i = R, I, F ) and set µ i = ξ i m T . Independent variations of the ξ i between 1/2 and 2 are restricted by keeping any ratio of the ξ i 's smaller than 2. We shall denote this choice of scales as the original prescription.
As already mentioned, this original scale choice does not provide a smooth transition to the FFNS at small p T . To achieve this we change the factorization scales to µ I = µ F = ξ 0 4m 2 Q + p 2 T with ξ 0 = 0.49. A similar choice with ξ 0 = 0.5 was used in a recent study of charm meson production [20]. In Ref. [21], using also ξ 0 = 0.49, good agreement was found with p-p collision data from the LHCb experiments [22][23][24] for p T values down to p T = 0. The choice of 4m 2 c + p 2 T in place of the transverse mass m T = m 2 c + p 2 T is motivated by the fact that the kinematic threshold for heavy-quark production is at 2m c in the FFNS approach. With the additional factor ξ 0 = 0.49 we can ensure that µ = m Q is reached already slightly above p T = 0. For m Q = m c = 1.3 GeV one has µ = m Q at p T = 0.528 GeV. We choose this value of m c to be consistent with the value used in the PDF set CTEQ14 from Ref. [35]; otherwise a smooth decoupling of the charm content of the proton PDF is not achieved. In our earlier calculations for larger values of p T [17] we had adopted m c = 1.5 GeV instead. We determine error bands for theoretical uncertainties from variations of the renormalization scale only, i.e., by varying ξ R between 1/2 and 2. We have to leave the factorization scales unchanged since otherwise the proper transition to the FFNS would be lost. This setting of scales will be called the modified scale in the following.
Before we apply this scale choice for a comparison with the ALICE data in p-Pb collisions at small p T [26], we have a look at the reference p-p cross section. The most precise data for the p T -differential cross section of prompt D 0 meson production at √ S = 7 TeV was obtained by a combination of measurements without decay-vertex reconstructed in the low-p T range, 0 < p T < 2 GeV, and an analysis using information from decay-vertex reconstruction at larger p T , 2 < p T < 16 GeV. In all cases, the rapidity is restricted to the range |y| < 0.5 and contributions from the b → D 0 transition have been subtracted. Data and results from the GM-VFNS are shown in Fig. 1 (left panel). The agreement with the default scale is very good in the large p T range, p T > 6 GeV, and for all p T values the data points lie inside the theoretical range obtained from the scale variation of µ R .
The ratio of data for dσ/dp T normalized to our prediction in the GM-VFNS with the modified scale choice is shown in the right panel of Fig. 1 (full-line histogram). For the larger p T ≥ 6 GeV the ratio is equal to one within the experimental accuracy. This is consistent with the prediction of the original scale choice, for which the same ratio was shown in Ref. [26] for p T ≥ 3 GeV. For the smaller p T , 1 < p T < 6 GeV, the ratio in Fig. 1 (right panel) increases to approximately 1.5. This is very similar to results based on the FONLL approach [29] and on the LO k T factorization calculation [39], which was also shown in [26]. The dashed-line histograms in the right panel of Fig. 1 show the ratios  Figure 1: Left panel: Differential production cross section dσ/dp T of prompt D 0 mesons in p-p collisions at √ S = 7 TeV with |y| < 0.5 in the p T interval 1 < p T < 16 GeV compared to ALICE data [26,28]. The data point for the bin 1 < p T < 2 GeV is from the analysis [26], while the data points for 2 < p T < 16 GeV are taken from [28]. The theoretical cross sections are calculated in the GM-VFNS with default scales µ R = 4m 2 c + p 2 T and µ I = µ F = 0.49 4m 2 c + p 2 T . The upper and lower dashed histograms are calculated with µ R changed by factors 1/2 and 2. Right panel: Ratios (see text) of the ALICE data over theory predictions. of the same data, but normalized to the GM-VFNS prediction with µ R varied by factors 1/2 and 2. In order to keep the plot readable, we do not display the error bars for the experimental uncertainties in this case. The band between the dashed histograms thus represents the scale uncertainty of the ratio dσ Data /dσ GM−VFNS . We observe that inside the scale variation this ratio is compatible with one. Now we continue with a comparison of theory predictions and ALICE data for p-Pb collisions. Theoretical predictions are obtained from the p-p cross section by multiplication with the mass number A = 208, Adσ/dp T . Data are available at √ S = 5.02 TeV in the rapidity region |y| < 0.5. Our results in the GM-VFNS with the modified scale choice are shown in Figs. 2, 3, 4, and 5 (left panels) for D 0 , D + , D * + and D + s production, in each case together with the data from [26] as a function of p T for bins in the range 1 < p T < 24 GeV. Except for two points at the largest p T (see Figs. 3 and 4) the error bars of the data points overlap with the uncertainty range due to scale variations. As for p-p collisions, the ALICE data shown in Fig. 2 are obtained for prompt D 0 production in the interval 0 < p T < 2 GeV (only data for p T > 1 GeV are shown) without decay-vertex reconstruc- (d /dp T ) pPb /A(d /dp T ) pp,GM-VFNS p Pb D 0 X GM-VFNS S = 5.02 TeV -0.5 y lab 0.5 ALICE Data     tion [26] and for p T > 2 GeV with decay-vertex reconstruction [25]. The data for the other three D-meson species D + , D * + and D + s are taken from Ref. [25]. Corresponding ratios for ALICE data normalized to our theoretical results for Adσ/dp T are presented in the right panels of Figs. 2, 3, 4, and the lower panels of 5. Again, we decide to present scale uncertainties by normalizing the data to varied theory predictions with µ R scaled up and down by factors of 1/2 and 2 (dashed-line histograms) and show the ratios R ± = dσ p−Pb,data /(Adσ p−p,GM−VFNS (µ ± )) where µ ± denotes the varied renormalization scale. The band enclosed by R ± should contain unity if there is a scale choice which leads to agreement between theory and experiment. This is indeed the case, except for D + and D * + production at the largest p T where the ratio falls slightly below one.
The shape of the p T dependence of the ratios looks rather similar for all cases, compare for example the case of D 0 production for p-Pb collisions at √ S = 5.02 TeV in Fig. 2 and for p-p collisions at √ S = 7 TeV in Fig. 1. The similarity between p-p and p-Pb collisions is even more clearly visible when we consider the ratios of the results shown in the right panels of Figs. 1 and 2. This is done in Fig. 6 where we show where µ i denotes the renormalization scale varied up and down by factos of 1/2 and 2 (d /dp T ) pPb /A(d /dp T ) pp,GM-VFNS normalized to ratio for p p p Pb D 0 X GM-VFNS S = 5.02 TeV -0.5 y lab 0.5 ALICE Data  around its central value µ 0 . The first factor in brackets is represented by the full histogram in Fig. 6. It is the ratio of p-Pb over p-p data, properly normalized to the same value of √ S using the GM-VFNS prediction. The error bars shown here represent the uncertainty of the p-Pb data only. The band of dashed-line histograms represents an estimate of the scale uncertainty, evaluated at √ s = 5.02 TeV (see the last factor in the definition of R i given above). Since R i ≡ 1 is contained inside this band we conclude that the data do not require corrections, for example due to initial-state interactions in the Pb nucleus.
For the other mesons, D + , D * + and D + s in Figs. 3, 4, and 5 the pattern of ratios looks also quite similar. For the larger p T bins the ratio is equal to one within errors, and for the smaller p T bins the ratio is close to 1.5. We remark that the nuclear modification factor R pPb is consistent with one for all four D meson species if the theoretical uncertainty due to scale variations is taken into account.
We can compare our results with the nuclear modification factor presented in Ref. [26]. The ratios R pPb for D 0 , D + and D * + given there are much closer to one than our calculated ratios shown in Figs. 2, 3, and 4. Note that the p-p cross sections used in Ref. [26] to obtain the ratios R pPb have been deduced from the measured cross sections at √ S = 7 TeV by extrapolation to √ S = 5.02 TeV. It would be premature to interpret the observed small deviations of the nuclear modification factors from one as a sign of initial-state interaction effects as long as we see similar deviations for p-p collisions as shown in Fig. 1, right panel. It has been shown in Ref. [26] that theoretical expectations for deviations of R pPb from one for several models existing in the literature are rather small at large p T . Only towards small values of p T model predictions start to deviate from one by more than 10 percent or so. Future higher-precision data may allow to exclude some of the theoretical approaches, but right now experimental uncertainties are still too large to draw any firm conclusion.
Finally we compare predictions from the GM-VFNS approach with most recent data from the LHCb collaboration [30]. For p-p collisions, a rather good agreement between theory predictions and LHCb data for the differential cross section dσ/dp T in various rapidity bins in the forward direction was already observed in Ref. [21]. The recent measurements of p-Pb cross sections at LHCb [30] have provided us with more information about the dependence on the rapidity y cm in the nucleon-nucleon centre-of-mass system and allow us to study the forward and backward regions separately. We note that experimental uncertainties are much smaller than for the other measurements described before. In Fig.  7 we show two sets of plots, one for the forward region, 1.5 ≤ y cm ≤ 4.0 (upper plots) and one for the backward region −5.0 ≤ y cm ≤ −2.5 (lower plots). All data points agree with theory within the scale uncertainty band. In the right plots of Fig. 7 we show ratios of data for p-Pb collisions normalized to A times theory predictions for p-p scattering. The deviation of these ratios from one are not very large in the forward region, but significantly above one for backward rapidities. We expect that this observation can be explained by using appropriately chosen nuclear PDFs. At present, nuclear PDFs have very large errors [31][32][33] and a direct comparison with the available nPDFs is not very instructive. However, one can conclude that these precise LHCb data will help to narrow down possible nPDF parametrizations. We note that the forward-backward ratio discussed in the LHCb  The p T distribution for D 0 +D 0 production in p-Pb collisions compared with data from the LHCb collaboration taken at √ S = 5 TeV. The left plots show the differential cross sections dσ/dp T , the right plots show the ratios of data over theory. Full and dashed lines are defined as in the previous figures (see also text). Data and ratios in the upper part are for the forward region 1.5 ≤ y cm ≤ 4.0 and in the lower part for the backward region −5.0 ≤ y cm ≤ −2.5.
publication will be particularly interesting for a study of nuclear PFFs since it is not affected by large scale uncertainties.

B meson production in p-Pb collisions
Up to now, cross section data of dσ/dp T for B-meson production (B + , B 0 and B 0 s ) in p-Pb collisions at √ S = 5.02 TeV are available only for larger p T values above 10 GeV [34], in the range 10 < p T < 60 GeV. In Ref. [34] data have been compared with A times the FONLL prediction for p-p collisions [29]. At √ S = 7 TeV the LHCb collaboration has measured the p-p cross section dσ/dp T down to p T = 0 for B + + B − , B 0 +B 0 and B 0 s +B 0 s production in the forward region 2 ≤ y ≤ 4.5 [13,40]. These data have been compared with our GM-VFNS predictions using the modified scale 0.5 m 2 b + p 2 T . The comparison between the LHCb data and our predictions showed reasonably good agreement for all three B meson species [18]. In this reference we compared the GM-VFNS predictions also for B + -meson production measured by the ATLAS collaboration [12] where data extend into the very large p T -range, 9 < p T < 120 GeV, for various rapidity intervals in the range 0 < |y| < 2.25. In this comparison we found agreement between data and theory except for the lowest p T bin, 9-13 GeV, where the data are slightly overestimated.
In the following we show the results for Adσ/dp T at √ S = 5.02 TeV in the rapidity interval |y| < 2.4, again obtained from the p-p cross section dσ/dp T by multiplication with the mass number A. We have done these calculations for the original scale choice µ o = m 2 b + p 2 T ; for the modified choice we decided to choose µ m = 0.5 m 2 b + p 2 T = 0.5µ o in order to allow for a direct comparison with the previous work [18] 1 . m b is the bottom quark mass, m b = 4.5 GeV. The FF for b → B was taken from [5] for all three B meson states. Cross sections for the different B meson species differ only by their respective constant fragmentation fractions. Our results are compared to the CMS data for p-Pb collisions [34] and are shown for B + , B 0 and B 0 s production, respectively, in the left panels of Figs. 8, 9, and 10 for µ = µ o and in the right panels of these figures for µ = µ m . As to be expected the results for the original scale choice µ o lie slightly higher than for the modified scale choice µ m , but the difference is decreasing towards larger p T . For all cases data and theory agrees within theoretical and experimental errors.
The comparison between the experimental cross section dσ/dp T for p-Pb scattering and the theoretical cross sections Adσ/dp T becomes more clear when presented in terms of the nuclear modification factors R pPb = (dσ/dp T ) pPb /A(dσ/dp T ) pp . We show these ratios for all three B meson species and for both scale choices, µ o and µ m , in Figs. 11, 12, and 13 (left and right panels). We notice that with the modified scale choice, the ratio R pPb agrees with one within experimental errors, even without taking into account the theory uncertainty 1 The value of µ I,F at p T = 0 is not very relevant here since we will compare with data at large p T . With µ I,F = 0.5 4m 2 b + p 2 T the cross section would increase by only 12 % in the first p T -bin (10 GeV ≤ p T ≤ 15 GeV) and by less than 2 % at higher p T .  [11][12][13]. For the modified scale choice our results agree also rather well with those presented in [34] where the p-p cross section used to obtain R pPb was calculated in the FONLL approach [29].
Our results for the nuclear modification factor R pPb compared with CMS data differ somewhat for the two scale choices µ = µ o and µ = µ m (compare left and right panels of Figs. 11,12,13). For µ = µ m the ratios R pPb are equal to one for all bins within the precision of the data. For the original scale choice µ = µ o deviations from one seem to occur already within present errors in some of the low-p T bins (see left panels of Figs. 11, 12, and 13). However, the observed deviations would become significant only if the experimental errors could be reduced, by at least a factor of two. It seems obvious to us that also theory uncertainties will have to be reduced before a conclusive interpretation of the data will be possible. This will require the calculation of higher-oder corrections which are expected to reduce the uncertainties due to the choice of renormalization and factorization scales.

Conclusions
We have studied D and B meson production in p-Pb collisions and made, for the first time, a comparison with predictions obtained at NLO in the GM-VFNS. Our main results are shown in the right panels of Figs. 2 -5, 7 for D-meson production and in Figs. 11 -13 for B-meson production. The comparison with data confirms our previous findings that a suitable choice of the factorization scale parameters can be found which brings the experimental data obtained by the LHC collaborations ALICE, CMS and LHCb into good agreement with predictions obtained in the general-mass variable-flavour-number scheme.
The ratio of data for p-Pb collisions over theory predictions for A times p-p cross sections is an important observable which could provide information about the nuclear modification of parton distribution functions, for example due to initial-state interaction effects. We found that for charmed meson production, the ratios of data over theory predictions at p T > 6 GeV are compatible with one within uncertainties and deviations are not larger than 40%. At small transverse momenta, p T < 6 GeV, the ratios for data from ALICE at mid-rapidity increase to values of about 1.5 and larger. The data from the LHCb collaboration for D-production in the forward region, however, do not show such a strong enhancement of the nuclear modification ratio. Interestingsly, the ratio of p-p data over theory show deviations from one of the same size and with a similar p T -dependence. It will be interesting to include forthcoming more precise data in our analysis, as for example from Ref. [41].
Experimental uncertainties are often still large, but data are steadily improving. In particular the most recent data from LHCb are promising and one can expect that updated fits of nuclear PDFs with smaller uncertainties than the existing parametrizations will be possible. At present, however, scale uncertainties are still very large and it is therefore doubtful whether the observed deviations can be interpreted as due to nuclear modification effects. Higher precision of the measurements as well as of theory predictions is needed in order to draw firm conclusions.