University of Birmingham Measurement of D-meson production versus multiplicity in pPb collisions at sNN = 5 . 02 TeV

The measurement of prompt D-meson production as a function of multiplicity in p–Pb collisions at √ sNN = 5.02 TeV with the ALICE detector at the LHC is reported. D0, D+ and D∗+ mesons are reconstructed via their hadronic decay channels in the centre-of-mass rapidity range −0.96 < ycms < 0.04 and transverse momentum interval 1 < pT < 24 GeV/c. The multiplicity dependence of D-meson production is examined by either comparing yields in p–Pb collisions in different event classes, selected based on the multiplicity of produced particles or zero-degree energy, with those in pp collisions, scaled by the number of binary nucleon-nucleon collisions (nuclear modification factor); as well as by evaluating the per-event yields in p–Pb collisions in different multiplicity intervals normalised to the multiplicity-integrated ones (relative yields). The nuclear modification factors for D0, D+ and D∗+ are consistent with one another. The D-meson nuclear modification factors as a function of the zero-degree energy are consistent with unity within uncertainties in the measured pT regions and event classes. The relative D-meson yields, calculated in various pT intervals, increase as a function of the charged-particle multiplicity. The results are compared with the equivalent pp measurements at √ s = 7 TeV as well as with EPOS 3 calculations.


Introduction
In high-energy hadronic collisions, heavy quarks (charm and beauty) are produced in hard parton scattering processes. Due to their large masses, their production cross sections can be calculated in the framework of perturbative Quantum Chromodynamics (pQCD) down to low transverse momenta. The differential cross section for heavy-flavour hadron production in nucleon-nucleon collisions can be calculated in the factorisation approach by the convolution of parton densities in the incoming nucleon, the short-distance partonic cross section of heavy quark production, and the fragmentation function that describes the transition of the heavy quark into a heavy-flavour hadron [1]. Thus, heavy-flavour production is sensitive to the gluon and the possible heavy-quark content in the nucleon and provides constraints on the parton distribution functions (PDFs) in the proton and in -1 -JHEP08(2016)078 the nucleus [2,3]. Measurements of heavy-flavour hadron production in hadronic collisions provide tests of pQCD and constitute a crucial baseline for the study of heavy-flavour production in heavy-ion collisions [4,5]. A suppression of heavy-flavour yields is observed in heavy-ion collisions at high transverse momentum (p T ), and is interpreted as being due to the formation of a Quark-Gluon Plasma (QGP).
Recently, the study of heavy-flavour production as a function of the multiplicity of charged particles produced in the collision has attracted growing interest. Such measurements probe the interplay between hard and soft mechanisms in particle production. At LHC energies, the multiplicity dependence of heavy-flavour production is likely to be affected by the larger amount of gluon radiation associated with short-distance production processes, as well as by the contribution of Multiple-Parton Interactions (MPI) [23][24][25]. It has also been argued that, due to the spatial distribution of partons in the transverse plane, the probability for MPI to occur in a pp collision increases towards smaller impact parameters [26][27][28]. This effect might be further enhanced by quantum-mechanical fluctuations of gluon densities at small Bjorken-x [29].
The measurements of prompt D mesons, inclusive and non-prompt J/ψ in pp collisions at √ s = 7 TeV [30, 31], and of the three Υ states in pp collisions at √ s = 2.76 TeV [32], provide evidence for a similar increase of open and hidden heavy-flavour yields as a function of charged-particle multiplicity. These results suggest that the enhancement probably originates in short-distance production processes, and is not influenced by hadronisation mechanisms. The enhancement is quantitatively described by calculations including MPI contributions, namely percolation model estimates [33,34], the EPOS 3 event generator [35,36] and PYTHIA 8.157 calculations [37].
In proton-nucleus collisions, several so-called 'Cold Nuclear Matter' (CNM) effects occur due to the presence of a nucleus in the colliding system, and, possibly, to the large density of produced particles. These CNM effects can affect the production of heavy-flavour hadrons at all the stages of their formation. In particular, the PDFs of nucleons bound in nuclei are modified with respect to those of free nucleons. This modification of the PDFs in the nucleus can be described by phenomenological parameterisations (nuclear PDFs, or nPDFs) [38][39][40]. Alternatively, when the production process is dominated by gluons at low Bjorken-x, the nucleus can be described by the Colour-Glass Condensate (CGC) effective theory as a coherent and saturated gluonic system [41][42][43][44]. The kinematics of the partons in the initial state can be affected by multiple scatterings (transverse momentum broadening, or k T broadening) [45][46][47] or by gluon radiation (energy loss) [48] before the heavy-quark pair is produced. Gluon radiation may also occur after the heavy-quark pair JHEP08(2016)078 is formed [49]. Other measurements in p-Pb collisions at √ s NN = 5.02 TeV, e.g. those of angular correlations between charged particles [50][51][52][53], of ψ(2S) suppression [54] and of the relative yields of the three Υ states [32], indicate that final-state effects also play an important role. The measured charm production cross section in minimum-bias p-Pb collisions at √ s NN = 5.02 TeV [55] is consistent within uncertainties with that in pp collisions at the same energy scaled by the atomic mass number of the Pb nucleus. The nuclear modification factor was also found to be consistent with calculations considering EPS09 nPDFs [38], CGC, or transverse momentum broadening and initial-state energy loss. The influence of cold nuclear matter effects on multiplicity-integrated D-meson production in p-Pb collisions is smaller than the measurement uncertainties.
Additional insight into CNM effects can be obtained by measuring the heavy-flavour hadron yields as a function of the multiplicity of charged particles produced in the p-Pb collision. The aim of these studies is to explore the dependence of heavy-flavour production on the collision geometry and on the density of final-state particles. Indeed, it is expected that the multiplicity of produced particles depends on the number of nucleons overlapping in the collision region, and therefore on the geometry of the collision (i.e. on the collision centrality).
Most of the aforementioned models of CNM effects consider a dependence on the collision geometry, usually expressed through the impact parameter of the collision, the number of participant nucleons (N part ), or the number of nucleon-nucleon collisions (N coll ). In general, CNM effects are expected to be more pronounced in central collisions, i.e. those having a small impact parameter. Some of the parameterisations of the nPDFs have studied the influence of the local nucleon density [56][57][58][59]. The spatially dependent EPS09 and EKS98 nPDF sets, EPS09s and EKS98s, are formulated as a function of the nuclear thickness [56]. The leading twist nuclear shadowing calculation [60] assumes the Glauber-Gribov approach of the collision geometry and predicts the dependence of the nPDF on the collision impact parameter. The estimates of the initial-state k T broadening due to multiple soft collisions also consider a dependence on the collision impact parameter [46,47]. Initial-state parton energy loss is also expected to evolve with the collision geometry as a consequence of the different nuclear density, though detailed calculations including this effect are not yet available. Finally, if final-state effects were to affect heavy-flavour production in p-Pb collisions, their influence would also vary with the density of produced particles.
In this paper, we report the p T -differential measurements of D 0 , D + and D * + production as a function of multiplicity in p-Pb collisions at √ s NN = 5.02 TeV. The experimental setup and the data sample are described in section 2. The determination of the multiplicity and the estimation of the collision centrality and of the number of nucleon-nucleon collisions are discussed in section 3. The D-meson reconstruction strategy is explained in section 4. The results are reported in the form of the D-meson nuclear modification factor in different centrality classes (section 5), and the relative D-meson yields as a function of the relative charged-particle multiplicity at central and backward rapidity (section 6).

