Measurement of prompt open-charm production cross sections in proton-proton collisions at $\sqrt{s} = $ 13 TeV

The production cross sections for prompt open-charm mesons in proton-proton collisions at a center-of-mass energy of 13 TeV are reported. The measurement is performed using a data sample collected by the CMS experiment corresponding to an integrated luminosity of 29 nb$^{-1}$. The differential production cross sections of the D$^{*\pm}$, D$^\pm$, and D$^0$ ($\overline{\mathrm{D}}^{0}$) mesons are presented in ranges of transverse momentum and pseudorapidity 4 $\lt$ $p_\mathrm{T}$ $\lt$ 100 GeV and $\lvert\eta\rvert$ $\lt$ 2.1, respectively. The results are compared to several theoretical calculations and to previous measurements.


Introduction
Measurements of the production cross sections for open-charm mesons, i.e., mesons containing a single charm quark, in hadronic collisions at LHC center-of-mass energies provide an important test of the theory of strong interactions known as quantum chromodynamics (QCD). Due in part to the presence of several competing scales (charm mass, transverse momentum) that are close to the threshold for the validity of the perturbative expansion, the theoretical uncertainties are rather large. Therefore, experimental constraints on heavy-quark production cross sections are relevant for physics phenomena directly involving heavy quarks or for which they constitute a background. An improved understanding of charm production is furthermore relevant for the exploration of physics processes arising from charm hadron decays, such as neutrinos and other feebly interacting particles [1][2][3], and as a reference for the study of the properties of QCD extreme media [4]. Several studies have been carried out at the LHC. These include measurements at √ s = 7 TeV by the ATLAS Collaboration [5], at √ s = 5 [6,7], 7 [8, 9], and 13 TeV [10] by the ALICE experiment, and at √ s = 5 [11], 7 [12], and 13 TeV [13] by the LHCb Collaboration. The CMS experiment has produced results on open-charm production, analyzing heavy-ion and proton-proton (pp) collisions at √ s NN = 5.02 TeV [4]. In this paper, a study of charm meson production by CMS in pp collisions at √ s = 13 TeV is presented. The analysis is focused on the measurement of cross sections for the prompt production of D * + , D 0 , and D + mesons. The charm mesons are identified via their exclusive decays: where X corresponds to any set of possible particles. Charge conjugation is implied throughout the paper, unless specified otherwise. The designation π + s indicates that, because of the specific kinematics of the D * + decay, this "slow" pion has significantly lower momentum than the kaon and pion decay products of the D 0 . The D * + , D 0 , and D + mesons are reconstructed in the range of transverse momentum 4 < p T < 100 GeV and pseudorapidity |η| < 2.1.
Charm mesons arising from the pp collision point, either directly or as decay products of excited charm resonances (e.g., D 0 coming from D * + decay), are referred to as promptly produced. Contrarily, charm mesons coming from the decays of b hadrons, e.g., B → DX, where X denotes any additional particles, are referred to as nonprompt and are considered as a background process.
Since the kinematic regions covered by the existing measurements differ from the one presented in this paper and are fully complementary in the case of LHCb, the new measurements described in this paper and their comparison to theoretical predictions provide an important contribution to a deeper understanding of the charm production mechanism.
The single-differential cross sections for prompt charm meson production are measured as a function of the transverse momentum p T and the absolute value of pseudorapidity |η|. In principle, rapidity (y) is the proper kinematic variable to study cross sections of massive particles. Here, pseudorapidity is used instead of rapidity to facilitate the comparison with the ATLAS measurement [5], which is the closest in kinematic range to the results presented in this paper. In the kinematic phase space of this measurement, the maximum difference between η and y is less than 5%. Tabulated results are provided in HEPData [14]. The measured cross sections are compared to the predictions from the Monte Carlo (MC) event generators PYTHIA 6.4 [15] and 8. 1 [16], to the theoretical calculations from FONLL [17,18], and to the previous LHC results [4-6, 8, 9, 11-13].

