Inclusive J/$\psi$ production at forward and backward rapidity in p-Pb collisions at $\sqrt{s_{\rm NN}}$ = 8.16 TeV

Inclusive J/$\psi$ production is studied in p-Pb interactions at a centre-of-mass energy per nucleon-nucleon collision $\sqrt{s_{\rm NN}}=8.16$ TeV, using the ALICE detector at the CERN LHC. The J/$\psi$ meson is reconstructed, via its decay to a muon pair, in the centre-of-mass rapidity intervals $2.03<y_{\rm {cms}}<3.53$ and $-4.46<y_{\rm {cms}}<-2.96$, where positive and negative $y_{\rm {cms}}$ refer to the p-going and Pb-going direction, respectively. The transverse momentum coverage is $p_{\rm T}<20$ GeV/$c$. In this paper, $y_{\rm cms}$- and $p_{\rm T}$-differential cross sections for inclusive J/$\psi$ production are presented, and the corresponding nuclear modification factors $R_{\rm pPb}$ are shown. Forward results show a suppression of the J/$\psi$ yield with respect to pp collisions, concentrated in the region $p_{\rm T}\lesssim 5$ GeV/$c$. At backward rapidity no significant suppression is observed. The results are compared to previous measurements by ALICE in p-Pb collisions at $\sqrt{s_{\rm NN}}=5.02$ TeV and to theoretical calculations. Finally, the ratios $R_{\rm FB}$ between forward- and backward-$y_{\rm {cms}}$ $R_{\rm pPb}$ values are shown and discussed.