JHEP08(2016)078
2 Experimental apparatus and data sample The ALICE apparatus is described in detail in [61] and its performance in [62]. It is composed of a series of detectors in the central barrel for tracking and particle identification; the Muon Spectrometer in the forward direction for muon tracking and identification; and a further set of detectors at forward rapidity for triggering and event characterisation. The central barrel detectors are located inside a large solenoid magnet that provides a 0.5 T field parallel to the beam direction, which corresponds to the z-axis of the ALICE coordinate system. In this section, the detectors used for the D-meson analysis are briefly described. The Inner Tracking System (ITS), the Time Projection Chamber (TPC) and the Time Of Flight detector (TOF) allow the reconstruction and identification of charged particles in the central pseudorapidity region. The V0 detector, composed of two scintillator arrays located in the forward and backward pseudorapidity regions, is used for online event triggering and multiplicity determination. The Zero Degree Calorimeters (ZDC) are used for event selection and to estimate the collision centrality via the zero-degree energy.
The ITS is composed of six cylindrical layers of silicon detectors, located at radii between 3.9 cm (about 1 cm from the beam vacuum tube) and 43.0 cm. The two innermost layers, which respectively cover |η| < 2.0 and |η| < 1.4, comprise the Silicon Pixel Detectors (SPD); the two intermediate layers, within |η| < 0.9, consist of Silicon Drift Detectors (SDD); and the two outer layers, also covering |η| < 0.9, consist of double-sided Silicon Strip Detectors (SSD). The low material budget, high spatial resolution, and position of the detector setup surrounding the beam vacuum tube and close to the interaction point allow it to provide a measurement of the charged-particle impact parameter in the transverse plane (d 0 ), i.e. the distance of closest approach between the track and the primary vertex along rφ, with a resolution better than 75 µm for transverse momenta p T > 1 GeV/c [63].
The TPC is a large cylindrical drift detector, extending from 85 cm to 247 cm in the radial direction and covering the range −250 < z < +250 cm along the beam axis [64]. It provides charged-particle trajectory reconstruction with up to 159 space points per track in the pseudorapidity range |η| < 0.9 and in the full azimuth. The primary interaction vertex position and covariance matrix are determined from tracks reconstructed from hits in the TPC and the ITS via a χ 2 analytic minimisation method.
The TOF detector is equipped with Multi-gap Resistive Plate Chambers (MRPCs) [62]. It is placed at radii between 377 cm and 399 cm, and has the same pseudorapidity and azimuthal coverage as the TPC. The TOF measures the flight times of charged particles from the interaction point to the detector with an overall resolution of about 85 ps. For events with the 20% lowest multiplicities, the resolution decreases to about 120 ps due to a worse start-time (collision-time) resolution. The start-time of the event is determined by combining the time estimated using the particle arrival times at the TOF and the time measured by the T0 detector, an array of Cherenkov counters located at +350 cm and −70 cm along the beamline. Particle identification (PID) is performed by comparing the measurement of the specific energy deposition dE/dx in the TPC and the time-of-flight information from the TOF with the respective expected values for each mass hypothesis.
The V0 detector consists of two arrays of scintillator tiles covering the pseudorapidity regions −3.7 < η < −1.7 (V0C) and 2.8 < η < 5.1 (V0A) [65]. The data sample analysed -4 -JHEP08(2016)078 in this paper was collected with a minimum-bias interaction trigger requiring at least one hit in both V0A and V0C counters coincident with the arrival time of the proton and lead bunches. The ZDC is composed of two sets of neutron (ZNA and ZNC) and proton (ZPA and ZPC) calorimeters positioned on either side of the interaction point at z = ±112.5 m. Contamination from beam-background interactions was removed via offline selections based on the timing information provided by the V0 and the ZNA. The signals registered by the SPD and V0 detectors were used to determine the event charged-particle multiplicity; the SPD, V0 and ZDC detectors were also exploited to classify the events in centrality classes, as will be described in section 3.
The data sample used in this paper was recorded in January 2013, during the p-Pb LHC run. Protons with an energy of 4 TeV were collided with Pb ions with an energy of 1.58 TeV per nucleon, resulting in collisions at a centre-of-mass energy per nucleon pair, √ s NN , of 5.02 TeV. With this beam configuration, the centre-of-mass system moves with a rapidity of ∆y cms = 0.465 in the direction of the proton beam, due to the different energies per nucleon of the proton and the lead beams. In the case of the D-meson analyses presented here, performed in the laboratory reference interval |y lab | < 0.5, this leads to a shifted centre-of-mass rapidity coverage of −0.96 < y cms < 0.04. In the following, we will use the notation η and y lab to refer to the pseudorapidity and rapidity values in the laboratory reference frame, and η cms and y cms for the values evaluated in the centre-of-mass reference frame. A total of 10 8 minimum-bias triggered events, corresponding to an integrated luminosity of L int = 48.6 ± 1.6 µb −1 , passed the selection criteria and were analysed.

Multiplicity determination
The production of D mesons in p-Pb collisions has been studied as a function of chargedparticle multiplicity using two different observables.
One observable is the p T -differential nuclear modification factor, which is defined as the ratio of the p T -differential yields measured in p-Pb collisions in centrality intervals to those in pp collisions, scaled by the number of binary nucleon-nucleon collisions. The centrality intervals were defined using three different estimators based on the multiplicity in the SPD and V0A detectors and the energy deposited in the zero-degree neutron calorimeter in the Pb-going side (ZNA). The procedure used to determine the number of binary nucleonnucleon collisions for each event class is described in section 3.1 and [66].
The other observable, referred to as the relative yield, is defined as the ratio of the per-event D-meson yields in p-Pb collisions in different multiplicity intervals normalised to the multiplicity-integrated yields. Details on the evaluation of the charged-particle multiplicity are discussed in section 3.2. In this analysis, the values of multiplicity measured in two different pseudorapidity intervals, namely at mid-rapidity with the SPD and at large rapidity in the Pb-going direction with the V0A, were considered.