The CMS detector
The central feature of the CMS apparatus is a superconducting solenoid of 6 m internal diameter, providing a magnetic field of 3.8 T. Within the solenoid volume are a silicon pixel and strip tracker, a lead tungstate crystal electromagnetic calorimeter, and a brass and scintillator hadron calorimeter, each composed of a barrel and two endcap sections. Forward calorimeters extend the pseudorapidity coverage provided by the barrel and endcap detectors. Muons are detected in gas-ionization chambers embedded in the steel flux-return yoke outside the solenoid. The most relevant subdetector for this analysis is the silicon tracker. The silicon tracker measures charged particles within the pseudorapidity range |η| < 2.5. It consists of 1440 silicon pixel and 15 148 silicon strip detector modules. For nonisolated particles of 1 < p T < 10 GeV and |η| < 1.4, the track resolutions are typically 1.5% in p T and 25-90  µm in the transverse (longitudinal) impact parameters [19].
Events of interest are selected using a two-tiered trigger system [20]. The first level, composed of custom hardware processors, uses information from the calorimeters and muon detectors to select events at a rate of around 100 kHz within a time interval of less than 4 µs. The second level, known as the high-level trigger, consists of a farm of processors running a version of the full event reconstruction software optimized for fast processing, and reduces the event rate to around 1 kHz before data storage. A more detailed description of the CMS detector, together with a definition of the coordinate system used and the relevant kinematic variables, can be found in Ref. [21].

Data samples, simulation, and event selection
The data sample was acquired in 2016 and corresponds to an integrated luminosity of 29 nb −1 , out of the 36.8 fb −1 collected in that year [22]. The average number of simultaneous pp collisions in the same or nearby bunch crossings, referred to as pileup (PU), for the subset of data used for this analysis is 14. To avoid any bias from a trigger requirement that would aim to select events and require an efficiency correction, the data used in this analysis were collected with an unbiased trigger that only required the presence of crossing beams. This trigger was heavily prescaled, which explains the low effective integrated luminosity for this analysis with respect to the total one collected in 2016.
The candidate vertex with the largest value of summed physics object p 2 T is taken to be the primary pp interaction vertex (PV). The physics objects are the jets, clustered using the jetfinding algorithm [23,24] with the tracks assigned to candidate vertices as inputs, and the associated missing transverse momentum, taken as the negative vector sum of the p T of those jets.
The effects arising from the detector acceptance, reconstruction efficiency, and selection efficiency, whose combination is referred to as the total reconstruction efficiency, are determined from simulated events. These events are generated with PYTHIA 6.4 [15], the heavy-flavor hadrons are decayed with EVTGEN 1.3.0 [25], and the final-state particles are propagated through a simulation of the CMS detector based on GEANT4 v10.00.p02 [26]. The simulated events used to determine the efficiency demand a D * + meson with p T > 3.9 GeV, which is required to decay via D * + → D 0 π + s → K − π + π + s . The p T threshold does not bias the p T spectrum for the D * + and D 0 events as this measurement is only for D mesons with p T > 4 GeV. This sample is also used to determine the D + efficiency, where the D + originates from the hadronization of the other charm quark in each event. In this case, the p T spectrum is biased by the p T threshold on the D * + , and therefore, the simulated D + p T distribution is reweighted to match the D * + spectrum. The effects of PU have been included by overlaying generated minimum-bias events with the main simulated collision. The distribution of the number of PU events is reweighted in the simulation to match the observed data distribution, separately for the 7 main data-taking periods. Following these corrections, the distributions of the kinematic and selection criteria variables obtained from simulation are found to agree with the data for all three mesons.