Introduction
Quarkonium production in nuclear collisions is sensitive to the temperature of the produced medium. In particular, the various quarkonium states are expected to melt in a Quark-Gluon Plasma (QGP), due to screening of the colour interaction in a deconfined state [1]. In addition, the abundant charm-quark production in the multi-TeV collisionenergy range can also lead to a (re)generation of charmonia during the QGP evolution and at the phase boundary [2,3]. A detailed investigation of these processes was carried out by ALICE, which has measured inclusive J/ψ production in Pb-Pb collisions down to zero transverse momentum (p T ). These results were reported at centre-of-mass energies per nucleon pair √ s NN = 2.76 and 5.02 TeV, at forward centre-of-mass rapidity y cms [4][5][6][7] for both energies, and at central y cms for √ s NN = 2.76 TeV [8]. The nuclear modification factor R AA was evaluated, corresponding to the ratio between the Pb-Pb and the pp production cross sections, normalised to the number of nucleon-nucleon collisions. A suppression of the J/ψ was observed, as indicated by values of R AA smaller than unity. However, the suppression was found to be systematically smaller with respect to results obtained at RHIC energies [9,10]. In addition, the suppression effects were less strong at low p T . These observations, together with the comparison to theoretical model calculations [11][12][13][14] and the measurement of a non-zero elliptic flow for the J/ψ [15], imply that a fraction of the J/ψ yield is produced via recombination of charm quarks, and that recombination is more prevalent at low p T , where the bulk of charm-quark production occurs. In addition to effects connected with the hot medium, cold nuclear matter (CNM) effects are expected to influence the charmonium yield in nuclear collisions. One of the most important is nuclear shadowing, i.e., the modification of the quark and gluon structure functions for nucleons inside nuclei (see e.g., refs. [16][17][18]). This effect modifies the probability for a quark or a gluon to carry a given fraction x of the momentum of the nucleon. It affects the elementary production cross section for the creation of the cc pair that will -1 -JHEP07(2018)160 eventually form a charmonium state. Modifications of the initial state of the nucleus are also addressed by calculations incorporating parton saturation, a coherent effect involving low-x quarks and gluons, described by the Colour Glass Condensate (CGC) effective theory [19]. In addition to these mechanisms, a coherent energy-loss effect involving partons in the initial and final state can also lead to a modification of the parton kinematics and consequently to a change in the quarkonium yields with respect to elementary nucleonnucleon collisions [20]. Finally, once produced, the charmonium state could be dissociated via inelastic interactions with the surrounding nucleons [21]. This process, which plays a dominant role among CNM effects at low collision energy [22,23], should become negligible at the LHC, where the crossing time of the two nuclei is much shorter than the formation time of the resonance [24][25][26].
The CNM effects introduced above are present in nucleus-nucleus collisions, but can be more directly investigated by studying proton-nucleus collisions, where the contribution of hot-matter effects are thought to be negligible. Previous results from p-Pb collisions at √ s NN = 5.02 TeV from ALICE [27][28][29], LHCb [30] and CMS [31] have shown a significant suppression of the J/ψ yield at forward rapidity (p-going direction) and low to intermediate p T ( 5 GeV/c). No significant effects, or at most a slight enhancement, were seen at high p T and at backward y cms (Pb-going direction). The results were compared to theoretical calculations that include various combinations of all the effects mentioned in the previous paragraph, except charmonium dissociation in cold nuclear matter [32][33][34][35][36][37]. A good agreement with the models was found, indicating on the one hand that mechanisms like shadowing, CGC-related effects and coherent energy loss can account for the observed nuclear effects, and on the other hand that final state break-up processes in nuclear matter have a negligible influence. It should be noted that the model of ref. [34] includes the effects of the interaction of charmonia with a dense hadronic medium possibly created in p-Pb collisions. However, such a medium may be expected to dissociate the weakly bound ψ(2S) state [26], but should have little or no effect on the strongly bound J/ψ meson. In 2016, p-Pb collisions at √ s NN = 8.16 TeV were delivered by the LHC. The interest in J/ψ studies at this energy is threefold: first, a significantly larger integrated luminosity with respect to studies performed at √ s NN = 5.02 TeV [27][28][29] has become available in ALICE, allowing a more detailed comparison to model calculations and an extended p T reach. Second, by varying the collision energy, it is possible to extend the investigations of shadowing and other CNM effects to a partly different x range. Finally, studies of various physics observables in p-Pb and high-muliplicity pp collisions at the LHC have shown effects such as long-range two-particle correlations [38][39][40][41][42][43] and an enhancement of strange and multi-strange hadron production [44], already seen in Pb-Pb collisions. These effects are usually connected with the formation of an extended system of strongly interacting particles. Concerning the specific case of charmonium production, in addition to the observations discussed above, long-range correlation structures in J/ψ production were recently observed in p-Pb collisions at √ s NN = 8.16 TeV [45]. Furthermore, for the weakly bound ψ(2S) a suppression signal, on top of the CNM effects discussed in the previous paragraphs, was seen in p-Pb and related to the resonance break-up in the medium created in such collisions [26,46]. As mentioned above, no extra suppression needs to be -2 -