Centrality estimators and T pPb determination
A centrality-dependent measurement of the nuclear modification factor requires the p-Pb data sample to be sliced into classes according to an experimental observable related to -5 -JHEP08(2016)078 the collision centrality, as well as a determination of the average nuclear overlap function T pPb , which is proportional to the number of nucleon-nucleon collisions N coll , for each centrality class.
The minimum-bias p-Pb data sample was divided into four centrality classes by exploiting the information from: (i) V0A, the amplitude of the signal measured by the V0 scintillator array located in the Pb-going side, covering 2.8 < η < 5.1, which is proportional to the number of charged particles produced in this pseudorapidity interval; (ii) CL1, the number of clusters in the outer layer of the SPD, covering |η| < 1.4, which is proportional to the number of charged particles at mid-rapidity; and (iii) ZNA, the energy deposited in the Zero Degree Neutron Calorimeter positioned in the Pb-going side by the slow nucleons produced in the interaction by nuclear de-excitation processes, or knocked out by wounded nucleons. The multiplicity of these neutrons is expected to grow monotonically with the number of binary collisions, N coll .
Centrality classes were defined as percentiles of the visible cross section, which was measured to be (2.09 ± 0.07) b [67]. For the centrality classes defined using the CL1 and V0A multiplicities, a Glauber Monte Carlo was used to calculate the relevant geometrical quantities, namely the average numbers of participant nucleons N Glauber part , of binary collisions N Glauber coll , and the average nuclear overlap function T Glauber pPb [66]. For the case where the ZNA information was used, the values of N part , N coll and T pPb were obtained using the so-called hybrid method [66]. In this approach, the determination of T pPb in a given ZNA-energy class relies on the assumption that the charged-particle multiplicity measured at mid-rapidity (−1 < η cms < 0) scales with the number of participant nucleons, N part .
1) where N MB part = 7.9 is the average number of participants in minimum-bias collisions and σ NN = (70 ± 5) mb is the interpolated inelastic nucleon-nucleon cross section at √ s NN = 5.02 TeV [66]. The values of T pPb obtained with the three estimators in the four multiplicity (zero-degree energy) classes used for the analysis are reported in table 1.
It was demonstrated by the studies of charged-particle production reported in [66] that when centrality classes are defined in p-Pb collisions, some biases are present. Firstly, there is a multiplicity selection bias due to the large multiplicity fluctuations for p-Pb interactions at a given impact parameter, which are comparable in magnitude to the full dynamic range of the minimum-bias multiplicity distribution. In addition, there is a jet-veto bias due to the contribution to the overall multiplicity from particles arising from the fragmentation of partons produced in hard-scattering processes. This causes low-(high-) multiplicity p-Pb collisions to correspond to a lower (higher) number of hard scatterings per nucleonnucleon collision. Furthermore, a purely geometrical bias was suspected to affect peripheral collisions for all centrality estimators, due to the fact that the mean impact parameter between the proton and each nucleon of the Pb nucleus, calculated from a Monte Carlo Glauber simulation, rises significantly for N part < 6, thus reducing the average number of multi-parton interactions for peripheral collisions. 60-100 0.037 0.037 23 0.046 6.2 Table 1. T pPb values in p-Pb collisions at √ s NN = 5.02 TeV obtained with a Glauber-model based approach for V0A and CL1, and from the hybrid method for ZNA, as described in [66].
These biases cause the nuclear modification factor of charged particles to differ from unity in the centrality classes even in the absence of nuclear effects. These biases decrease with increasing rapidity separation between the centrality estimator and the region where the nuclear modification factor is measured. A strong selection bias is observed for the CL1 estimator, due to the full overlap with the tracking region, which is reduced with the V0A estimator. By contrast, the selection based on the energy deposited in the ZNA is expected to be free from the biases related to the event selection, and is only affected by the geometrical bias.
For these reasons, the results based on the ZNA selection, which is the least biased [66], provide insight into possible centrality-dependent nuclear effects on charm production in p-Pb collisions. Moreover, the measurements of the D-meson nuclear modification factor in centrality intervals defined with the three estimators described above offer the possibility to study these biases based on heavy-flavour production, which, due to the large mass of the charm quarks, is expected to scale with the number of binary collisions over the whole p T range, provided that cold nuclear matter effects are negligible. This is in contrast to the charged-particle yield, where a scaling with N coll is expected to occur only in the high-p T region.

Relative event multiplicity determination
The charged-particle multiplicity, N ch , was estimated at mid-rapidity by measuring the number of tracklets, N tracklets , reconstructed in the SPD. A tracklet is defined as a track segment that joins a pair of space points on the two SPD layers and is aligned with the reconstructed primary vertex. N tracklets was counted within |η| < 1.0.
The pseudorapidity acceptance of the SPD depends on the position of the interaction vertex along the beam line z vtx , both due to the asymmetry of the collision system and the limited coverage of the detector. In addition, the overall SPD acceptance varies as a function of time due to a varying number of active channels. A data-driven correction was applied to the N tracklets distributions on an event-by-event basis to account for these two effects. This was done by renormalising the N tracklets distributions to the overall minimum with a Poissonian smearing to account for the fluctuations. Multiplicity classes were then defined based on the percentiles of analysed events in each N tracklets range.

JHEP08(2016)078
The conversion of N tracklets to N ch was performed using minimum-bias Monte Carlo simulations. The distribution of the measured N tracklets as a function of the number of generated "physical primaries" (N ch ) in the simulation was considered for this purpose. Physical primaries are defined as prompt particles produced in the collision and their decay products, excluding those from weak decays of strange particles. The proportionality factor was evaluated from a linear fit to the distribution, and was then applied to the mean N tracklets in each interval to give the estimated N ch values. These values were then divided by the width of the considered η range, ∆η = 2, to give an estimated dN ch /dη. The uncertainty of the N tracklets to N ch conversion was estimated by testing its deviation from linearity. A linear fit to the distribution was performed in each multiplicity interval to evaluate the possible changing slope of the distribution between intervals. From these fits, a series of scaling factors were obtained and compared to the multiplicity-integrated one, resulting in a 5% uncertainty.
The results are given as a function of the relative charged-particle multiplicity, (dN ch /dη)/ dN ch /dη , where dN ch /dη = 17.64 ± 0.01 (stat.) ± 0.68 (syst.) was measured by ALICE for inelastic p-Pb collisions at √ s NN = 5.02 TeV with at least one charged particle within |η| < 1.0 [68]. The N tracklets ranges considered in this analysis, and the corresponding relative multiplicity values, are given in table 2.
The production of D mesons was also studied as a function of charged-particle multiplicity in the region 2.8 < η < 5.1, as measured with the signal amplitude in the V0A detector, N V0A , reported in units of the minimum-ionising-particle charge. This estimator allows the multiplicity and the D-meson yields to be evaluated in two different pseudorapidity intervals (backward and central η), avoiding possible auto-correlations.
The average N V0A depends on z vtx , due to the varying distance between the primary vertex and the detector array. This effect was corrected with the same method used for the N tracklets case, leading to an overall average N V0A of 82.7. In this case, the results are considered as a function of the V0A multiplicity relative to the mean multiplicity in the same rapidity region, rather than performing a conversion to dN ch /dη. The N V0A intervals considered, and the corresponding relative multiplicity intervals, are reported in table 3.
It should be noted that the analyses performed as a function of centrality examine the events in samples populated by 20% of the analysed events (40% for the most peripheral events, see table 1), whereas those performed as a function of charged-particle multiplicity explore events from low to extremely high multiplicities, corresponding to about 60% and 5% of the analysed events, respectively (see tables 2 and 3). For the latter analyses, the event classes were defined to study the D-meson yield at extreme multiplicities.