Charm meson reconstruction
The first step in the reconstruction of the charm mesons is the selection of tracks corresponding to the final-state objects. The criteria used to select the tracks include a minimum p T , a maximum χ 2 of the track fit divided by the number of degrees of freedom (dof), a minimum number of hits in the pixel detector and in the full tracker (pixel and strip detectors), and maximum impact parameters with respect to the PV in the transverse plane (IP xy ) and longitudinal direction (IP z ). The track requirements are summarized in Table 1.
The D 0 (D + ) mesons are reconstructed by combining two (three) tracks with total charge 0 (1) and having an invariant mass M cand within 100 MeV of the nominal meson mass M PDG [27]. When the D 0 candidates are reconstructed in the D * + decay chain, a mass window of 23 MeV on the D 0 mass is used, instead. In the p T range relevant for this analysis, charged pions and kaons cannot be identified efficiently in the CMS detector. A kaon or pion mass hypothesis is thus assumed for the tracks, according to the charge and the specific decay channel. Three topological requirements, whose values are given in Table 1, are also used to reduce the background, which is primarily from random combinations of tracks. First, the secondary vertex (SV) fit χ 2 probability from fitting the two (D 0 ) or three (D + ) tracks to a common vertex is used to ensure the tracks originate from a common point. Second, the cosine of the angle (cos α) between the charm candidate momentum and the vector pointing from the PV to the SV is used to ensure the D meson is consistent with originating from the PV, which reduces background from b hadron decays as well as from random track combinations. Third, the SV significance is the distance between the PV and SV divided by its uncertainty. This is a crucial requirement in the analysis, which provides a considerable reduction in the combinatorial background that is dominated by tracks originating from the PV.
To complete the D * + meson reconstruction, a third track, corresponding to the slow pion, has a pion mass assigned and is kinematically combined with the D 0 candidate. Looser requirements on the p T , χ 2 /dof, and total number of hits are used for this track, as detailed in Table  1. In addition, the impact parameter requirements are changed to require that IP xy and IP z be less than three times their respective uncertainties. To improve the mass resolution, the mass difference ∆M = m(Kππ + s ) − m(Kπ) is used in the analysis instead of the invariant mass of the three-track combination.
For each event, we require there be no more than one candidate for each of the six decay modes (three mesons and two charge-conjugate states). For events in which there is more than one candidate in a particular decay mode, the candidate whose invariant mass is closest to M PDG (D) [27] is chosen for D 0 and D + candidates and the smallest ∆M for the D * + candidates. This arbitration is required for 2, 3, and 11% of the events with a D * + , D 0 , and D + candidate, respectively. By comparing with a random arbitration scheme, it was verified that this method does not introduce a statistically significant bias in terms of the signal yield or signal invariant mass distribution. However, a bias in the background shape has been identified for the D + mesons, which is described in Section 5.

Signal yield determination
The prompt charm meson differential cross section dσ/dp T is measured in 9 bins of p T between 4 and 100 GeV in the range |η| < 2.1; the differential cross section dσ/d|η| is measured in 10 bins of |η|, for |η| < 2.1 and 4 < p T < 100 GeV.
The signal yields, including both prompt and nonprompt decays, are determined using unbinned maximum-likelihood fits to the invariant mass distributions for the various decay modes (the ∆M distribution is used for the D * + ) in each p T and |η| bin. The signal components are modeled by the sum of two Gaussian functions to account for the nonuniform resolution over the detector acceptance. The two Gaussian function means are constrained to be the same. The mean, widths, and normalizations are treated as free parameters. An additional Gaussian function is used to describe the invariant mass shape of D 0 candidates with incorrect pion and kaon mass assignments. The width of this wide Gaussian is taken from simulation bin by bin. The normalization of the wide Gaussian function contribution is fixed to be the same as that of the sum of the two narrow signal Gaussian functions in each bin, reflecting the fact that the number of correct and swapped K/π D 0 signal candidates is the same by construction.
The combinatorial background is described with different functions, according to the decay mode. For the D * + meson, the background is described by a phenomenological threshold function [28] given by where M 0 is the endpoint, taken to be the pion mass, and p 0,1,2 are free parameters. For the D 0 and D + mesons, the combinatorial background component is modeled by a third-degree polynomial function, where the four parameters, including the normalization, are unconstrained in the fits.
In Fig. 1, the fitted invariant mass distributions are reported for two example p T bins, low (5-6GeV) and high (16-24GeV); while Fig. 2 shows the fitted invariant mass distributions for the bins |η| < 0.2 and 1.6 < |η| < 1.8. As expected, at low p T and high |η|, the track momentum and position resolutions are worse, which affect the reconstructed mass width and the distribution shapes, resulting in an increase in the combinatorial background under the peak. Despite the different kinematic regions, it was found that the same functions with different parameter values reproduce all the distributions well.
, and M(K − π + π + ) (lower); charge conjugation is implied. Plots in the left column show the 5 < p T < 6 GeV bin, while the 16 < p T < 24 GeV bin is shown in the right column. The vertical bars on the points represent the statistical uncertainties in the data. The overall result from the fit is shown by the solid line, the fit to the combinatorial background by the dotted line, and, in the middle plots, the fit to the K/π swapped candidates by the red dot-dashed line.
The signal yields and the statistical uncertainties returned by the fits are reported in Tables 2  and 3 for the p T and |η| bins, respectively.