JHEP07(2018)160
introduced for the strongly bound J/ψ in order to reproduce the experimental observations at √ s NN = 5.02 TeV. However, higher energy p-Pb collisions may create a more extended and longer-lived medium, which might lead to a suppression effect also on the J/ψ. In this paper, we report ALICE results on cross sections and nuclear modification factors for inclusive J/ψ production in p-Pb collisions at √ s NN = 8.16 TeV, in the rapidity regions 2.03 < y cms < 3.53 and −4.46 < y cms < −2.96, and for p T < 20 GeV/c. In section 2, the experimental apparatus, the data sample and the event selection criteria are presented. section 3 contains a description of the analysis procedure, including a discussion of the evaluation of the systematic uncertainties. The results and their comparison to theoretical models, to recent LHCb results [47] and to √ s NN = 5.02 TeV data are shown in section 4, while conclusions are drawn in section 5.
2 Experimental apparatus, data sample and event selection The ALICE detector design and performance are extensively described in [48,49]. The analysis presented here is based on the detection of muons in the ALICE forward muon spectrometer [50], which includes five tracking stations (Cathode Pad Chamber detectors), followed by two triggering stations (Resistive Plate Chamber detectors). An absorber, 10 interaction-length (λ I ) thick and made of carbon, concrete and steel, positioned in front of the tracking system, filters out most hadrons produced in the collision. A second (7.2 λ I thick) iron absorber, positioned between the tracking and the triggering system, absorbs secondary hadrons escaping the first absorber and low-momentum muons. Finally, a 3 T·m dipole magnet, positioned in the region of the third tracking station, provides the track bending for momentum evaluation. Particles are detected in the pseudo-rapidity range −4 < η < −2.5 in the laboratory system and muon triggering is performed with a programmable transverse momentum threshold, set to p µ,T = 0.5 GeV/c for the data sample analysed in this paper. The trigger threshold is not sharp, and the single muon trigger efficiency reaches its plateau value (∼ 96%) at p µ,T ∼ 1.5 GeV/c. In addition to the muon spectrometer, four other sets of detectors play an important role for this analysis. The Silicon Pixel Detector (SPD) [51], with its two layers covering the pseudo-rapidity intervals |η| < 2 and |η| < 1.4, is part of the ALICE central barrel and is used to reconstruct the primary vertex. A coincidence of a signal in the two V0 scintillator detectors [52], covering 2.8 < η < 5.1 and −3.7 < η < −1.7, provides a minimum-bias (MB) trigger. The luminosity determination is obtained from the V0 information and, independently, using the T0 Cherenkov detectors [53], which cover 4.6 < η < 4.9 and −3.3 < η < 3.0. Finally, the timing information from the V0 and the Zero Degree Calorimeters (ZDC) [54] is used to remove beam-induced background.
The trigger condition used in the analysis is a µµ − MB trigger formed by the coincidence of the MB trigger and an unlike-sign dimuon trigger. By taking data in two configurations of the beams corresponding to either protons or Pb ions going towards the muon spectrometer, it was possible to cover the dimuon rapidity ranges 2.03 < y cms < 3.53 and −4.46 < y cms < −2.96, respectively. The two configurations are also referred to as p-Pb and Pb-p in the following.

JHEP07(2018)160
The data samples used in this analysis correspond to an integrated luminosity L pPb int = 8.4 ± 0.2 nb −1 for p-Pb, and L Pbp int = 12.8 ± 0.3 nb −1 for Pb-p collisions [55]. These values are larger by about a factor 2 with respect to √ s NN = 5.02 p-Pb collision data [27].
The selection criteria used by ALICE in previous J/ψ analyses [27,28] have been applied. Namely, both muons belonging to the pair must have −4 < η µ < −2.5, to reject tracks at the edges of the acceptance. In addition, each muon must have 17.6 < R abs < 89.5 cm, where R abs is the radial transverse position of the muon tracks at the end of the absorber, to remove tracks crossing its thicker region, where energy loss and multiple scattering effects are more important. Finally, each track reconstructed in the tracking chambers of the muon spectrometer has to match a trigger track reconstructed in the trigger system.

