Measurement of prompt D 0 and D 0 meson azimuthal anisotropy and search for strong electric ﬁelds in PbPb collisions at

The strong Coulomb ﬁeld created in ultrarelativistic heavy ion collisions is expected to produce a rapidity- dependent difference ( (cid:3) v 2 ) in the second Fourier coeﬃcient of the azimuthal distribution (elliptic ﬂow, v 2 ) between D 0 (uc) and D 0 (uc) mesons. Motivated by the search for evidence of this ﬁeld, the CMS detector at the LHC is used to perform the ﬁrst measurement of (cid:3) v 2 . The rapidity-averaged value is found to be (cid:3) (cid:3) v 2 (cid:4) = 0 . 001 ± 0 . 001 (stat) ± 0 . 003 (syst) in PbPb collisions at √ s NN = 5 . 02 TeV. In addition, the inﬂuence of the collision geometry is explored by measuring the D 0 and D 0 mesons v 2 and triangular ﬂow coeﬃcient ( v 3 ) as functions of rapidity, transverse momentum ( p T ), and event centrality (a measure of the overlap of the two Pb nuclei). A clear centrality dependence of prompt D 0 meson v 2 values is observed, while the v 3 is largely independent of centrality. These trends are consistent with expectations of ﬂow driven by the initial-state geometry. ©


Introduction
The observation of a strongly-coupled quark-gluon plasma (QGP), a state of matter composed of deconfined quarks and gluons, was established by experiments investigating ultrarelativistic heavy ion collisions at the BNL RHIC [1][2][3][4] and CERN LHC [5,6]. The azimuthal particle correlations constitute an effective tool to probe the properties of the QGP [1][2][3][4][5][6][7][8][9]. These correlations are parameterized by a Fourier expansion [10][11][12], with the magnitude of the Fourier coefficients, v n , providing information about the initial collision geometry and its fluctuations [12]. The second-(v 2 ) and third-(v 3 ) order Fourier coefficients are referred to as "elliptic" and "triangular" flow harmonics, respectively. Measuring these coefficients for particle species with different quark composition provides additional information about this hot and dense medium [13]. Because of their large mass, charm and bottom quarks are produced earlier in the collisions than the light quarks (up and down) [14,15]. In addition, the charm and bottom quarks have masses many times larger than the typical temperatures in the QGP [16]. These heavy quarks experience the full evolution of the medium until the hadronization phase. As a consequence, the v n of charmed D 0 (uc) and D 0 (uc) mesons (henceforth referred to as D 0 mesons, except where explicitly stated otherwise) are ex-E-mail address: cms -publication -committee -chair @cern .ch.
pected to receive important contributions from medium energy loss and coalescence effects [17,18]. In ultrarelativistic heavy ion collisions, very strong and transient (∼10 −1 fm/c) magnetic and electric fields are expected to be induced by the collision spectators and participants [19]. Such electromagnetic (EM) fields are predicted to produce a difference in the v n harmonics for positively and negatively charged particles [19]. In such a picture, the magnetic field is mainly responsible for splitting the rapidity ( y)-odd directed flow (v 1 ) [19,20]. The electric field is predicted to induce a charge-dependent splitting in the v 2 coefficient and in the average transverse momentum ( p T ) values of the emitted particles [19]. As charm quarks are expected to be created very early in the collision, they have a higher probability of interacting with this strong EM field than the light flavor quarks [20,21].
In this letter, measurements of the v 2 and v 3 coefficients as functions of D 0 meson rapidity, p T , and lead-lead (PbPb) collision centrality are presented. The collision centrality bins are given in percentage ranges of the total inelastic hadronic cross section, with the 0-10% centrality bin corresponding to the 10% of collisions having the largest overlap of the two nuclei. The flow harmonics are measured using the scalar product method [22,23]. In this analysis, the selection of D 0 mesons uses multivariate methods [24] for selecting D 0 candidates and their antiparticles. The contamination from nonprompt D 0 candidates, arising from B meson decay, is considered as a systematic uncertainty. Using the data recorded in 0.58 nb −1 of integrated luminosity, the flow coefficients are measured within the rapidity range | y| < 2, which is twice as large as achieved in previous CMS measurements [25]. The extension of the measurements to this larger rapidity range, together with smaller statistical uncertainties provided by a larger data set, furnish important inputs for a better understanding of the three-dimensional evolution of the QGP formed in heavy ion collisions. Measurements of the v 2 difference between D 0 and D 0 mesons, v 2 , as a function of rapidity are presented as a method to probe possible effects originating from the Coulomb fields.