Efficiency estimation
The efficiency is estimated using the signal MC sample and is defined as the fraction of charm signal decays, generated in the kinematic region 4 < p T < 100 GeV and |η| < 2.1, that is reconstructed and survives the selection criteria described in Section 4.1. The efficiency is thus determined for each p T and |η| bin and for both the charge-conjugate states. Taking the D * + channel as an example, these values range from 0.6% for 4 < p T < 5 GeV to 30% for 40 , and M(K − π + π + ) (lower); charge conjugation is implied. Plots in the left column show the |η| < 0.2 bin, while the 1.6 < |η| < 1.8 bin is shown in the right column. The vertical bars on the points represent the statistical uncertainties in the data. The overall result from the fit is shown by the solid line, and the fit to the combinatorial background by the dotted line, and, in the middle plots, the fit to the K/π swapped candidates by the red dot-dashed line.

Contamination from nonprompt decays
The aim of this work is to measure the prompt open-charm production cross sections. Thus, it is important to evaluate and subtract the contribution coming from nonprompt charm mesons arising from b hadron decays. Since consistency with the PV is already part of the selection requirements, the prompt signals and secondary-decay components have similar kinematic variable distributions. The nonprompt-background fraction is estimated from simulation, using minimum-bias events generated with the PYTHIA 8 tune CUETP8M1 [29]. From the generatorlevel information, two subsamples are identified as being representative of the prompt and nonprompt charm meson contributions. The same reconstruction strategy as the one described in Section 4.1 is applied to each of them, and the yields are computed following the method for   yield evaluation reported in Section 4.2, and are labeled N prompt and N nonprompt , respectively. The contamination is then evaluated as the ratio of N nonprompt to the sum (N prompt + N nonprompt ) for each p T and |η| bin. The contamination is nonnegligible, ranging from 5 to 17%, depending on the p T and |η| bin and on the reconstructed meson. This is expected because the requirement on the decay length significance to reject combinatorial backgrounds tends to enhance the contribution from long-lived hadrons. In Fig. 3, the nonprompt-background fractions for the three mesons are shown as a function of p T (left) and |η| (right).
The fraction of nonprompt events in the signal samples is thus taken from simulation, after ensuring a good description of data events in all the relevant quantities. The fractions are found to be consistent using different generator settings and with CMS measurements at 5 TeV [4,30], extrapolated to the final states and kinematic region of this measurement. The nonprompt contribution obtained with this method is subtracted from the measured values of the visible event rate for each p T and |η| bin to get the final prompt results.