Data analysis
The analysis procedure is the same for the two data sets discussed in this paper, and very similar to the one reported in refs. [27,28] for the √ s NN = 5.02 TeV p-Pb sample. The inclusive J/ψ production cross section was obtained from is the corresponding product of acceptance times reconstruction efficiency, B.R.(J/ψ → µ + µ − ) = 5.961 ± 0.033% is the branching ratio for the decay to a muon pair [56] and L pPb int is the integrated luminosity for the data sample under study.
The quantities N J/ψ (∆y cms , ∆p T ) were obtained through fits to the invariant mass spectra of the opposite-sign muon pairs. The fitting functions are the sum of two resonance contributions (J/ψ and ψ(2S)) and a continuum background. For the resonances [57], an "extended" Crystal Ball (CB2) function was adopted, which accommodates a non-Gaussian tail both on the right and on the left side of the resonance peak. Alternatively, a pseudo-Gaussian function was used, corresponding to a resonance Gaussian core around the J/ψ pole and tails on the right and left side of it, parameterised by varying the width of the Gaussian as a function of the mass. The background was described by empirical functions, either with a Gaussian with a mass-dependent width or with an exponential function times a fourth-order polynomial [57]. Fits were performed using all the combinations of the signal and background functions, and varying the fitting ranges (2.2 < m µµ < 4.5 GeV/c 2 or 2 < m µµ < 5 GeV/c 2 ). Figure 1 shows an example of fits to the invariant mass distributions of the p-Pb and Pb-p data samples, for opposite-sign dimuons in the region p T < 20 GeV/c.
When fitting the mass spectra, the value of the J/ψ mass and its width (σ) at the pole position are free parameters of the fit. The contribution of the ψ(2S) was found to have a negligible impact on the evaluation of N J/ψ .
A study of the influence of the non-Gaussian tails of the shapes of the reconstructed resonance spectra was also performed. The corresponding fit parameters were extracted either from the MC or directly from data. In the latter case, the tail parameters were evaluated either by leaving them as free parameters in the fit to the p-Pb and Pb-p samples, or using values obtained from the corresponding pp data samples at √ s = 8 TeV [58] (about the same energy of the collisions under study) or √ s = 13 TeV [59] (largest data sample collected by ALICE).
The N J/ψ values were finally obtained as the average of the results of all the fits performed. The statistical uncertainties were obtained as the average of the statistical uncertainties over the various fits, while the standard deviations of the N J/ψ distributions were taken as the systematic uncertainties. Typical values of the signal over background ratio in a 3σ window around the J/ψ peak range from 1.4 (0.7) to 2.8 (1.4) from low to high p T in p-Pb (Pb-p) collisions. For the p T -integrated data samples (p T < 20 GeV/c), N J/ψ amounts to (1.67 ± 0.01 ± 0.05) · 10 5 and (2.52 ± 0.01 ± 0.08) · 10 5 for p-Pb and Pb-p respectively, where the first uncertainty is statistical and the second systematic. The latter quantity (which amounts in percentage terms to ∼ 3%) is dominated by the choice of the J/ψ tail parameters. When extracting N J/ψ in narrower p T and y cms ranges, the systematic uncertainties turn out to be similar (from 3% up to 4% in the highest p T bins).
The quantity (A × ε) (∆ycms,∆p T ) was evaluated by means of MC simulations, performed separately for each data taking run, in order to follow the evolution of the detector conditions. The input p T and y cms distributions for the J/ψ were tuned directly to the data by means of an iterative procedure. In detail, a first set of differential distributions, corresponding to the results of the measurements performed at √ s NN = 5.02 TeV [28], was taken as an input to the calculation, and the resulting (A × ε) (∆ycms,∆p T ) values were then used to correct the raw J/ψ distributions obtained from the fits of the invariant mass spectra. The corrected differential distributions were then used as an input for another (A × ε) (∆ycms,∆p T ) calculation, and so on. Convergence was reached at the second iteration.