Experimental apparatus and data sample
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, there are four primary subdetectors including a silicon pixel and strip tracker detector, a lead tungstate crystal electromagnetic calorimeter, and a brass and scintillator hadron calorimeter, each composed of a barrel and two endcap sections. Iron and quartz-fiber Cherenkov hadron forward (HF) calorimeters cover the pseudorapidity range 2.9 < |η| < 5.2. The HF calorimeters are segmented to form 0.175×0.175 ( η× φ) towers. Muons are measured in gas-ionization detectors embedded in the steel flux-return yoke outside the solenoid. The silicon tracker measures charged particles within the range |η| < 2.5. A 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. [26].
The analysis presented in this letter uses approximately 4.27 × 10 9 minimum bias (MB) PbPb collision events collected by the CMS experiment during the 2018 LHC run. The MB events are triggered by requiring signals in both forward and backward sides of the HF calorimeters [27]. Further selections are applied offline to reject events from background processes (beam-gas interactions and nonhadronic collisions), see Ref. [28] for details. Events are required to have at least one interaction vertex, reconstructed based on two tracks or more, and with a distance of less than 15 cm from the center of the nominal interaction point along the beam axis. The primary interaction vertex is defined as the one with the highest track multiplicity in the event. The shapes of the clusters in the pixel detector have to be compatible with those expected from particles produced at the primary vertex location. The PbPb collision events are also required to have at least two calorimeter towers in each HF detector with energy deposits of more than 4 GeV per tower. These criteria select (99 ± 2)% of inelastic hadronic PbPb collisions. The possibility to have values higher than 100% reflects the possible presence of ultra-peripheral (nonhadronic) collisions in the selected event sample.
Events from Monte Carlo (MC) simulations are used to study both prompt and nonprompt D 0 meson processes. The events are generated using an embedding procedure, in which D 0 mesons generated by pythia 8.212 [29] (tune CP5 [30]) are embedded into MB events from hydjet 1.9 [31]. A full simulation of the CMS detector is performed using Geant4 [32]. The prompt D 0 meson MC simulation is employed to define signal selections and measure efficiency corrections, while the nonprompt D 0 meson MC sample is used to estimate systematic uncertainties coming from nonprompt D 0 contamination.

Reconstruction and selection of D 0 mesons
Prompt D 0 mesons are reconstructed from the decay D 0 → π + + K − and D 0 → π − + K + with a branching fraction of (3.94 ± 0.04)%, using selected tracks with p T > 1.0 GeV/c and within the acceptance of |η| < 2.4. Candidates are formed by combining pairs of tracks from oppositely charged particles and requiring an invariant mass (m inv ) within a ±200 MeV/c 2 window of the worldaverage D 0 meson mass of (1864.83 ± 0.05) MeV/c 2 [33]. For each pair of selected tracks, two possible candidates for D 0 and D 0 mesons are considered by assuming one of the tracks has the pion mass, while the other track has the kaon mass, and vice versa. Kinematic vertex fits are performed to reconstruct the secondary vertices of D 0 candidate decays. After the D 0 candidate reconstruction, a selection using a boosted decision tree (BDT) algorithm from the tmva package [24] is employed. For the BDT training, misidentified D 0 candidates in data events, where pion and kaon have the same charge, are used to mimic the combinatorial background. The signal candidates are taken from MC simulations of prompt D 0 mesons and are required to match D 0 particles at the generator level. The variables related to D 0 mesons used to discriminate the signal from the background are: χ 2 probability for the D 0 vertex fit, 3D distance between the secondary and primary vertices and its significance, the decay length significance projected in the xy-plane, and the angle in two and three dimensions between the momentum of the D 0 meson candidate and the line connecting the primary and the secondary vertices (pointing angle). Related to the decay products of the D 0 meson candidate, the variables used are: the uncertainty in p T returned by the track fitting procedure, the significance of the z and the xy distances of closest approach to the primary vertex, and the number of hits in the tracker detector. These variables are chosen by analyzing their BDT ranking (variables more frequently used in the decision tree) and correlation matrix among all variables. Different BDT boost algorithms are tested, choosing the adaptive boost algorithm [24] as default. Overtraining checks are done for all analysis bins by comparing the BDT distributions from training and testing D 0 meson samples. In addition, a BDT cut optimization is performed in bins of centrality, p T , and rapidity, doing a scan in different BDT scores and finding the one resulting in maximal D 0 mesons signal significance for each analysis bin. Compared to a cutoff-based procedure, this BDT selection almost doubles the signal significance for D 0 mesons in 1 < |y| < 2, and increases the signal significance by 30% for D 0 mesons in |y| < 1, for events with collision centrality in the range 0-30%. For the remaining analysis bins a similar performance of BDT and cutoffbased methods is observed.