Systematic uncertainties
Several systematic uncertainty sources are considered in the measurement of the charm meson cross sections. The dominant effects come from the uncertainties related to the tracking efficiency and the modeling of the invariant mass distributions used in the fit for both the signal and background components. The uncertainties considered in this analysis can be organized into three different categories: decay mode and kinematic bin dependent, only decay mode dependent, and independent of both decay mode and bin.
The first category includes the uncertainty in the efficiency coming from the finite number of MC simulation events, resulting in systematic uncertainties of 0.3, 0.3, and 3.5%, respectively, for the D * + , D 0 , and D + meson cross sections. The last uncertainty is larger than the other two since the sample was enriched in D * + → D 0 π + s → K − π + π + s decays. The uncertainties in the nonprompt event contamination have been treated similarly, considering the limited number of events in the MC simulation. The uncertainties in the CMS nonprompt charm measurement [30] and its extrapolation are propagated as an additional uncertainty of 2%, applied to the three mesons. The resulting systematic uncertainties are 3.5, 2.2, and 2.4% for the D * + , D 0 , and D + cross sections, respectively.
The mass arbitration requirement can produce a peaking background for the small number of events affected by the arbitration, which contributes to the first uncertainty category. This was studied in simulation by selecting a control region that does not contain events with charm mesons and was found to have a nonnegligible contribution only for the D + meson. The effect for the D + was evaluated in each bin and found to contribute a systematic uncertainty of 8% for p T < 12 GeV, but to be negligible at higher p T values. This contribution is independent of |η|, and an uncertainty of 6% is assigned for all |η| bins. This is considered as part of the background modeling systematic uncertainty.
Another source of systematic uncertainty that is included in the first category comes from the p T selection criterion applied to the π + s in the D * + decay chain. The p T spectrum of the slow pion peaks below 0.5 GeV and the selection requirement of p T > 0.3 GeV affects the reconstruction efficiency of the D * + in the first p T bin (4-5 GeV). A systematic uncertainty of 9% is assigned for this bin, which reflects the variation of the simulated event efficiency calculation in the two p T bins of the π + s : 0.2 < p T < 0.3 GeV and 0.3 < p T < 0.4 GeV. Since a reweighting is applied to simulated events in order to reproduce the data PU distribution, a systematic uncertainty is evaluated for each final state. The systematic uncertainty is estimated using the weight calculated for each bin. The cross sections are re-evaluated using the weights raised and lowered by their statistical uncertainties. The largest bin-by-bin change with respect to the cross section calculated with the central value of the weight is taken as the corresponding systematic uncertainty. The total effect is 1% for D * + and D 0 , and 2% for D + .
The uncertainties associated with the branching fraction values and the track reconstruction efficiency depend on the decay mode but not on the candidate p T or |η|. The first one is taken from Ref. [27] and has values of 1.1, 0.8, and 1.7% for D * + , D 0 , and D + , respectively. An uncertainty is assigned to the track reconstruction efficiency according to Ref. [31]. A different procedure is needed for the slow pion from the D * + decay. Because of the soft p T spectrum, a lower tracking efficiency is expected. The uncertainty related to the slow pion is computed by comparing the yield in data and MC simulated events when varying the p T and |η| of the slow pion. This results in a systematic uncertainty of 5.2%. Combining the uncertainties for each track, the total systematic uncertainty from the tracking efficiency is 9.4, 4.2, and 6.1% for the D * + , D 0 , and D + meson cross sections, respectively.
The systematic uncertainty due to the modeling of the invariant mass distribution also falls into the second category. As described in Section 4.2, the signal yields are computed by modeling the resonance peaks with a double-Gaussian function in order to take into account the different resolution effects in various kinematic regions. The uncertainty is estimated by using instead a single-Gaussian function, the sum of three Gaussian functions, and a Crystal Ball function [32,33]. The largest deviation with respect to the default choice is then taken as the systematic uncertainty, yielding 3.6, 5.0, and 4.2% for the D * + , D 0 , and D + meson cross sections, respectively. For the combinatorial background description, the systematic uncertainty is evaluated by replacing the baseline function with a fourth-degree polynomial, resulting in 1.2, 4.8, and 5.3% for the D * + , D 0 , and D + meson cross sections, respectively.
The last category, containing uncertainties independent of both decay modes and kinematic variables, includes those due to data-taking conditions and variation in detector performance. The systematic uncertainty in the integrated luminosity for 2016 is 2.5% [34]. During the 2016 run, the CMS tracker suffered some time-dependent inefficiencies that resulted in a nonnegligible change in charm meson yields. The average PU rate also varied during the data taking. Both these effects are taken into account by correcting the different runs for the tracker inefficiency, which was determined from simulation after the PU distribution in the simulated events was reweighted to match the data. The resulting systematic uncertainty in the correction is estimated to be 1.4%.
All the systematic uncertainties, except the one applied only to the first p T bin of the D * + , are summarized in Table 4. For the bin-dependent systematic uncertainties, the value given in Table 4 is an average that is weighted by the number of signal events. The total uncertainty is evaluated as the sum in quadrature of the individual contributions.