JHEP07(2018)160
The systematic uncertainties are related to the corresponding uncertainties on the trigger and tracking efficiencies, as well as to the choice of the input distributions. Concerning the efficiencies, for the muon trigger the procedure already used for the analysis of p-Pb data at √ s NN = 5.02 TeV was followed [27,28]. The response function of the muon trigger obtained in MC and in data was used for the J/ψ (A × ε) (∆ycms,∆p T ) calculation. Integrating over p T , a difference of 2.4% (2.9%) on the trigger efficiency for J/ψ was estimated in p-Pb (Pb-p) collisions. The difference can become as high as 4% for low-p T J/ψ. A 1% contribution due to the uncertainty on the intrinsic efficiency of the muon-trigger detectors was then added in quadrature to the quoted uncertainties. For the tracking efficiency, the corresponding systematic uncertainty was calculated by comparing the efficiencies evaluated in data and MC. The efficiency of each tracking plane was obtained using the redundancy of the tracking system (two independent planes per station). Then, the single muon tracking efficiencies were calculated according to the tracking algorithm, and finally combined, in order to get the dimuon tracking efficiency. The estimated value of the systematic uncertainty on the tracking efficiency is 1% (2%) for p T -integrated J/ψ production in p-Pb (Pb-p) collisions, and shows no appreciable dependence on the dimuon kinematics. A further systematic uncertainty, related to the choice of the χ 2 cut applied to the matching of tracks reconstructed in the muon tracking and triggering systems, was also included. Its value is 1%, independent of p T and y cms . Finally, the choice of the MC input distributions was found to induce a 0.5% systematic uncertainty on the acceptance calculation for the p T -integrated data samples. This effect is due to the statistical uncertainty on the measured y cms and p T distributions that were used for the calculation, and to possible correlations between the distributions in the two kinematic variables. The maximum value of this uncertainty becomes 3% at very low p T (<1 GeV/c). The integrated luminosities for the two data samples were obtained from L int = N MB /σ MB where N MB is the number of MB events and σ MB the cross section corresponding to the MB trigger condition. The latter quantity was evaluated from a van der Meer scan, obtaining 2.09 ± 0.03 b for p-Pb and 2.10 ± 0.04 b for Pb-p [55]. The N MB quantity was estimated as N µµ−MB · F norm , where N µµ−MB is the number of analysed dimuon triggers and F norm is the inverse of the probability of having a triggered dimuon in a MB event. F norm was calculated using the event trigger information, as the ratio between the number of collected MB triggers and the number of times the MB condition is verified together with the dimuon trigger condition, with the latter information obtained from the level-0 trigger mask. The F norm values were evaluated, and corrected for the small pile-up contribution to the MB sample (∼ 3% on average), for each run and finally averaged using as a weight the number of µµ − MB triggers. In this way one obtains F pPb norm = 679 ± 7 and F Pbp norm = 371 ± 4. The quoted uncertainties (1%) are systematic and were obtained by comparing the results of the evaluation described above with an alternative method based on the information of the trigger scalers [27]. Statistical uncertainties on F norm are negligible.
The nuclear effects on J/ψ production in p-Pb collisions were estimated via the nuclear modification factor, defined as: R pPb (y cms , p T ) = d 2 σ J/ψ pPb /dy cms dp T A Pb · d 2 σ J/ψ pp /dy cms dp T where the p-Pb production cross section is normalised to the corresponding quantity for pp collisions times the atomic mass number of the Pb nucleus (A Pb = 208). The reference pp cross section was evaluated starting from the available results for forward-y cms inclusive J/ψ production at √ s = 8 TeV from ALICE [58] and LHCb [60]. These results are in fair agreement, as their maximum difference is 1.4σ, in the region close to y cms = 2.5. Since the ALICE pp data cover a different y cms -range (2.5 < y cms < 4) with respect to those accessible in p-Pb and Pb-p collisions, a rapidity extrapolation by ∼ ±0.5 y-units was performed to match the kinematic window of the various samples, following the procedure described in [61]. In addition, a √ s-interpolation [59] was performed to account for the small difference in the centre-of-mass energy between pp and proton-nucleus collisions.
The rapidity extrapolation was performed on the ALICE data using three different functions (Gaussian, 2 nd and 4 th degree polynomials) and taking the weighted average of the extrapolated values. The associated systematic uncertainty was calculated as the maximum difference between the results obtained with the different functions. Typical values are ∼ 2-3%, reaching a maximum of ∼ 25% at the very edge of the extrapolation region. For LHCb, the same procedure was used in order to match the rapidity binning of the p-Pb and Pb-p data. The procedure corresponds in this case to an interpolation, because of the larger rapidity acceptance (2 < y cms < 4.5) of LHCb. Finally, the weighted average of the ALICE/LHCb based extrapolations/interpolations was calculated, and a small correction factor (1.5%), obtained via a √ s-interpolation of data at various centreof-mass energies, was introduced to account for the slight centre-of-mass energy difference between p-Pb ( √ s NN = 8.16 TeV) and pp data ( √ s = 8 TeV).
For the p T -differential studies, the reference pp cross section was obtained as a weighted average of the ALICE and LHCb p T -differential cross sections at √ s = 8 TeV [58,60], extrapolated/interpolated to the proton-nucleus rapidity domains. The ALICE values, which are extrapolated beyond the measured pp rapidity range, were also corrected by p Tdependent factors, which account for the softening/hardening of the p T -differential cross section when y cms increases/decreases, and were calculated from the LHCb pp results on the p T -differential inclusive J/ψ cross section in narrow y cms bins [60]. Since the p T coverage of pp data at √ s = 8 TeV by LHCb extends only up to p T = 14 GeV/c, a linear extrapolation of the correction factors up to p T = 20 GeV/c was performed. The size of this correction is < 10% for p T 6 GeV/c and increases up to ∼ 40% in the highest p T bin. The uncertainty associated with this correction factor is small (1 − 2%), thanks to the very good accuracy of the LHCb results. Finally, the effect of the slight centre-of-mass energy difference between proton-nucleus and pp data sets ranges from 1% to 3.5% when increasing p T . Table 1 summarises the systematic uncertainties on the various contributions entering the cross section and the nuclear modification factor determination. The uncertainty on the integrated luminosity is the sum in quadrature of the uncertainties on σ MB [55] and F norm . The fractions correlated/uncorrelated between p-Pb and Pb-p measurements are separately quoted.   Table 1. Summary of systematic uncertainties on the calculation of cross sections and nuclear modification factors. Uncertainties on signal extraction, MC input and efficiencies are considered as uncorrelated over p T and y cms . The uncertainties on the luminosity and on the pp reference result from the combination of two contributions, one uncorrelated and the other correlated, which are separately quoted. The uncorrelated uncertainty on luminosity includes the contribution of the systematic uncertainty on F norm as well as a 1.1% (0.6%) contribution due to the difference between the luminosities obtained with the V0 and T0 detectors.