Analysis technique
The elliptic and triangular flow coefficients of D 0 mesons are extracted using the scalar product (SP) method, similarly to what was done in a previous CMS publication [25]. In this method, the v n coefficients of D 0 candidates (including backgrounds) are measured using v n {SP} ≡ with the Q -vectors expressed as Q n ≡ M j=1 w j e inφ j , where the sum is over the total number (M) of HF towers above a certain energy threshold (with the weights w j taken as the energy deposited in the HF tower at azimuthal angle φ j ), of tracks with p T above a certain threshold (with w j taken as track p T in φ j angle), or of selected D 0 meson candidates (with w j taken equal to 1).
The Q -vectors related to HF and the tracker are measured and corrected for detector irregularities by applying a flattening and a recentering procedure [12,34]. The Q nA and Q nB are defined using the event-plane measurements from the negative (−5 < η < −3, HF−) and the positive (3 < η < 5, HF+) sides of HF, and Q nC is measured using the tracker information in the region of |η| < 0.75, allowing to minimize the correlations among the three regions, with a gap of more than two units of rapidity. The Q D 0 n vector is defined for each D 0 meson candidate. The averages Q nA Q * nB , Q nA Q * nC , and Q nB Q * nC are made considering all selected events, while the average Q D 0 n Q * nA is made considering all D 0 meson candidates in all selected events. To avoid autocorrelations, the terms Q D 0 n Q * nA and Q nA Q * nB use A = HF− (HF+) when the D 0 meson candidate is at positive (negative) pseudorapidity. One goal of this analysis is to measure the difference ( v n ) between D 0 and D 0 meson flow coefficients, v n , as a function of rapidity, to probe effects from EM fields. The difference v n is measured as: The The α(m inv ) parameter, which characterizes the signal fraction as a function of mass, is defined as follows: For extracting the difference v sig n , the following expression is em- while the constant parameter const in Eq. (5) is added to account for possible fluctuations in the background v n component. The relevance of this const parameter was investigated by redoing v n measurements in MC simulation (without azimuthal correlations or effects from EM fields), indicating that this parameter improves the fit quality and does not introduce artificial signals. A crosscheck is performed by redoing the measurements using a linear function instead of a constant. No significant changes in the central values of v 2 and on their uncertainties are observed. Fig. 1 shows an example of a simultaneous fit for v 2 and v 2 .
After performing the fits for extracting the signal v n , there is still a sizable fraction of nonprompt D 0 mesons embedded in v sig n . The extracted v n can be written as The nonprompt D 0 meson contamination is taken into account as a systematic uncertainty, by checking that the nonprompt D 0 meson fraction is always smaller than 12% (i.e., comparable to the uncertainties in the reconstructed D 0 meson yield). This implies that the central values of v n will not be considerably affected by this component, being compatible within statistical uncertainties. Such a low fraction arises from the use of prompt D 0 meson signals in the BDT training, together with variables that are highly correlated with the distance of closest approach (DCA) to the primary vertex, which is defined as the flight distance of the D 0 particle times the sine of the pointing angle in three dimensions. Additional DCA selection and dedicated training, involving prompt and nonprompt D 0 meson signals, do not bring considerable improvements in performance. The prompt and nonprompt D 0 meson fractions are obtained using the DCA variable. For prompt D 0 mesons, the nonzero DCA corresponds to the detector resolution, and is expected to be concentrated around zero. For nonprompt D 0 mesons, larger values of DCA result from the B meson decay. To extract the prompt and nonprompt D 0 meson fractions, a fit to the DCA distributions is performed in data considering DCA shapes from MC simulations for prompt and nonprompt D 0 meson components. The nonprompt D 0 meson v n is estimated by considering two regions in the DCA: one with very low fraction (2.7-8.0%) of nonprompt D 0 particles (DCA < 0.012 cm), and one with a high fraction (62.0-88.0%) of nonprompt D 0 particles (DCA > 0.012 cm). Using this information together with Eq. (6), it is possible to estimate v nonprompt n by solving a system of two equations from the two DCA regions. In the current analysis this procedure can only be done in wide p T , centrality, and rapidity bins, because of the limited amount of data available in the region with DCA > 0.012 cm.