D meson reconstruction
The D 0 , D + , and D * + mesons were reconstructed via their hadronic decay channels D 0 → K − π + (with a branching ratio, BR, of 3.88±0.05%), D + → K − π + π + (BR of 9.13±0.19%), and D * + → D 0 π + (BR of 67.7 ± 0.05%) followed by D 0 → K − π + , and their corresponding charge conjugates [69]. The D 0 and D + weak decays, with mean proper decay lengths (cτ ) of about 123 and 312 µm, respectively, were selected from reconstructed secondary  Table 2. Summary of the multiplicity intervals at central rapidity used for the analyses. The number of reconstructed tracklets N tracklets , the average charged-particle multiplicity dN ch /dη (uncertainty of 5% not quoted), and the relative charged-particle multiplicity (dN ch /dη)/ dN ch /dη (uncertainty of 6.3% not quoted) are listed (see section 6.1 for the uncertainties description). The number of events analysed for the D 0 -meson analysis is also reported for each multiplicity range.  Table 3. Summary of the multiplicity intervals at backward rapidity used for the analyses. The V0A signal N V0A intervals and the relative multiplicity (N V0A )/ N V0A (uncertainty of 5% not quoted) are listed (see section 6.1 for the uncertainties description). The number of events analysed for the D 0 -meson analysis is also reported for each multiplicity range.
vertices separated by a few hundred microns from the interaction point. The D * + meson decays strongly at the primary vertex, and the decay topology of the produced D 0 was reconstructed along with a soft pion originating at the primary vertex. Events were selected by requiring a primary vertex within ±10 cm from the centre of the detector along the beamline. An algorithm to detect multiple interaction vertices was used to reduce the pile-up contribution. D 0 and D + candidates were defined using pairs or triplets of tracks with the proper charge sign combination, within the fiducial acceptance |η| < 0.8 and with transverse momentum p T > 0.3 GeV/c. Only good quality tracks were considered in the combinatorics by requiring selection criteria as described in [19,20,55]. The selection of tracks with |η| < 0.8 reduces the D-meson acceptance, which drops steeply to zero for |y lab | > 0.5 at low p T and for |y lab | > 0.8 at p T > 5 GeV/c. Therefore, a p T -dependent fiducial acceptance region was defined, as reported in [19,20,55].
The selection strategy of the D-meson decay topology was based on the displacement of the decay tracks from the interaction vertex, the separation between the secondary and primary vertices, and the pointing angle, defined as the angle between the reconstructed -9 -JHEP08(2016)078 D-meson momentum and its flight line (the vector between the primary and the secondary vertices). The cuts on the selection variables were chosen in order to obtain a large statistical significance of the D-meson signals, as well as an as large as possible selection efficiency. Therefore, the cut values depend on the D-meson p T and species. In the case of the analysis of the relative yields as a function of multiplicity, the same selections were used in all multiplicity intervals in order to minimise the effect of the efficiency corrections on the ratio of the yields in the multiplicity intervals to the multiplicity-integrated ones. On the other hand, for the analysis of the nuclear modification factor in different centrality classes, the cut values were optimised in each centrality class. Particle identification criteria were applied on the decay tracks, based on the TPC and TOF detector responses, in order to obtain a further reduction of the combinatorial background as explained in [19,20,55].
The raw D-meson yields, both multiplicity-integrated and in each multiplicity or centrality class, were extracted in the considered p T intervals by means of a fit to the invariant mass (M ) distributions of the selected candidates (for the D * + meson the mass difference distributions ∆M = M (Kππ) − M (Kπ) were used). The fit function is the sum of a Gaussian to describe the signal and a function describing the background shape, which is an exponential for D 0 and D + and a threshold function multiplied by an exponential where M π is the pion mass and a and b are free parameters) for the D * + . The centroids and the widths of the Gaussian functions were found to be in agreement with the world average D-meson masses and the values obtained in simulations, respectively, in all multiplicity, centrality and p T intervals. In particular, the widths of the Gaussian functions are independent of multiplicity (or centrality) and increase with increasing D-meson p T . In the relative yield analysis, in order to reduce the effect of the statistical fluctuations, the fits were performed by fixing the Gaussian centroids to the world average D-meson masses, and the widths to the values obtained from a fit to the invariant mass distribution in minimum-bias events, where the signal statistical significance is larger. Figure 1 shows the D 0 and D + invariant mass, and D * + mass difference distributions in the 2 < p T < 4 GeV/c, 4 < p T < 6 GeV/c, 6 < p T < 8 GeV/c intervals, respectively, for the 0-20% and 60-100% centrality classes defined with the ZNA estimator. The fits to the invariant mass distributions were repeated under different conditions and the raw yields were extracted by using alternative methods in order to determine the systematic uncertainties related to the extraction of the raw D-meson counts. The fits were performed by varying the invariant mass ranges and bin widths of the histograms, and considering different functions to describe the background, namely parabolic or linear functions. The raw yields were also obtained by counting the entries of the histograms within a 3σ interval centred on the peak position, after the subtraction of the background estimated from a fit to the side bands, far away from the D-meson peaks.
The raw counts of D mesons extracted in each p T and multiplicity interval were corrected for the acceptance and the reconstruction and selection efficiency. The correction factor for each D-meson species was obtained by using Monte Carlo simulations. Events containing a cc or bb pair were generated by using the PYTHIA v6.4.21 event generator [70] with the Perugia-0 tune [71] and adding an underlying event generated with -10 -JHEP08(2016)078      candidates and of the mass difference for D * + candidates (right column) in two centrality classes defined with the ZNA estimator: 0-20% and 60-100%. The red lines in each plot represent the fit to the background, and the blue lines represent the sum of signal and background. One p T interval is shown for each meson species: 2 < p T < 4 GeV/c for D 0 , 4 < p T < 6 GeV/c for D + , and 6 < p T < 8 GeV/c for D * + .
HIJING v.1.36 [72]. Detailed descriptions of the detector response, the geometry of the apparatus and the conditions of the luminous region were included in the simulation. The generated D-meson p T distribution was tuned in order to reproduce the FONLL [16] spectrum at √ s = 5.02 TeV. The reconstruction and selection efficiency depends on the multiplicity of charged particles produced in the collision, since the primary vertex resolution and the resolution on the topological selection variables improve at high multiplicity. The generated events were weighted on the basis of their charged-particle multiplicity in order to match the multiplicity distribution observed in the data. The reconstruction and selection efficiency depends on the D-meson species and on p T . For prompt D 0 mesons it is about 1-2% in the 1 < p T < 2 GeV/c interval, where the selection criteria are more stringent due to the higher combinatorial background, and it increases to 20% in 12 < p T < 24 GeV/c. The efficiency for D mesons from B decays is higher because the decay vertices of feed-down D mesons are more displaced from the primary vertex and they are more efficiently selected by the topological selections. The efficiencies are slightly larger at high multiplicity, by about 4-10%.
The D-meson raw yields have two components: the prompt D-meson contribution (produced in the charm quark fragmentation, either directly or through strong decays of excited open charm states) and the feed-down contribution originating from B-meson decays. The yield of D mesons from B decays was subtracted from the raw counts by applying a correction factor, f prompt , which represents the fraction of promptly produced -11 -JHEP08(2016)078 D mesons. The f prompt factor was evaluated using the B-hadron production cross section obtained from the FONLL pQCD calculation [16][17][18], the B → D + X kinematics from the EvtGen package [73], and the acceptance times efficiency for D mesons from B decays obtained from the Monte Carlo simulations [19]. The value of f prompt depends on the nuclear modification factor, R feed-down pPb , of the feed-down D mesons. This quantity is related to the nuclear modification of beauty production, which has not been measured in the p T interval of these analyses. Therefore, the nuclear modification factor of feed-down D mesons was assumed to be equal to that of prompt D mesons, R feed-down pPb = R prompt pPb , and a systematic uncertainty was assigned considering the variation 0.9 < R feed-down pPb /R prompt pPb < 1.3. These assumptions were based on the study of the possible modification of the Bhadron production due to the modification of the PDFs in the nucleus through either CGC or pQCD calculations with the EPS09 parameterisation of the nPDFs [38,41].