Results
The differential cross sections for prompt charm meson production as a function of p T and |η| are determined using the equations: where N i (D → f ) is the number of prompt charm mesons reconstructed in the selected final state f (including the charge-conjugate final state) for each bin i, after subtracting the nonprompt backgrounds, ∆p T and ∆η = 2∆|η| are the bin widths, B(D → f ) is the branching fraction of the reconstructed decay, ε i (D → f ) is the total reconstruction efficiency of the decay chain evaluated using simulated events, and L is the integrated luminosity.
Here the universality of charm fragmentation is implicitly still assumed, although recent results in a different kinematic range seem to demonstrate that universality is violated in a p T -dependent way [36][37][38][39][40]. This might in principle require an evaluation of p T -dependent downward corrections to the predicted D meson yields of order 5-20% in the p T range measured here, while the measurements are not affected. Since a scheme to consistently evaluate and apply such corrections to the FONLL predictions does not yet exist, and the difference is still expected to be subdominant in particular compared to the large QCD scale uncertainties, no uncertainties were assigned for this source here.
• leading-order (LO) plus parton shower (PS) simulations using PYTHIA 6.4 [15] with the tune Z2* [41], which is based on the CMS Z1 tune [42], but adopting the CTEQ6L PDF set instead of the previous CTEQ5L. The tune is the result of optimizing two parameters that refer to the regularization scale: p T⊥0 , for multiple interactions at a reference energy and the power exponent of the energy rescaling used to determine the value of p T⊥0 as it goes to zero at scales different from the reference scale; • LO plus PS simulations by PYTHIA 8.202 [16] with the tunes: • A2, which is an ATLAS minimum-bias tune [43] validated using their kinematic distributions and based on the tune 4C [44], using the MSTW PDF; • Monash [45], which was developed by re-evaluating the constraints imposed by LEP and SLD on hadronization, in particular with regard to heavy-quark fragmentation and strangeness production; it is a PYTHIA 8 tune using the NNPDF2.3 LO PDF. • CUETP8M1 [29], which is a CMS-specific tune and stands for "CMS Underlying Event Tune PYTHIA 8". It is based on Monash (M1), but the two multiple-parton-interaction (MPI), energy-dependent parameters, which are the MPI cutoff value and the exponent of the √ s dependence, are determined by fitting underlying events in CMS data at √ s = 0.9, 1.96, and 7 TeV. In Ref. [46] it was shown that neither the Monash nor the CUETP8M1 tunes describe well the central value of underlying events in data at √ s = 13 TeV. This suggests that the tune CUETP8M1 does not produce enough charged particles at 13 TeV.
The PYTHIA predictions are from samples of various sizes with corresponding varying statistical uncertainties, which are not shown in the figures. The lower panels in Figs. 4-6 display the ratio of the FONLL and PYTHIA predictions to the data, for which the statistical and total uncertainties are shown by the inner and outer bars, respectively.
The agreement with the different predictions is fair in the wide kinematic range analyzed. However, none of the MC event generators models the data well over the entire measurement region. The measurements tend to favor a higher cross section than predicted by FONLL and a lower one than estimated by PYTHIA, although the predictions from the PYTHIA generator are sensitive to the different tunes used. The cross section predictions from the different PYTHIA tunes differ in both normalization and shape, which confirms that the description of the data provided by the models is sensitive to further model improvements. Overall, the best description of the data is given by the upper edge of the FONLL uncertainty band.