Systematic uncertainties
The sources of systematic uncertainties include the D 0 identification requirements (BDT selection); the probability distribution function (PDF) for modeling the background in the invariant mass fit; the impact of acceptance and efficiency of the D 0 meson yield; the variation of the PDF for modeling the background v n ; and the remaining nonprompt D 0 contamination. With the exception of the last component, the uncertainties are quoted as absolute values of v n and v n after comparing the default analysis configuration with the variations. To diminish the influence of statistical fluctuations, after observing no special trends in the deviations from the default measurements, the systematic uncertainties are evaluated by averaging the deviations with a constant fit as a function of the analysis bins.
In order to take into account the systematic uncertainty associated with the BDT selection, the BDT cut is varied up and down by the maximal deviation between the BDT optimized selection based on MC simulations and data. The BDT cuts (and variations for systematic uncertainties) are defined in bins of collision centrality, p T , and rapidity, ranging from 0.28 to 0.47 (±0.02-0.03). Regarding the effect of the background mass modeling, either an exponential function together with a second order polynomial, or just a second order polynomial, are considered instead of the default fit function using a third-order polynomial. To fit v n as a function of mass, the default configuration using a linear function is replaced by either a constant or a second order polynomial. Although the D 0 meson selection efficiency essentially cancels in v n measurements, a systematic uncertainty is assigned by comparing the results with and without applying corrections based on MC simulations in bins of p T and rapidity. The D 0 meson selection efficiency times acceptance varies from 0.5 to 12.5% in the p T range of 1.  meson v n in wide bins of p T , rapidity, and centrality. A relative systematic uncertainty is obtained by comparing v n from mixed prompt and nonprompt D 0 mesons to the v n derived from nonprompt D 0 mesons. Table 1 summarizes the estimates of systematic uncertainties in absolute values for v 2 , v 3 , and v 2 . The ranges of variation of the uncertainties are presented for each binning.