Nuclear modification factor as a function of centrality
The nuclear modification factor of prompt D 0 , D + and D * + mesons was studied as a function of p T using the three different centrality estimators introduced in section 3.1, based on different measurements of the centrality in terms of multiplicity (CL1 and V0A estimators) or zero-degree energy (ZNA estimator). For each estimator, the analysis of D-meson production was carried out in four event classes, and the nuclear modification factor was calculated as: where (dN D /dp T ) cent pPb is the yield of prompt D mesons in p-Pb collisions in a given centrality class, (dσ D /dp T ) pp is the cross section of prompt D mesons in pp collisions at the same √ s, and T pPb is the average nuclear overlap function in a given centrality class, which was estimated with the Glauber-model approach for the CL1 and V0A estimators (T Glauber pPb ) and with the hybrid method for the ZNA estimator (T mult pPb ) (see section 3.1). In contrast to the multiplicity-integrated R pPb = (dσ D /dp T ) pPb / A · (dσ D /dp T ) pp , Q pPb is influenced by potential biases in the centrality estimation that are not related to nuclear effects, as explained in section 3.1. Hence, Q pPb may be different from unity even in the absence of nuclear effects, in particular if measured with respect to the CL1 and V0A estimators. Complementary to this, the measurement of Q pPb with the ZNA estimator allows the least biased estimation of the possible centrality-dependent modification of the p T -differential yields in p-Pb collisions with respect to the binary-scaled yields in pp collisions.
The cross sections of prompt D-meson production in pp collisions at √ s = 5.02 TeV were obtained by a pQCD-based energy scaling of the p T -differential cross sections measured at √ s = 7 TeV with the scaling factor evaluated by the ratio of the FONLL [16][17][18] calculations at 5.02 and 7 TeV [74]. The scaling procedure was validated by comparing the D-meson p T -differential cross sections at 2.76 TeV with the 7 TeV data scaled down to 2.76 TeV [20]. In the case of D 0 mesons, some refinements were considered for the lowest and highest p T intervals. For 1 < p T < 2 GeV/c, where the D 0 cross section was -12 -JHEP08(2016)078 measured at both 7 and 2.76 TeV [4,19], both measurements were scaled to 5.02 TeV and averaged using the inverse squared of their relative statistical and systematic uncertainties as weights. Since the ALICE measurements of the D 0 cross section in pp data are limited to p T < 16 GeV/c, the estimate for 16 < p T < 24 GeV/c was determined by extrapolating the 7 TeV cross section to higher p T using the FONLL p T -differential spectrum normalised to the measurement in 5 < p T < 16 GeV/c, and scaling it down to 5.02 TeV.
The raw numbers of D mesons in each p T and centrality interval were extracted and corrected by the acceptance and efficiency obtained from Monte Carlo simulations, as described in section 4. The feed-down from B-hadron decays was subtracted from the extracted yields by calculating f prompt in each centrality class independently, as described in section 4.

Systematic uncertainties
The systematic uncertainties (yield extraction, reconstruction and selection efficiency determination and feed-down subtraction) do not depend on the estimator used to define the centrality classes. A mild dependence of the uncertainty on the multiplicity that populates the different centrality classes was observed, resulting in slightly larger uncertainties in the event class with the lowest multiplicity.
The systematic uncertainty of the yield extraction procedure was estimated by varying the fit conditions and by using the bin counting method as introduced in section 4. It is about 3-4% at intermediate p T (2 < p T < 6 GeV/c) and increases to 8-10% at p T < 2 GeV/c and p T > 6 GeV/c. For the D 0 meson, the yield extraction systematic uncertainty includes the contribution to the raw yield of signal candidates reconstructed by assigning the wrong mass to the final state hadrons (about 3-4% for all p T intervals) [55].
The influence of the tracking efficiency was estimated by varying the track selection criteria. The corresponding uncertainty was found to be about 3% per track, resulting in a total uncertainty of 6% (9%) for a two-(three-)particle decay. The uncertainty due to the Dmeson candidate selection criteria was evaluated by varying the topological selections used. It was estimated to be 10% for the interval 1 < p T < 2 GeV/c and 5% for p T > 2 GeV/c.
The effect of the generated D-meson p T shape used to compute the efficiency was estimated by comparing the efficiency values obtained with the PYTHIA and the FONLL p T spectra. A systematic uncertainty of 2-3% was applied only in the interval 1 < p T < 2 GeV/c due to this. The uncertainty due to the multiplicity dependence of the reconstruction and selection efficiency was evaluated changing the weight functions used to reproduce the measured charged-particle multiplicity in the simulations. The multiplicity weights were determined by the ratio of the distribution of the number of tracklets within |η| < 1 in data and Monte Carlo. The weights were computed for: (i) all events selected in the analysis, (ii) events with a D-meson candidate within approximately ±10σ of the invariant mass peak, and (iii) events with a D-meson candidate in the ±3σ invariant mass region. A deviation of about 10% is observed for D mesons at low p T . For high-p T D mesons (p T > 12 GeV/c), the weights have a smaller effect on the efficiency determination, introducing a difference of only 4%.

JHEP08(2016)078
The analysis was repeated without applying the particle identification selections to the D-meson decay hadrons. The corrected yields were consistent, within statistical fluctuations, with those calculated considering particle identification selections. Therefore, no corresponding uncertainty was assigned.
The systematic uncertainty due to the subtraction of feed-down D mesons from B decays was estimated by considering the FONLL uncertainties on the normalisation and factorisation scales and using a second subtraction method based on the ratio of FONLL calculations for D-and B-meson cross sections [19]. The magnitude of this systematic uncertainty depends on the meson species and on the p T interval considered in the measurement, since it is related to the topological selections applied in each analysis. As explained in section 4, a variation of the feed-down D-meson nuclear modification factor was also taken into account as part of the systematics. The quadratic sum of the two contributions to the Q pPb was found to range from a few percent up to 30%.
The denominator of the Q pPb has an uncertainty on the T pPb , which is reported in table 1, and an uncertainty on the pp reference. The latter has a contribution coming from the 7 TeV measurement (ranging from 15% up to 25%) and one from the scaling factor ranging from +17% −4% at p T = 1 GeV/c to ±3% for p T > 8 GeV/c. The uncertainty on the energy scaling factor was estimated by varying the calculation parameters as described in [74]. A larger uncertainty for D 0 in 16 < p T < 24 GeV/c was quantified due to the extrapolation procedure explained above; in that case the uncertainty is +17.5% −4% . The global Q pPb uncertainties were determined by adding the pp and p-Pb uncertainties in quadrature, except for the branching ratio uncertainty, which cancels out in the ratio, and the feed-down contribution, which partially cancels out.