Comparison with other experiments
While there are several previous LHC measurements of charm meson production cross sections, none of them cover the same kinematic region at the same center-of-mass energy as the results presented here. However, since the previous measurements can also be compared to FONLL predictions, it is useful to see if a consistent picture emerges. In the comparisons that follow, no scaling is performed on the previous measurements to account for different kinematic regions, center-of-mass energies, or cross section definitions, but the FONLL predictions were adjusted to match the conditions of the various results.
The measurements from the ATLAS experiment [5], although performed with data collected at √ s = 7 TeV, are the closest to the ones presented in this paper in terms of the acceptance and kinematic regime. Figure 7 shows both the ATLAS and CMS results, compared to the respective FONLL predictions at √ s = 7 and 13 TeV for both D * + (left) and D ± (right) mesons. Since the ATLAS measurement includes both prompt and nonprompt charm mesons, the corresponding FONLL predictions include both components as well. The lower two panels in the figure give the ratio of the FONLL predictions to the respective measurements. The results show good agreement in terms of shape, and the comparison between the data and the theoretical predictions is very similar for the two experiments. The central value of the FONLL predictions tends to underestimate the data, but the upper edge of the FONLL band agrees reasonably well. Figure 8 shows the comparison with the ALICE results [8,9] for the D * + , D 0 , and D + cross sections at √ s = 7 TeV in the range 1 < p T < 24 GeV (0 < p T < 36 GeV for the D 0 ) and for the rapidity region |y| < 0.5. Between the two ALICE measurements, the more recent one [9] is chosen for the comparison. It should be noted that the cross section definition by ALICE includes a factor of 1/2 that accounts for the fact that the measured yields include particles and antiparticles, while the cross sections are given for particles only. The same is true for the corresponding FONLL predictions, as well. To provide a relevant comparison, the CMS measurements are given for p T < 24 GeV (p T < 40 GeV for the D 0 ), for a better comparison with the ALICE points. Both sets of results are consistent with the respective FONLL predictions and close to their upper edge, as shown in the lower two panels.
The CMS experiment has also measured the D 0 cross section in pp collisions at √ s = 5.02 TeV for |y| < 1 [4].  Figure 4: Differential cross sections dσ/dp T (upper) and dσ/d|η| (lower) for prompt D * ± meson production. Black markers represent the data and are compared with several MC simulation models and theoretical predictions. The statistical and total uncertainties are shown by the inner and outer vertical lines, respectively. The FONLL band represents the standard uncertainties in the prediction as detailed in the text. The lower panel gives the ratios of the predictions to the data.  Figure 5: Differential cross section dσ/dp T (upper) and dσ/d|η| (lower) for prompt D 0 (D 0 ) meson production. Black markers represent the data and are compared with several MC simulation models and theoretical predictions. The statistical and total uncertainties are shown by the inner and outer vertical lines, respectively. The FONLL band represents the standard uncertainties in the prediction as detailed in the text. The lower panel gives the ratios of the predictions to the data.  Figure 6: Differential cross section dσ/dp T (upper) and dσ/d|η| (lower) for prompt D ± meson production. Black markers represent the data and are compared with several MC simulation models and theoretical predictions. The statistical and total uncertainties are shown by the inner and outer vertical lines, respectively. The FONLL band represents the standard uncertainties in the prediction as detailed in the text. The lower panel gives the ratios of the predictions to the data.  : Differential cross section dσ/dp T for D * ± (left) and D ± (right) meson production, comparing the production from CMS (black circles, prompt, this paper) at √ s = 13 TeV and ATLAS (red squares, prompt + nonprompt) at √ s = 7 TeV [5]. The corresponding predictions from FONLL are shown by the unfilled and filled boxes, respectively. The vertical lines on the points give the total uncertainties in the data, and the horizontal lines show the bin widths. The two lower panels in each plot give the ratios of the FONLL predictions to the CMS and ATLAS data, shown by circles and squares, respectively. the LHCb Collaboration [13]. Since the η acceptances of the CMS and LHCb experiments differ, the two measurements are complementary and the results presented in this paper extend the reconstruction to a rapidity region not covered by LHCb. The two measurements are compared in Fig. 10 for the LHCb rapidity bin closest to the CMS fiducial region, together with the FONLL predictions. The CMS measurements are shown for p T < 16 GeV to allow a better comparison with the LHCb results. Again, both sets of measurements are in reasonable agreement with the FONLL predictions.  Figure 8: Differential cross section dσ/dp T for prompt D * ± (upper left), D 0 +D 0 (upper right) and D ± (lower) meson production with p T < 24 GeV from CMS (black circles, this paper) at √ s = 13 TeV and ALICE [9] (magenta squares) at √ s = 7 TeV and |y| < 0.5. The corresponding predictions from FONLL are shown by the unfilled and filled boxes, respectively. The cross section definition by ALICE includes a factor of 1/2 that accounts for the fact that the measured yields include particles and antiparticles while the cross sections are given for particles only. The same is true for the corresponding FONLL predictions, as well. The vertical lines on the points give the total uncertainties in the data, and the horizontal lines show the bin widths. The two lower panels in each plot give the ratios of the FONLL predictions to the CMS and ALICE data, shown by circles and squares, respectively.   Figure 9: Differential cross section dσ/dp T for the prompt D 0 + D 0 meson production from CMS at √ s = 13 TeV (black circles, this paper) and 5.02 TeV [4] (light blue squares) for |y| < 1. The corresponding FONLL predictions are shown by the unfilled and filled boxes, respectively. The vertical lines on the points give the total uncertainty in the data, and the horizontal lines show the bin widths. The two lower panels give the ratios of the FONLL predictions to the CMS data at √ s = 13 TeV and 5.02 TeV, shown by circles and squares, respectively.   Figure 10: Differential cross section dσ/dp T for prompt D * ± (upper left), D 0 +D 0 (upper right) and D ± (lower) meson production at √ s = 13 TeV with p T < 16 GeV for CMS (black circles, this paper) for |η| < 2.1 and LHCb [13] (blue squares) for 2 < |y| < 2.5. The corresponding FONLL predictions are shown by the unfilled and filled boxes, respectively. To simplify the results representation, the equivalence between dσ/dp T for 2 < |y| < 2.5 and d 2 σ/dp T dy for 2 < y < 2.5, as in the original publication, has been used. The vertical lines on the points give the total uncertainty in the data, and the horizontal lines show the bin widths. The two lower panels in each plot give the ratios of the FONLL predictions to the CMS and LHCb data, shown by circles and squares, respectively.