Results
Results for prompt D 0 meson v 2 and v 3 anisotropic flow coefficients, obtained with 2018 PbPb data, as functions of p T and for |y| < 1, are shown in Fig. 2 for three centrality ranges: 0-10%, 10-30%, and 30-50%. The results extend previously published data from CMS [25], by extending the high-p T coverage to ∼60.0 GeV/c and by providing finer p T bins. These high-precision data are compatible with previous measurements from Ref. [25], and a clear trend of rise and fall from low to high p T is observed for both v 2 and v 3 across the full centrality range. This behavior is similar to that observed for inclusive charged particles [35] for |η| < Calculations from theoretical models at midrapidity (| y| < 1) are also presented. These models use different assumptions of the QGP properties, for example in the thermal evolution of the collision system and in the initial-state conditions before the formation of the QGP. In addition, different mechanisms are considered regarding the interaction of heavy quarks with the medium and for the hadronization process. Results from the models LBT [36], CUJET 3.0 [37], and SUBATECH [38] include collisional and radiative energy losses, while those from the models TAMU [39], PHSD [15], and TAMU SMCs [40] include only collisional energy loss. Initial-state fluctuations are included in the calculations by LBT, SUBATECH, and PHSD, and calculations for the v 3 coefficient are only available from these three models. Coalescence mechanisms are also included in LBT, SUBATECH, TAMU, PHSD, and TAMU SMCs. While most models seem to capture the qualitative trend of the data (except for the v 2 description provided by TAMU in the 10-50% centrality range), most of the models do not provide a quantitative description over the full range, except for TAMU SMCs. The TAMU SMCs version improves the TAMU model by implementing event-by-event space-momentum correlations (SMCs) between charm quarks and the high-flow partons in the QGP medium [40]. Since it does not include initial-state fluctuations, TAMU SMCs does not provide v 2 calculations for centrality values between 0-10%. This puts more stringent constraints on the development of the collective flow for charm quarks in the QGP medium, giving  [15]. The TAMU SMCs model [40] is available only in the 10-50% centrality bins. further inputs for understanding heavy-quark interactions with the medium (for example, energy loss and coalescence mechanisms).
Results for the rapidity dependence of heavy-flavor collective flow are presented for the first time for prompt D 0 meson v 2 and v 3 as functions of p T , both at midrapidity (| y| < 1) and in the forward (1 < |y| < 2) region, as shown in Fig. 3. No clear rapidity dependence is observed for both v 2 and v 3 as functions of p T . This observation is similar to that for inclusive charged-hadron measurements [41].
In Fig. 4 (left), results for prompt D 0 mesons v 2 and v 3 , averaged over 2.0 < p T < 8.0 GeV/c, for |y| < 1 and 1 < |y| < 2, are presented as a function of collision centrality. This p T range is chosen in order to cover the widest possible p T range, while maximizing the D 0 meson signal yield significance. These p T -and rapidity- integrated results include an additional centrality bin (50-70%), which has an insufficient number of events for the full differential analysis. For both mid-and forward-rapidity regions, the v 2 results show a clear increase from the most central to mid-central events, and then a declining trend toward the most peripheral events. This trend is similar to that observed for inclusive charged particles (also shown in Fig. 4), and can be understood in terms of collision geometry and viscosity effects. In particular, a faster increase of v 2 is observed from central to peripheral collisions for charged particles compared to prompt D 0 mesons. This feature was also observed when comparing v 2 of low-p T J/ψ with charged pions [42], where it is claimed that this could be understood in terms of two phenomena: one, associated with transport models predicting an increasing fraction of regenerated J/ψ at low-p T , when going from peripheral to central collisions; the other, not related to regeneration, is associated with a possible partial or later thermalization of charm quarks compared to light quarks [42]. The v 3 shows no centrality dependence, which is also consistent with expectations from collision geometry fluctuations [43]. are not yet available, they are expected to be much larger than those for inclusive charged particles. In the case of v 1 , the ALICE collaboration reported results about three orders of magnitude larger than measurements for charged hadrons [44], although the uncertainties prevent a clear conclusion. The main reason is that heavy-flavor quarks are usually produced much earlier than lightflavor quarks, the former being predominantly produced soon after the collision takes place, when the EM field strength is several or-

Summary
Measurements of the elliptic (v 2 ) and triangular (v 3 ) flow coefficients of prompt D 0 mesons are presented as functions of transverse momentum (p T ), rapidity, and collision centrality, in PbPb collisions at √ s NN = 5.02 TeV. The results improve previously published CMS data by extending the p T and rapidity coverage and by providing more differential information in p T , rapidity, and centrality. A clear centrality dependence of prompt D 0 meson v 2 is observed, while v 3 is largely centrality independent. These trends are consistent with the expectation that v 2 and v 3 are driven by initial-state geometry. A weak rapidity dependence of prompt D 0 meson v 2 and v 3 is observed. When comparing various theoretical calculations to the data at midrapidity, no model is able to describe the data over the full centrality and p T ranges. Motivated by the search for evidence of the strong electric field expected in PbPb collisions, a first measurement of the v 2 flow coefficient difference ( v 2 ) between D 0 and D 0 mesons as a function of rapidity is presented. The rapidity-averaged v 2 difference is measured to be v 2 = 0.001 ± 0.001 (stat) ± 0.003 (syst). This indicates that there is no evidence that charm hadron collective flow is affected by the strong Coulomb field created in ultrarelativistic heavy ion collisions. Future comparisons of theoretical models with these results may provide constraints on the electric conductivity of the quark-gluon plasma.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgements
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 for delivering so effectively the computing infrastructure essential to our analyses.