Results
The nuclear modification factors of D 0 , D + and D * + mesons were calculated according to eq. (5.1) in four centrality classes (0-20%, 20-40%, 40-60% and 60-100%) defined with the ZNA estimator, and applying the hybrid method to obtain the T pPb in each class. Figure 2 illustrates these results for 0-20% and 40-60% centrality classes. The Q pPb of the three D-meson species were found to be consistent with one another within the statistical and systematic uncertainties for each p T and centrality class considered. Therefore, the average of the D 0 , D + and D * + meson results was evaluated in each centrality class considering the inverse square of the relative statistical uncertainties as weights. The systematic uncertainties on the averages were computed considering the tracking efficiency, the B feed-down subtraction and the scaling of the pp reference as correlated uncertainty sources among the three mesons. The averages of the D 0 , D + and D * + p T -differential nuclear modification factors in different centrality classes obtained with the ZNA estimator are presented in figure 3 and table 4. The D-meson Q pPb results in the different centrality classes are consistent with unity within the uncertainties in the measurement p T interval. Typical values of the Q pPb uncertainties are of 7% (stat.) and 16% (syst.) for 2 < p T < 4 GeV/c. It should be noted that with this centrality estimator no bias is expected due to the event selection, and only a small bias in peripheral events, due to the geometrical bias in the determination of the number of hard scatterings, was observed in the studies with charged   particles [66]. Therefore, with the least biased centrality estimator, the D-meson Q pPb results are consistent within statistical and systematic uncertainties with binary collision scaling of the yield in pp collisions, independent of the geometry of the collision.

Q pPb with CL1 and V 0A estimators
As explained in section 3.1, the D 0 , D + and D * + Q pPb were also calculated with the CL1 and V0A estimators in four centrality classes to study the centrality selection biases based on heavy-flavour production from low to high p T . The Q pPb results for the three D-meson species were found to be consistent with one another within the statistical and systematic uncertainties for each p T and centrality class considered. Therefore, the averages of the D 0 , D + and D * + meson results and the systematic uncertainties were evaluated as explained before. The averages of the p T -differential D 0 , D + and D * + nuclear modification factors in different centrality classes with CL1 and V0A estimators are presented in figure 4 (see also tables 5 and 6).
The centrality estimation from the CL1 multiplicity suffers from a large bias introduced by multiplicity fluctuations in the central rapidity region caused by fluctuations of the number of hard scatterings per nucleon collision, which affect the T pPb determination [66]. The Q CL1 pPb results show an ordering from low (60-100%) to high (0-20%) multiplicity, with a difference larger than a factor of two between the most central and most peripheral classes, induced by the bias on the centrality estimator.
The V0A estimator classifies the events as a function of the multiplicity in the backward rapidity region. The rapidity gap with respect to the central rapidity D-meson analyses      . Average D 0 , D + and D * + meson Q pPb as a function of centrality with the CL1, the V0A and the ZNA estimators for (a) 2 < p T < 4 GeV/c and (b) 8 < p T < 12 GeV/c. The average D-meson Q pPb in 8 < p T < 12 GeV/c is compared with the charged-particle Q pPb calculated for p T > 10 GeV/c [66]. The vertical error bars and the empty boxes represent, respectively, the statistical and systematic uncertainties on the D-meson results. The filled boxes at Q pPb = 1 indicate the correlated systematic uncertainties: the grey-filled box represents the uncertainty on the pp reference and the p-Pb analysis PID and track selection uncertainties, common to all estimators for a given p T interval; the red-filled box represents the correlated systematic uncertainty on N coll determination for the ZNA energy estimator. similar qualitative behaviour to the Q CL1 pPb ones, with a smaller difference between centrality classes. This is consistent with the expectation of a smaller bias when there is a rapidity gap between the regions where the centrality and the D-meson yield are studied.

Comparison with charged-particle Q pPb
The average D-meson Q pPb results obtained with the three estimators, for 2 < p T < 4 GeV/c and 8 < p T < 12 GeV/c, are displayed as a function of centrality in figure 5. The D-meson Q pPb for 8 < p T < 12 GeV/c is compared with the analogous measurement for charged hadrons with p T > 10 GeV/c [66]. In this transverse momentum region also the production of charged hadrons is expected to scale with the number of binary nucleon-nucleon collisions [66]. The measured trends of charged-particle Q pPb at high p T in all the CL1 and V0A centrality classes were found to be reasonably described by an incoherent superposition of N coll pp collisions generated with PYTHIA, after defining the event centrality from the charged-particle multiplicity in the rapidity region covered by each estimator in the same way as in data (|η| < 1.4 for CL1, 2.8 < η < 5.1 for V0A) [66].
The Q pPb results for D mesons and charged hadrons with p T > 10 GeV/c show a similar trend as a function of centrality and estimator due to the bias in the centrality determination, as observed in [66] based on high-p T particle production in the light flavour -17 -

JHEP08(2016)078
sector. The results presented in this paper allow these studies to be extended into the charm sector and down to low p T .
6 Relative yields as a function of multiplicity D 0 , D + and D * + meson yields were also studied as a function of the charged-particle multiplicity in two pseudorapidity intervals, see section 3.2. The D-meson yields were evaluated for various multiplicity and p T intervals and the results are reported in terms of corrected per-event yields normalised to the multiplicity-integrated values where the index j identifies the multiplicity interval, N j raw D is the raw yield extracted from the fit to the invariant mass distribution in each multiplicity interval, j prompt D represents the reconstruction and selection efficiencies for prompt D mesons, and N j events is the number of events analysed in each multiplicity interval. The efficiencies were estimated with Monte Carlo simulations (see section 4). Equation (6.1) holds under the assumption that the relative contribution to the raw D-meson yield due to the feed-down from beauty-hadron decays does not depend on the multiplicity of the event, and therefore cancels out in the ratio to the multiplicity-integrated values. This assumption is justified by the beauty production measurements as a function of multiplicity in pp collisions, and also by PYTHIA simulations [31]. The acceptance correction, defined as the fraction of D mesons within a given rapidity and p T interval that decay into pairs or triplets of particles within the detector coverage, cancels out in this ratio. The number of events used for the normalisation of the multiplicity-integrated yield must be corrected for the fraction of non-single diffractive events that are not accepted by the minimum-bias trigger condition, expressed as N MB trigger / MB trigger with MB trigger = (96.4±3.1)% [67]. It was verified with PYTHIA 6.4.21 Monte Carlo simulations that the minimum-bias trigger is 100% efficient for D mesons in the kinematic range of the measurement, meaning that the number of D mesons in the minimum-bias triggered events is the same as in the sample of non-single diffractive events.

Systematic uncertainties
In this section the systematic uncertainties estimated for the D-meson measurements as a function of N tracklets and as a function of the N V0A multiplicity are outlined.
The most significant source of systematic uncertainty is the one related to the signal extraction procedure. The raw D-meson yields were obtained by fixing the position of the Gaussian signal peak to the world averages of the D-meson masses, and the widths to the values obtained from the fit to the multiplicity integrated invariant mass distributions. To estimate the yield extraction uncertainty the fit parameters were varied as described in section 4. In addition to the variations listed in section 4, the fits were performed also allowing the position and the width of the Gaussian terms to remain free in the individual JHEP08(2016)078 multiplicity intervals. The yield extraction uncertainty was estimated based on the stability of the ratio of the raw yields N j raw D / N raw D , where the same raw yield extraction method was used in the multiplicity interval j and for the multiplicity-integrated result. The magnitude of this uncertainty depends on p T and meson species. The contribution of the yield extraction procedure to the systematic uncertainties varied between 4-10%.
The influence of D-meson selections, due to the PID and the topological selections, were examined and found to have no significant effect on the final result, since they enter equally into the numerator and denominator of eq. (6.1).
As mentioned in section 4, the contribution of feed-down from B decays to the raw yield was estimated based on FONLL calculations [18]. In this case, it was assumed that the fraction of D mesons that are not from feed-down decays, f prompt , remains constant as a function of multiplicity, causing it to cancel out in the numerator and denominator of the ratio in eq. (6.1). The feed-down contribution was therefore not explicitly subtracted from the final result. A systematic uncertainty related to this hypothesis was assigned by assuming that the fraction f j B / f B , where f B = 1 − f prompt , increases linearly from 1/2 to 2 from the lowest to the highest multiplicity intervals. The resulting uncertainty depends on multiplicity, p T and meson species, and ranges from +4 −0 % to +10 −0 % at low multiplicity and from +0 −4 % to +0 −20 % at high multiplicity. In the analyses as a function of N tracklets , the relative average values of N j tracklets / N tracklets for each interval were corrected to give relative (dN ch /dη) j / dN ch /dη values, as described in section 3.2. The systematic uncertainty due to this correction was estimated in the simulations based on the resolution and the linearity of the correlation between the number of tracklets, N tracklets , and the number of generated charged primary particles, N ch . The deviation from linearity was found to contribute by roughly 5% to the uncertainty on the relative multiplicity. Finally, the uncertainty on the measured dN ch /dη in inelastic p-Pb collisions measured in [68] was considered. This contributed an uncertainty of approximately 4%. The total systematic uncertainty on the relative charged-particle density per N tracklets interval was found to be 6.3%.
In the analyses as a function of N V0A , the measurements are reported as a function of the relative multiplicity N V0A N V0A . The uncertainty on the mean multiplicity values, N V0A , was determined by comparing the mean and median values of the distributions. It was found to be below 5% for each multiplicity interval, and about 30% for the multiplicityintegrated value.