Summary
The differential cross sections dσ/dp T and dσ/d|η| for prompt charm meson (D * ± , D 0 (D 0 ), and D ± ) production are measured in the transverse momentum range 4 < p T < 100 GeV and pseudorapidity |η| < 2.1, using data collected by the CMS experiment in proton-proton collisions in 2016 at √ s = 13 TeV, corresponding to an integrated luminosity of 29 nb −1 . The charm mesons were identified with signal invariant mass peaks of high statistical significance. The contamination arising from nonprompt D mesons originating from b hadron decays was removed using Monte Carlo event simulations, validated by measurements.
The measured cross section values are compared to predictions from a theoretical calculation and several different Monte Carlo generators. The agreement with the various models can be considered fair, but no single Monte Carlo simulation or theoretical prediction describes the data well over the entire kinematic range. The measurements tend to favor a higher cross section than predicted by the FONLL calculations [17, 18] and lower than estimated by the PYTHIA event generators [15,16]. The cross section predictions from the different PYTHIA tunes differ in both normalization and shape, which confirms that the description of the data provided by the models is sensitive to further model improvements. Overall, the best description is obtained by the upper edge of the FONLL uncertainty band, which could be taken as a reference prediction for background estimations for other processes, over the full kinematic range covered by all the LHC measurements. By confirming this finding in kinematic regions not previously covered, this measurement makes a contribution to the understanding of charm meson production in hadronic collisions, which is still dominated by large uncertainties in the present theoretical models.

Acknowledgments
We congratulate our colleagues in the CERN accelerator departments for the excellent performance of the LHC and thank the technical and administrative staffs at CERN and at other CMS institutes for their contributions to the success of the CMS effort. In addition, we gratefully acknowledge the computing centers and personnel of the Worldwide LHC Computing Grid and other centers for delivering so effectively the computing infrastructure essential to our analyses. Finally, we acknowledge the enduring support for the construction and operation of the LHC, the CMS detector, and the supporting computing infrastructure provided by the following funding agencies: [8] ALICE Collaboration, "Measurement of charm production at central rapidity in proton-proton collisions at √ s = 7 TeV", JHEP 01 (2012) 128, doi:10.1007/JHEP01(2012)128, arXiv:1111.1553.