Results
In figure 2 the differential cross sections are presented for inclusive J/ψ production as a function of rapidity in p-Pb and Pb-p collisions, integrated over the transverse momentum interval p T < 20 GeV/c. The same figure shows the reference cross sections for pp collisions, obtained through the interpolation procedure described in section 3 and scaled by A Pb . Figure 3 reports the p-Pb differential cross sections as a function of p T , separately for the forward (2.03 < y cms < 3.53) and backward (−4.46 < y cms < −2.96) rapidity regions, where the corresponding pp cross sections scaled by A Pb are also superimposed. The comparison of proton-nucleus and scaled pp cross sections shows that at forward y cms a suppression of the inclusive J/ψ production is visible, while no significant nuclear effects can be seen at backward y cms . Nuclear effects, already visible from the different behaviour of p-Pb and pp-scaled cross sections, are quantified through the nuclear modification factors, shown as a function of y cms in figure 4 and of p T in figure 5. The results are compared with the corresponding nuclear modification factors at √ s NN = 5.02 TeV [28]. Although the √ s NN = 8.16 TeV data are systematically lower, the difference is not significant given the uncertainties of the measurements. As a function of y cms , R pPb decreases when moving from the Pb-going to the p-going direction, showing a significant suppression at forward rapidity, while the negative rapidity measurements do not show any significant deviation from unity. As a function of p T , an increase is seen at forward y cms and the data become compatible with  unity for p T 5 GeV/c. At negative y cms an increasing trend is also likely to be present at low transverse momentum, as shown by a fit in the region p T < 4 GeV/c with a constant function, which gives χ 2 /ndf = 3.3. For p T > 4 GeV/c the nuclear modification factor is systematically larger than 1, but compatible with unity within 1.9σ. Concerning the compatibility of the results at the two energies, the integration over different p T ranges (p In figure 6 the ALICE results are compared to the corresponding LHCb values [47], which cover a slightly wider y cms range and are integrated up to p T = 14 GeV/c, showing a good agreement between the two measurements. The LHCb results refer to prompt J/ψ production, i.e., include decays of higher-mass charmonium states but do not include the contribution from b-hadron decays (non-prompt production). For the region p T ≤ 5 GeV/c, which dominates the p T -integrated results, the size of the latter contribution amounts to 10-15% of the inclusive production. An estimate of the difference between prompt and inclusive nuclear modification factors, based on LHCb results [47], gives a 3 − 4% (1 − 2%) effect at positive (negative) y cms .
In figure 6 a comparison with the results of several theoretical models for prompt J/ψ production is also presented. The results of two calculations based on a pure shadowing scenario (Vogt [62], Lansberg et al. [37,63]) show good agreement with data when the nCTEQ15 [17] or EPPS16 [18] set of nuclear parton distribution functions (nPDF) is adopted, while using the EPS09 [16] set of nPDF leads to a slightly worse agreement at forward y cms . Calculations based on a CGC approach coupled with various elementary production models are able to reproduce the data in their domain of validity, corresponding  to the forward-y cms region (Venugopalan et al. [64], Ducloue et al. [65]). The model of Arleo et al. [66], based on the calculation of the effects of parton coherent energy loss, gives a good description of backward-y cms results and reproduces the data at forward y cms fairly well. Finally, models including a contribution from the final state interaction of the cc pair with the partonic/hadronic system created in the collision (Ferreiro [34], Zhuang et al. [67]) can also reproduce the trend observed in the data. In such a class of models nuclear shadowing is included, and is anyway the process that plays a dominant role in determining the values of the nuclear modification factors. Figure 7 shows a comparison of the p T -dependence of R pPb and R Pbp with the calculations of the models discussed above. Thanks to the extended p T range, these data explore a wide x interval. At y cns = 2.78 (centre of the forward-y interval), the covered range for 0 < p T < 20 GeV/c is 2.3 · 10 −5 < x < 1.5 · 10 −4 while at y cms = −3.71 one has 1.5 · 10 −2 < x < 10 −1 . These values were calculated in the so-called 2 → 1 approach, where the production channel is based on the gluon fusion process gg → J/ψ. The agreement between data and models is rather good. It should be noted that for models that include uncertainty bands, such uncertainties are generally larger than those of the data, both as a function of y cms and p T .
By forming the ratio of the nuclear modification factors at forward and backward rapidity, it is possible to obtain a quantity, R FB , with smaller uncertainties, provided that the same absolute values of the y cms -ranges are chosen for the ratio. In this way, the   Figure 8. The ratio R FB between the inclusive J/ψ nuclear modification factors, as a function of y cms (left) and p T (right), relative to 2.96 < |y cms | < 3.53. The horizontal bars correspond to the bin size. The vertical error bars represent the statistical uncertainties, the boxes around the points the systematic uncertainties. The results are compared with those obtained at √ s NN = 5.02 TeV [27].
The latter have been plotted at a slightly shifted y cms and p T , for better visibility.
reference pp cross section, and the related uncertainties, cancel out. R FB is calculated in the rapidity range 2.96 < |y cms | < 3.53, which is covered by both p-Pb and Pb-p samples. In figure 8 the y cms -and p T -dependence of R FB are shown, and compared with the corresponding results at √ s NN = 5.02 TeV [27]. No appreciable dependence on y cms can be seen, while R FB steadily increases as a function of p T , reaching unity at p T ∼ 12 GeV/c. Results at √ s NN = 8.16 and 5.02 TeV are compatible within uncertainties.

Conclusions
Inclusive J/ψ production in p-Pb collisions at √ s NN = 8.16 TeV was measured by ALICE, with about twice the integrated luminosity of the corresponding data sample at √ s NN = 5.02 TeV [28]. Results on the cross sections and on the nuclear modification factors were shown, in six rapidity bins, for the p-going (2.03 < y cms < 3.53) and Pb-going (−4.46 < y cms < −2.96) directions. The corresponding results as a function of transverse momentum were also shown, separately for the two y cms regions, for p T < 20 GeV/c. A suppression of the J/ψ was observed at positive y cms , concentrated in the p T 5 GeV/c range. For negative y cms , an increasing trend in the nuclear modification factor is present at low p T , and the data are compatible with unity within 1.9σ for p T > 4 GeV/c. The ratios R FB between forward-and backward-y cms R pPb in the region 2.96 < |y cms | < 3.53 were also shown as a function of y cms and p T . The results on the nuclear modification factors and on R FB were found to be compatible with those obtained at √ s NN = 5.