Results
The relative D-meson yields were calculated for each p T and multiplicity interval according to eq. (6.1). The results are reported as a function of the relative charged-particle multiplicity at both backward and central rapidity. It is worth noting that the smaller number of reconstructed D mesons in the lowest and highest p T intervals 1 limited the number of multiplicity intervals of the measurement for those p T intervals. 1 The number of reconstructed D mesons in the lowest and highest pT intervals is smaller than in the other pT intervals. At low pT, the strategy employed to cope with the low signal-to-background ratio was (b) D mesons with 4 < pT < 8 GeV/c. Figure 6. Relative D 0 , D + and D * + meson yields for two selected p T intervals as a function of charged-particle multiplicity at central rapidity. The relative yields are presented in the top panels with their statistical (vertical bars) and systematic (empty boxes) uncertainties, apart from the feed-down fraction uncertainty, which is drawn separately in the bottom panels. The position of the points on the abscissa is the average value of (dN ch /dη) dN ch /dη . For D + and D * + mesons the points are shifted horizontally by 1.5% to improve the visibility. The diagonal (dashed) line is also shown to guide the eye.
The relative D 0 , D + and D * + yields were measured in five p T intervals from 1 to 24 GeV/c as a function of the charged-particle multiplicity at mid-rapidity. Figure 6 presents the measurements for selected p T intervals with their statistical (vertical bars) and systematic (boxes) uncertainties, apart from the feed-down fraction uncertainty, which is drawn separately in the bottom panels. The position of points on the abscissa is the average value of (dN ch /dη) dN ch /dη , but for some meson species they are shifted horizontally by 1.5% to improve the visibility. The relative yields of the three D-meson species are consistent with one another in all p T intervals within uncertainties.
The average of the relative D 0 , D + and D * + yields was evaluated considering the inverse square of their relative statistical uncertainties as weights. The yield extraction uncertainties were treated as uncorrelated systematic uncertainties, while the feed-down subtraction uncertainties were considered as correlated uncertainty sources. Figure 7a presents the average D-meson yields for each p T interval. The results are reported in table 7. The p T evolution of the yields was examined using the results in the 2 < p T < 4 GeV/c interval as reference and by computing the ratio between the average relative D-meson yields in the various p T intervals and those in 2 < p T < 4 GeV/c. The results are shown in figure 7b.
to apply tight topological selections, decreasing the selection efficiency and consequently the number of reconstructed D mesons. At high pT, the small number of candidates is the consequence of the steeply falling D-meson pT spectra.     Average of relative D 0 , D + and D * + yields as a function of the relative charged-particle multiplicity at central rapidity. (a) Average of relative D-meson yields in p T intervals. (b) Ratio of the average relative yields in all p T intervals with respect to that of the 2 < p T < 4 GeV/c interval. The results are presented in the top panels with their statistical (vertical bars) and systematic (boxes) uncertainties, apart from the feed-down fraction uncertainty, which is drawn separately in the bottom panels. The position of the points on the abscissa is the average value of (dN ch /dη) dN ch /dη . For some p T intervals the points are shifted horizontally by 1.5% to improve the visibility. The dashed lines are also shown to guide the eye, a diagonal on (a) and a constant on (b).
The yield increase is independent of transverse momentum within the uncertainties of the measurement. The D-meson yields show a faster-than-linear increase with the chargedparticle multiplicity at central rapidity. The yield increase is approximately a factor of 7 for multiplicities of 4.2 times dN ch /dη . These results are compared with the equivalent measurements in pp collisions, as well as with model calculations, in section 6.2.1.
The measurement of the relative D 0 , D + and D * + yields was also performed as a function of the relative charged-particle multiplicity at large rapidity in the Pb-going direction, thus introducing an η gap between the regions where the D mesons and the multiplicity are measured. The charge collected by the V0A detector, N V0A , was considered as a multiplicity estimator (see section 3.2). Simulations have shown that the collected charge is proportional to the charged-particle multiplicity in the measured η range, 2.8 < η < 5.1. The relative D-meson yields measured in p T and N V0A intervals are reported as a function of the relative multiplicity in the V0A detector, N V0A N V0A . The D 0 , D + and D * + yields are consistent with one another in all the measurement intervals, within uncertainties. The average D-meson yield was calculated with the same procedure used for the results as a function of charged-particle multiplicity at mid-rapidity.    multiplicity at backward rapidity. The yield increase is consistent with a linear growth as a function of multiplicity. The results as a function of V0A multiplicity indicate that the per-event D-meson yield increases as a function of multiplicity, regardless of the η range in which the multiplicity is measured. This remains the case even when the charged-particle yield is measured in a different η interval from the D mesons, which originate from the fragmentation of charm quarks produced in hard partonic scattering processes.
One notable effect to consider when comparing the trends of D-meson production as a function of multiplicity at central and large rapidity is that the charged-particle multiplicity was observed to scale differently with the number of nucleons involved in the p-A interaction depending on η [66,75]. In particular, at central rapidity the charged-particle multiplicity is found to scale with the number of participant nucleons, N part , while at large rapidities in the Pb-going direction (i.e. in the V0A acceptance) it scales with the number of participants of the Pb nucleus, which is equal to N part − 1 = N coll in p-Pb collisions.
It was verified that the results of the D-meson yields as a function of multiplicity are consistent with those of the Q pPb analysis (see section 5). In the Q pPb analysis, D-meson production is studied by dividing the events into centrality classes equally populated by 20% of the events, whereas in this section we examine events with extremely high multiplicity (see tables 2 and 3). Events with low (high) multiplicity correspond to interactions with a smaller (larger) number of hard scatterings per nucleon-nucleon collision, as well as to -22 -JHEP08(2016)078 negative (positive) multiplicity fluctuations which affect event classification and influence both measurements.

Comparison of p-Pb data with pp results and models
The relative D-meson yield (average of D 0 , D + and D * + ) as a function of charged-particle multiplicity at central rapidity in p-Pb collisions at √ s NN = 5.02 TeV is compared with the corresponding pp measurements at √ s = 7 TeV for 2 < p T < 4 GeV/c in figure 9a. A similar relative increase of charmed-meson yield with charged-particle multiplicity is observed in pp and p-Pb collisions. Note that the multiplicity is measured for both pp and p-Pb collisions in the same pseudorapidity range in the laboratory system, which corresponds to different ranges in the centre-of-mass frame for the two collision systems, due to the asymmetry of the beam energies in the p-Pb case.
The increasing yield in pp data can be described by calculations taking into account the contribution of Multiple-Parton Interactions (MPI) [23][24][25], by the influence of the interactions between colour sources in the percolation model [33,34], or by the effect of the initial conditions of the collision followed by a hydrodynamic evolution computed with the EPOS 3 event generator [35,36] where the individual scatterings are identified with parton ladders. In p-Pb collisions, the multiplicity dependence of heavy-flavour production is also affected by the presence of multiple binary nucleon-nucleon interactions, and the initial conditions of the collision are modified due to CNM effects.
Charmed-meson yields in pp and p-Pb collisions as a function of the relative multiplicity at large rapidity are compared in figure 9b for 2 < p T < 4 GeV/c. The multiplicity in p-Pb collisions is measured in 2.8 < η < 5.1 in the Pb-going direction, whereas in pp data the multiplicities at backward (2.8 < η < 5.1) and forward (−3.7 < η < −1.7) pseudorapidity were summed together. The D-meson yields increase faster in pp than p-Pb collisions as a function of the relative multiplicity at backward rapidity. The different pseudorapidity intervals of the multiplicity measurement may contribute to this observation. In addition, measurements in p-Pb collisions differ from those in pp interactions because the initial conditions of the collision are affected by the presence of the Pb nucleus, and because there are multiple binary nucleon-nucleon interactions per p-Pb collision. Figures 10 and 11 present comparisons of the D-meson results and EPOS 3.116 model estimates. The EPOS 3 event generator [35,36] imposes the same theoretical framework for various colliding systems: pp, p-A and A-A. The initial conditions are generated using the "Parton-based Gribov-Regge" formalism [35] of multiple scatterings. Each individual scattering is identified with a parton ladder, composed of a pQCD hard process with initialand final-state radiation. The non-linear effects of parton evolution are treated introducing a saturation scale below which those effects become important. With these initial conditions, a 3D+1 viscous hydrodynamical evolution is applied to the core of the collision [36]. The measurements agree with the EPOS 3 model calculations within uncertainties. The results at high multiplicity are better reproduced by the calculation including a viscous hydrodynamical evolution of the collision, which predicts a faster-than-linear increase of the charmed-meson yield with multiplicity at central rapidity. The same calculation evaluates an approximately linear increase of the charmed-meson yield with the multiplicity   . Average relative D-meson yields in |y lab | < 0.5 as a function of (a) the relative chargedparticle multiplicity at mid-rapidity |η| < 1.0, and (b) at backward-rapidity 2.8 < η < 5.1 (including also −3.7 < η < −1.7 in pp data) for 2 < p T < 4 GeV/c. The relative yields are presented in the top panels with their statistical (vertical bars) and systematic (boxes) uncertainties, apart from the uncertainty on the B feed-down fraction, which is drawn separately in the bottom panels. The positions of the points on the abscissa are the average values of (dN ch /dη) dN ch /dη or N V0A N V0A . A diagonal (dashed) line is also shown to guide the eye. measured at backward rapidity due to the reduced influence of flow on charged particles produced at large rapidity.

Summary
The production of D 0 , D + and D * + mesons as a function of multiplicity in p-Pb collisions at √ s NN = 5.02 TeV, measured with the ALICE detector, has been reported. D mesons were reconstructed in their hadronic decays in different transverse momentum intervals within 1 < p T < 24 GeV/c, in the centre-of-mass rapidity range −0.96 < y cms < 0.04. The multiplicity dependence of D-meson production was studied both by comparing their yields in p-Pb collisions for various centrality classes with those of binary scaled pp collisions at the same centre-of-mass energy via the nuclear modification factor, and by evaluating the relative yields sliced in multiplicity intervals with respect to the multiplicity-integrated ones.
The p T -differential nuclear modification factor, Q pPb , of the D mesons was evaluated with three centrality estimators according to the multiplicity measured in different pseudorapidity intervals: CL1 in |η| < 1.4, V0A in 2.8 < η < 5.1 in the Pb-going direction, and the energy of slow neutrons detected by the ZNA calorimeter at very large rapidity. For each estimator, the events were classified in four classes corresponding to percentiles of the cross section: 0-20%, 20-40%, 40 [35,36] are also shown. The coloured lines represent the calculation curves, whereas the shaded bands represent their statistical uncertainties at given values of (dN ch /dη) dN ch /dη . A diagonal (dashed) line is also shown to guide the eye. three D-meson species fluctuate around unity and are consistent in the measured p T and centrality intervals within uncertainties. The results with the CL1 estimator suggest an ordering from higher (> 1) to lower (< 1) Q pPb values from the 0-20% to the 60-100% centrality class. This disparity is reduced when Q pPb is calculated using the V0A estimator, and vanishes when it is determined with the ZNA estimator (Q pPb ≈ 1). These effects are understood to be due to the biases in the centrality determination in p-Pb collisions based on measurements of multiplicity. The ZNA estimator is the least affected by these sources of biases, and the Q pPb results obtained with this estimator indicate that there is no evidence of a centrality dependence of the D-meson production in p-Pb collisions with respect to that of pp collisions at the same centre-of-mass energy in the measured p T interval within the uncertainties.    [35,36] are also shown. The coloured lines represent the calculation curves, whereas the shaded bands represent their statistical uncertainties at given values of N V0A N V0A . A diagonal (dashed) line is also shown to guide the eye.
The D-meson yields were also studied in p-Pb collisions as a function of the relative charged-particle multiplicity at mid-rapidity, |η| < 1.0, and at large rapidity, 2.8 < η < 5.1, in the Pb-going direction. The relative yields, i.e. the yields in a given multiplicity interval divided by the multiplicity-integrated ones, were calculated differentially in transverse momentum. In contrast to Q pPb , which examines particle production in samples of 20% of the analysed events, this observable explores events from low to extremely high multiplicities corresponding to only 5% (1%) of the analysed events in p-Pb (pp) collisions. The measurements of the relative yields for D 0 , D + and D * + mesons are consistent within the uncertainties. The D-meson yields increase with charged-particle multiplicity, and the increase is independent of p T within the measurement uncertainties. The yield increases with a faster-than-linear trend as a function of the charged-particle multiplicity at mid-rapidity. This behaviour is similar to that of the corresponding measurements in pp collisions at -26 -JHEP08(2016)078 √ s = 7 TeV. Possible interpretations include short-distance gluon radiation, contributions from Multiple-Parton Interactions, the influence of initial conditions followed by a hydrodynamic expansion (EPOS 3 event generator), or the percolation model scenario. In addition, the contribution from multiple binary nucleon-nucleon collisions must be considered in p-Pb collisions. By contrast, the increase of the charmed-meson yields as a function of charged-particle multiplicity at large rapidity in the Pb-going direction is consistent with a linear growth as a function of multiplicity. EPOS   -31 -

JHEP08(2016)078
Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY 4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited.