Prompt and non-prompt $J/\psi$ elliptic flow in Pb+Pb collisions at $\sqrt{s_{_\text {NN}}} = 5.02$ Tev with the ATLAS detector 

Abstract

The elliptic flow of prompt and non-prompt $J/\psi$ was measured in the dimuon decay channel in Pb+Pb collisions at $\sqrt{s_{_\text {NN}}}=5.02$ $\text {TeV}$ with an integrated luminosity of $0.42~\mathrm {nb}^{-1}$ with the ATLAS detector at the LHC. The prompt and non-prompt signals are separated using a two-dimensional simultaneous fit of the invariant mass and pseudo-proper decay time of the dimuon system from the $J/\psi$ decay. The measurement is performed in the kinematic range of dimuon transverse momentum and rapidity $9<p_\mathrm {T}<30$ $\text {GeV}$, $|y|<2$, and 0–60% collision centrality. The elliptic flow coefficient, $v_2$, is evaluated relative to the event plane and the results are presented as a function of transverse momentum, rapidity and centrality. It is found that prompt and non-prompt $J/\psi$ mesons have non-zero elliptic flow. Prompt $J/\psi$ $v_2$ decreases as a function of $p_\mathrm {T}$, while for non-prompt $J/\psi$ it is, with limited statistical significance, consistent with a flat behaviour over the studied kinematic region. There is no observed dependence on rapidity or centrality.

1 Introduction

With the advent of lead–lead collisions at the centre-of-mass energy of $5.02~\text {TeV}$ per nucleon–nucleon pair, new opportunities open up for understanding the detailed properties of the hot dense plasma produced in such collisions. A special advantage of studies with quarkonia as hard probes of the plasma properties is that the comparison of prompt and non-prompt $J/\psi$ mesons elucidates the differences between the responses of the produced c-quark and b-quark systems. This is because the prompt $J/\psi$ mesons are $c\bar{c}$ systems produced soon after the collision whereas the non-prompt $J/\psi$ mesons come from decays of b-hadrons that are formed outside the medium [1]. Thus, the comparison of these two classes of $J/\psi$ mesons probes the flavour dependence of the mechanisms of the interactions of heavy quarks with the medium. ATLAS measurements of the attenuation of both the prompt and non-prompt $J/\psi$ meson yields indicate very strong medium effects that are surprisingly similar in magnitude at this collision energy [2].

A complementary and powerful probe into the heavy-quark flavour dependence of interaction mechanisms can be obtained by studying the azimuthal asymmetries of prompt and non-prompt quarkonia. Such studies [3,4,5] are especially useful in the dimuon transverse momentum range $9<p_{\text {T}} <30~\text {GeV}$, investigated in this paper, since this range represents the transition between the lower $p_{\text {T}}$ region, in which recombination processes are believed to play an important role [6,7,8], and the higher $p_{\text {T}}$ region in which other processes are expected to dominate, such as absorption due to colour-exchange interactions with the medium [9,10,11,12]. The naive expectation in this range is that recombination processes will partially couple the produced $J/\psi$ to the hydrodynamic flow of the hot medium, resulting in an enhancement of the observed azimuthal asymmetry at lower $p_{\text {T}}$ relative to higher $p_{\text {T}}$  [7, 11]. In this picture, a flavour-dependent enhancement of the azimuthal asymmetry of $J/\psi$ at low $p_{\text {T}}$ can be interpreted as a difference in the degree of recombination between c- and b-quarks and it is expected that any flavour dependence will vanish at higher values of $p_{\text {T}}$, which are accessible by this measurement. A recent transport model study suggests a sensitivity of charm quarks to hydrodynamic flow [13]. In this model, additionally, a strong suppression of the prompt $J/\psi$ yield in the final-state medium should lead to an azimuthal asymmetry even in the high $p_{\text {T}}$ region [11].

In non-central collisions, the overlap region of the colliding ions has an elliptic shape. The particle yield is influenced by this matter distribution, leading to the observation of an azimuthal anisotropy relative to the reaction plane as observed for charged hadrons [14,15,16,17,18]. The azimuthal distribution of particles is characterised by a Fourier expansion of the particle yield:

$$\begin{aligned} \frac{\mathrm {d}N}{\mathrm {d}\phi } \propto 1 + \sum _{n=1}^{\infty } 2v_{n}\cos [n(\phi - \Psi _{n})], \end{aligned}$$

where $\phi$ is the azimuthal angle of the particle relative to the detector frame of reference, and $\Psi _{n}$ is the nth harmonic of the event-plane angle, which can be estimated using the event-plane method [19]. The second-order coefficient, $v_2$, is referred to as elliptic flow and its magnitude quantifies the yield modulation relative to the elliptical shape of the initial matter distribution.

Interestingly the observed azimuthal asymmetry for prompt $J/\psi$ is the same in central collisions as in non-central collisions [5], although for inclusive $J/\psi$ some indication of the centrality dependence is reported in the lower $p_{\text {T}}$ region at forward rapidities at 5.02 $\text {TeV}$ [4] and at 2.76 $\text {TeV}$ [3]. This is in contradiction with the expected hydrodynamic behaviour, which is confirmed by the results for charged hadrons where the anisotropies are more significant in semi-central collisions than in peripheral and central collisions [14,15,16,17,18]. This intriguing observation may ultimately provide more insight into the origins of azimuthal asymmetries beyond a simple hydrodynamic picture. Further, there is evidence of a surprising universality among many different probes, such as D-mesons and jets [20, 21], for which the $v_2$ values are very similar at high $p_{\text {T}}$. This paper provides $v_2$ measurements as a function of transverse momentum, rapidity and collision centrality for both prompt and non-prompt $J/\psi$ in the dimuon decay channel, extending the kinematic range covered by recent results from other LHC experiments [3,4,5].

2 ATLAS detector

The ATLAS detector [22] at the LHC covers nearly the entire solid angle around the collision point. It consists of an inner tracking detector (ID) surrounded by a thin superconducting solenoid, electromagnetic and hadronic calorimeters, and a muon spectrometer (MS) incorporating three large superconducting toroid magnets.¹

A high-granularity silicon pixel detector surrounds the interaction region and typically provides four measurements per track. It is followed by a silicon microstrip tracker, which provides around eight two-dimensional measurement points per track. These silicon detectors are complemented by a transition radiation tracker (TRT), which enables radially extended track reconstruction up to $|\eta | = 2.0$. The ID system is immersed in a 2 T axial magnetic field and provides charged-particle tracking in the range $|\eta | < 2.5$.

The calorimeter system covers the pseudorapidity range $|\eta | < 4.9$. Within the region $|\eta |< 3.2$, electromagnetic calorimetry is provided by barrel and endcap high-granularity lead/liquid-argon (LAr) electromagnetic calorimeters, with an additional thin LAr presampler covering $|\eta | < 1.8$ to correct for energy loss in material upstream of the calorimeters. Hadronic calorimetry is provided by a steel/scintillator-tile calorimeter, segmented into three barrel structures within $|\eta | < 1.7$, and two copper/LAr hadronic endcap calorimeters covering $1.5< |\eta | < 3.2$. The high $|\eta |$ region, $3.2< |\eta | < 4.9$, is covered by forward copper/LAr and tungsten/LAr calorimeter (FCal) modules optimised for electromagnetic and hadronic measurements respectively.

The MS comprises separate trigger and high-precision tracking chambers measuring the deflection of muons in a magnetic field generated by the superconducting air-core toroid magnets. The precision chamber system covers the region $|\eta | < 2.7$ with three layers of monitored drift tubes, complemented by cathode strip chambers in the forward region, where the background is the highest. The muon trigger system covers the range $|\eta | < 2.4$ with resistive plate chambers in the barrel, and thin gap chambers in the endcap regions.

In addition to the muon trigger, two triggers are used in Pb+Pb collisions to select minimum-bias events. These are based on the presence of a minimum amount of transverse energy, of at least 50 GeV, in all sections of the calorimeter system with $|\eta | < 3.2$ or, for events which do not meet this condition, on the presence of energy deposits in both zero-degree calorimeters, which are primarily sensitive to spectator neutrons in the region $|\eta | > 8.3$.

A two-level trigger system is used to select events [23]. The first-level trigger is implemented in hardware and uses a subset of detector information to reduce the event rate to a design value of at most 100 kHz. This is followed by a software-based high-level trigger which reduces the event rate to about 1 kHz.

3 Analysis

3.1 Data, event selection and centrality definition

Data from Pb+Pb collisions at $\sqrt{s_{_\text {NN}}} =5.02~\text {TeV}$ were recorded by the ATLAS experiment in 2015. Events were collected using a trigger requiring at least two muons, both with $p_{\text {T}} ^\mu > 4~\text {GeV}$. This muon triggered dataset has an integrated luminosity of 0.42 $\text{ nb }^{-1}$. In the offline analysis, reconstructed muons are required to satisfy the tight muon working point ignoring the TRT requirements [24], have $p_{\text {T}} ^\mu > 4~\text {GeV}$, $|\eta | < 2.4$, and be matched to the muons reconstructed at the trigger level. In addition, muon pairs are required to have pair $p_{\text {T}} > 9~\text {GeV}$, rapidity $|y| < 2$ and be in the invariant mass range $2.6< m_{\mu \mu } < 3.5~\text {GeV}$. In addition to the muon triggered event sample, a minimum-bias triggered event sample and Monte Carlo (MC) simulated event samples were used for studies of the detector performance. Prompt ($pp \rightarrow J/\psi \rightarrow \mu \mu$) and non-prompt ($pp \rightarrow b\bar{b} \rightarrow J/\psi \rightarrow \mu \mu$) samples of $J/\psi$ were produced using Pythia 8.212 [25] for event generation and Photos [26] for electromagnetic radiation corrections. The A14 set of tuned parameters [27] is used together with CTEQ6L1 parton distribution function set [28]. The response of the ATLAS detector was simulated using Geant4 [29, 30]. Simulated events are overlaid on a sample of minimum-bias Pb+Pb events produced with HIJING [31] to replicate the high-multiplicity environment of heavy-ion collisions.

To characterise the Pb+Pb collision geometry, events are classified into centrality intervals determined by the summed transverse energy deposited in the FCal, $\sum E_{\text {T}} ^{\mathrm {FCal}}$, in each event. Centrality intervals are defined according to successive percentiles of the $\sum E_{\text {T}} ^{\mathrm {FCal}}$ distribution ordered from the most central (the highest $\sum E_{\text {T}} ^{\mathrm {FCal}}$, the smallest impact parameter) to the most peripheral collisions (the lowest $\sum E_{\text {T}} ^{\mathrm {FCal}}$, the largest impact parameter). The average number of nucleons participating in the collision, $\langle N_{\text {part}} \rangle$, is calculated using a Glauber model analysis of the $\sum E_{\text {T}} ^{\mathrm {FCal}}$ distribution [32, 33]. The centrality intervals used in this analysis are quoted together with $\langle N_{\text {part}} \rangle$ in Table 1. Only events in the 0–60% centrality interval are used. Events in the 60–80% are disregarded from the study, since they represent only a small fraction of events, below $3\%$ of the full $J/\psi$ sample, which is not significant for the measurement.

3.2 Signal extraction and observable determination

The $J/\psi$ $v_2$, the second-order coefficient of the Fourier decomposition of the azimuthal angle distribution, is measured using the event-plane method [19]. The event-plane angle is estimated by its second-order harmonic, $\Psi _2$, using the distribution of transverse energy deposited in the forward calorimeters. Similar methods are described in detail for previous azimuthal anisotropy analyses of charged hadrons with the ATLAS detector [14, 15, 34]. To reduce autocorrelations in the event-plane analysis, $v_2$ is measured by correlating $J/\psi$ with positive (negative) pseudorapidity with the event-plane angle measured using the FCal in the negative (positive) $\eta$-region. The prompt and non-prompt $J/\psi$ yields are obtained from two-dimensional fits of the reconstructed $J/\psi$ invariant mass and pseudo-proper decay time distributions. The azimuthal distributions of the prompt and non-prompt yields are fitted simultaneously to obtain the elliptic flow coefficients.

The second-order harmonic of the event-plane angle is determined using measurements of transverse energy deposits in each FCal system positioned at $\eta >3.2$ and $\eta <-3.2$. The flow vector is $\pmb {q}_2=\sum _{i} w_i(\cos (2\phi _i)\pmb {\hat{x}}+\sin (2\phi _i)\pmb {\hat{y}})$, where $\phi _{i}$ is the azimuthal coordinate of the ith calorimeter tower, $w_{i}$ is a weight that equals the transverse energy deposited in the calorimeter tower, and the sum is over all the FCal towers. The FCal towers consist of calorimeter cells grouped into regions in $\Delta \eta \times \Delta \phi$ of $0.1\times 0.1$. The flow vector is determined separately in the positive and negative rapidity regions. The event-plane angle is then calculated as $2\Psi _{2} = \tan ^{-1}(\pmb {q}_{2}\cdot \pmb {{\hat{y}}}/\pmb {q}_{2}\cdot \pmb {{\hat{x}}})$.

To ensure the uniformity of the event-plane angle distribution, the raw flow vector is corrected by subtracting its mean value, obtained by averaging over all events in a given centrality interval, so that $\pmb {q}_2=\pmb {q}_2^\text {raw}-\langle \pmb {q}_2^\text {raw}\rangle$. Since the mean values $\langle \pmb {q}_2^\text {raw}\rangle$ are found to be independent of the collision centrality, the $\langle \pmb {q}_2^\text {raw}\rangle$ averaged over all analysed events (centrality 0–60%) is used. The remaining modulations of the event-plane angle distribution are removed by including a shift $\delta \Psi _{2} = \sum _{k=1}^{k_{\text {max}}}(1/k)\left[ A_2 \cos (2k\Psi _2) + B_2 \sin (2k\Psi _2)\right]$ [19]. The calculated coefficients $A_2=-\langle \sin (2k\Psi _2) \rangle$ and $B_2=\langle \cos (2k\Psi _2) \rangle$ are found to be centrality dependent up to $k=2$, and the sum is performed up to a conservative choice of $k_\text {max} = 12$. After these recentring and flattening corrections the event-plane angle distribution follows a uniform distribution.

The event-plane resolution, $\mathcal {R}$, is determined using the two-sub-event method [19] and the minimum-bias event sample. Values of $\Psi _2$ are determined on both sides of the detector and are used to calculate $\mathcal {R} = \sqrt{\langle \cos (2\Delta \Psi _2)\rangle }$, where $\Delta \Psi _2$ is the difference between the values of $\Psi _2$ computed using the FCal modules in the positive and negative $\eta$-region of the detector. The resulting event-plane resolution depends strongly on centrality: it is poorer at low and high $\sum E_{\text {T}} ^{\mathrm {FCal}}$ and better at middle values of $\sum E_{\text {T}} ^{\mathrm {FCal}}$. The resolution is calculated in minimum-bias events in very fine bins of transverse energy. The average value of the resolution in wider bins must account for the different $\sum E_{\text {T}} ^{\mathrm {FCal}}$ distribution of the sample of events containing $J/\psi$ candidates. Thus, the event-plane resolution is weighted by the number of $J/\psi$ candidates in a given centrality interval relative to the number of minimum-bias events in the same interval. The values of $\mathcal {R}$ for the centrality intervals used in this analysis are shown in Table 1.

Table 1 The average number of participating nucleons, $\langle N_{\text {part}} \rangle$, and the event-plane resolution, $\mathcal {R}$, with their total uncertainties in each centrality interval

To account for detector effects, each muon pair is corrected for trigger efficiency, $\epsilon _\text {trig}$, reconstruction efficiency, $\epsilon _\text {reco}$, and detector acceptance, A. These three quantities form a per-dimuon weight:

$$\begin{aligned} w^{-1} = A\times \epsilon _\text {trig}\times \epsilon _\text {reco}. \end{aligned}$$

Trigger and reconstruction efficiencies are studied using the tag-and-probe method in data and in MC simulations as a function of the muon $p_{\text {T}} ^\mu$ and $\eta ^\mu$. The reconstruction efficiency increases from low to high $p_{\text {T}} ^\mu$ and decreases from central to forward pseudorapidity, becoming constant at $p_{\text {T}} ^\mu > 6~\text {GeV}$ with a maximum efficiency of about 90%. Trigger efficiency increases from low to high $p_{\text {T}} ^\mu$ and from central to forward pseudorapidity, increasing from 50% to 85% between the lowest and highest $p_{\text {T}} ^\mu$. The acceptance is studied from MC simulations. It is defined as the probability that the $J/\psi$ decay products fall within the fiducial volume $p_{\text {T}} ^\mu > 4~\text {GeV}$ and $|\eta ^\mu | < 2.4$ and assuming unpolarised $J/\psi$ production [35,36,37]. A detailed description of the performance studies is presented in Ref. [2].

The separation of the prompt and non-prompt $J/\psi$ signals is performed using the pseudo-proper decay time of the $J/\psi$ candidate, $\tau _{\mu \mu } = L_{xy}m_{J/\psi }/p_{\text {T}}$, where $L_{xy}$ is the distance between the position of the dimuon vertex and the primary vertex projected onto the transverse plane, $m_{J/\psi } = 3.096~\text {GeV}$ is the value of the $J/\psi$ mass [38], and $p_{\text {T}}$ is the transverse momentum of the dimuon system.

The corrected two-dimensional distribution of the number of events as a function of pseudo-proper decay time and dimuon invariant mass is used to determine the prompt and non-prompt $J/\psi$ yields. The probability distribution function (PDF) for the fit is defined as a sum of five terms, where each term is the product of functions that depend on the dimuon invariant mass or pseudo-proper decay time. The PDF is written in a compact form as:

$$\begin{aligned} P(m_{\mu \mu },\tau _{\mu \mu }) = \sum ^5_{i=1}N_i f_i(m_{\mu \mu }) \cdot h_i(\tau _{\mu \mu }) \otimes g(\tau _{\mu \mu }), \end{aligned}$$

where $N_i$ is the normalisation factor of each component, $f_i(m_{\mu \mu })$ and $h_i(\tau _{\mu \mu })$ are distribution functions of the invariant mass, $m_{\mu \mu }$, and the pseudo-proper decay time, $\tau _{\mu \mu }$, respectively; $g(\tau _{\mu \mu })$ is the resolution function described by a double Gaussian distribution, and the $\otimes$-symbol denotes a convolution. The PDF terms are defined by Crystal Ball (CB) [39], Gaussian (G), Dirac delta ($\delta$), and exponential (E) distributions as specified in Table 2.

Table 2 Individual components of the probability distribution function in the default fit model used to extract the prompt and non-prompt contribution for $J/\psi$ signal and background. $F_\mathrm {CB}$ and $F_\mathrm {G}$ are the Crystal Ball (CB) and Gaussian (G) distribution functions respectively, $\omega$ is the relative fraction of the CB and G terms, $F_\mathrm {E}$ is an exponential (E) function, and $\delta (\tau _{\mu \mu })$ is the Dirac delta function

The signal invariant mass shapes are described by the sum of a CB function and a single Gaussian function with a common mean. The width term in the CB function is equal to the Gaussian standard deviation times a scaling term, fixed from MC simulation studies. The CB left-tail and height parameters are also fixed from MC studies and variations of the two parameters are considered as part of the fit model’s systematic uncertainties. The relative fraction of the CB and Gaussian functions, $\omega$, is free in the fit. The prompt background contribution to the invariant mass spectrum follows a nearly flat distribution, and is modelled by an exponential function, denoted $F_{\mathrm {E}_2}(m_{\mu \mu })$ in Table 2. The non-prompt contribution to the background requires two exponential functions, denoted $F_{\mathrm {E}_3}(m_{\mu \mu })$ and $F_{\mathrm {E}_5}(m_{\mu \mu })$ in Table 2, respectively.

The pseudo-proper decay time of the prompt signal is modelled with a Dirac delta function, while the non-prompt signal is described by a single-sided exponential, denoted $F_{\mathrm {E}_1}(\tau _{\mu \mu })$ in Table 2. The backgrounds are represented by the sum of one prompt component and two non-prompt components. The prompt background component is described by a Dirac delta function. One of the non-prompt background contributions is described by a single-sided decay model (for positive $\tau _{\mu \mu }$ only), and the other is described by a double-sided decay model, denoted $F_{\mathrm {E}_4}(\tau _{\mu \mu })$ and $F_{\mathrm {E}_6}(|\tau _{\mu \mu }|)$ in Table 2, accounting for candidates of mis-reconstructed or non-coherent muon pairs resulting from other Drell–Yan muons and combinatorial background. A double Gaussian resolution function, $g(\tau _{\mu \mu })$, is used in convolution with the background and signal terms. These resolution functions have a fixed mean at $\tau _{\mu \mu } = 0$ and free widths.

The free parameters in the fit are the number of signal candidates, the number of background candidates, the non-prompt fraction of signal candidates, the non-prompt fraction of background candidates, the non-prompt fraction of mis-identified candidates, the mean and width of the $J/\psi$ mass peak, the slopes of the exponential distribution functions, and the widths of the pseudo-proper decay time resolution functions.

The relevant quantities extracted from the fit are: the number of signal candidates, $N_\text {signal}$, and the fraction of the signal that is non-prompt, $f_\text {NP}$. These are used to build azimuthal distributions of the prompt and non-prompt yields, as the fits are done in intervals of relative azimuthal angle $2|\phi -\Psi _2|$, $p_{\text {T}}$, y and the collision centrality. Example plots of fit projections are shown in Fig. 1. The prompt and non-prompt signals are obtained from the fit as:

$$\begin{aligned} N_\text {prompt}= & {} N_\text {signal}(1 - f_\text {NP})\text {,} \\ N_\text {non-prompt}= & {} N_\text {signal}f_\text {NP}\text {.} \end{aligned}$$

Fig. 1

Fit projections of the two-dimensional invariant mass ($m_{\mu \mu }$) and pseudo-proper decay time ($\tau _{\mu \mu }$) for the signal extraction for the azimuthal bin $0< 2|\phi -\Psi _2| < \pi /4$ in the kinematic range $9< p_{\text {T}} < 11~\text {GeV}$, $0<|y|<2$ and 0–60% centrality

While the total signal and the non-prompt fraction are weakly correlated, approximately less than 1%, an artificial correlation is introduced when transforming these variables to the prompt and non-prompt yields. The sum of the two yields is constrained to the total number of signal candidates. To compute the correlation factor a toy Monte Carlo model is implemented using the same fit model and the output of the fits to data in bins of relative azimuthal angle, $p_{\text {T}}$, rapidity and centrality. The correlation varies with $p_{\text {T}}$ from $-18\%$ to $-24\%$; in rapidity from $-22\%$ to $-16\%$; and is approximately constant as a function of centrality. The average correlation coefficient is $-20\%$ for all slices and azimuthal bins with a standard deviation of about 2%. It is important to note that this correlation is merely due to the procedure used to extract the yields from the invariant mass and pseudo-proper decay time distributions.

The elliptic flow coefficient is computed by fitting the prompt and non-prompt yields simultaneously to:

$$\begin{aligned} \frac{\mathrm {d} N}{\mathrm {d} (2|\phi -\Psi _2|)} = N_0 \left( 1 + 2 v_2^{\text {fit}} \cos (2|\phi -\Psi _2|)\right) , \end{aligned}$$

(1)

in order to account for the anti-correlation between the two signals. This is achieved by minimising the $\chi ^2$ function:

$$\begin{aligned} \chi ^2(\pmb \theta ) = \left( \pmb {y} - \pmb {\mu }(\pmb \theta )\right) ^{T} V^{-1} \left( \pmb {y} - \pmb \mu (\pmb \theta )\right) , \end{aligned}$$

where $\pmb {y}$ is the vector of measurements, $\pmb {\mu }(\pmb \theta )$ is the vector of predicted values with parameters $\pmb \theta$, and V is the error matrix. The two elements of the vector of measurements are the prompt and non-prompt yields; the vector of predicted values is given by Eq. (1) with the set of free parameters $\{N_{0},v_{2}^\text {fit}\}^\text {prompt}$ and $\{N_{0},v_{2}^\text {fit}\}^\text {non-prompt}$ for the modelling of the prompt and non-prompt yields respectively. The elements in the diagonal of V are the yield uncertainties and the off-diagonal terms are the correlation terms between the prompt and non-prompt yields.

An example of the prompt and non-prompt $J/\psi$ yields normalised by the inclusive $J/\psi$ yield and the projection of the fit result are shown in Fig. 2. The simultaneous fit of the prompt and non-prompt yields correctly accounts for the correlation between the two signals that arose from the modelling used for the signal extraction. The correlation between the fit parameters obtained from the simultaneous fit is shown in Fig. 3.

Fig. 2

The azimuthal distribution of prompt (left) and non-prompt (right) $J/\psi$ yields for the lowest $p_{\text {T}}$ bin studied. The yields are normalised by the inclusive $J/\psi$ signal and the error bars are fit uncertainties associated with the signal extraction. The dotted red line is the result of the simultaneous fit used to compute $v_2$

Fig. 3

Results of the error analysis for the fitted values of the prompt and non-prompt $J/\psi$ $v_2$. The contour lines correspond to the $n\sigma$ fit uncertainties. For this bin, prompt $J/\psi$ $v_2$ has a significance of $3\sigma$ and non-prompt $J/\psi$ has a significance of $1\sigma$

In the final step, the fitted value of $v_{2}^\text {fit}$, is corrected for the event-plane resolution:

$$\begin{aligned} v_2 = v_2^\text {fit}/\mathcal {R}. \end{aligned}$$

3.3 Systematic uncertainties

The systematic uncertainties of this measurement are classified into three groups: (a) related to the centrality definition, (b) related to the estimation of the event-plane method, and (c) related to extraction of the signal. The assigned systematic uncertainty from each source is defined in each bin of $p_{\text {T}}$, rapidity or centrality as the root mean square of the difference between the nominal and varied values of the elliptic flow coefficient. All the uncertainties affecting the extraction of the signal are bin-to-bin uncorrelated, while the uncertainties related to the event-plane method can be correlated or uncorrelated, depending on the studied dependence.

The centrality intervals are defined by values of $\sum E_{\text {T}} ^\text {FCal}$. These intervals have an uncertainty associated primarily to the effect of trigger and event selection inefficiencies as well as backgrounds in the most peripheral $\sum E_{\text {T}} ^\text {FCal}$ intervals [15, 40]. To test the sensitivity of the results to this uncertainty, modified centrality intervals are used, for which the $\sum E_{\text {T}} ^\text {FCal}$ cuts involved in the definition of the centrality intervals are shifted upward and downward, and the analysis is repeated. These changes affect the number of muon pairs entering the signal fitting procedure and thus have an impact on the final value of $v_2$. For the $v_2$ measurements as a function of $p_{\text {T}}$ or rapidity the uncertainty is about $2\%$ for both prompt and non-prompt $J/\psi$ $v_2$, while for the centrality dependence this source contributes a $10\%$ systematic uncertainty to both $v_2$ measurements.

For the estimation of the event-plane angle, the calibration coefficients for the recentring of the flow vector are calculated using a narrower centrality interval (20–60%) instead of the full centrality range (0–60%). For the evaluation of $\delta \Psi _2$ the sum limit is changed to $k_\text {max}=4$. No significant differences are observed, so a systematic uncertainty is not assigned. For the event-plane resolution, the three-sub-event method [19] is used as an alternative to compute $\mathcal {R}$ for the event-plane angle calculated with FCal in $\eta <0$ and $\eta >0$ independently. By using the electromagnetic and hadronic calorimeters, the event-plane angle is calibrated and determined in two different sections with $0.5<\eta <2$ and $-1.5<\eta <0$ and compared with the event-plane angle as measured in the FCal in $\eta <0$ to obtain its resolution. For the resolution of the FCal in the opposite side ($\eta >0$) a reflection of this selection is performed. Both the two-sub-event and three-sub-event methods give consistent results for collisions in the 0–60% centrality interval. To account for the different $\sum E_{\text {T}} ^\text {FCal}$ distributions in minimum-bias and $J/\psi$ triggered events, the difference between the resolutions computed in the two datasets is assigned as a systematic uncertainty and it is the dominant source of uncertainty. The total uncertainty for the average event-plane resolution adds a 4% correlated uncertainty to the measurements integrated over centrality, while for the centrality dependence each point has an approximate $1.5\%$ uncertainty due to resolution. In the centrality interval considered in this analysis (0–60%), it is found that the uncertainty related to the centrality definition has no effect on the event-plane resolution.

Many variations of the PDF defined for the signal extraction are considered. For the mass fit, the most important variations are the release of the fixed parameters of the CB function [39], the substitution of the Gaussian + CB function by a single CB function, and variations of the Gaussian standard deviation and CB width scaling parameter. For the time dependence, the notable variations include that a single exponential function is replaced by short and long lifetime exponential decays, and one Gaussian function instead of two is considered for time resolution. Among all of these variations the biggest contribution is the release of the parameters of the CB function, which contributes between $10\%$ and $15\%$ uncertainty for the $p_{\text {T}}$ and rapidity dependence of $v_2$ and up to $20\%$ for the centrality dependence. Deviations from the case of unpolarised $J/\psi$ production are studied for different spin-alignment scenarios, corresponding to the extreme cases, as explained in Ref. [41]. These alternative scenarios are covered by a theoretical uncertainty of 3% in $v_2$ for prompt $J/\psi$ and 4% for non-prompt $J/\psi$.

The choice of mass window used for the study is changed to analyse potential biases from the mass peak of the $\psi (2\mathrm {S})$. At high rapidity its width increases and the fit response for the background estimation changes. The mass window range is narrowed to $2.7< m_{\mu \mu } < 3.4~\text {GeV}$ and its impact is between $5\%$ and $10\%$ for the $p_{\text {T}}$, rapidity and centrality dependencies

The correlation between the prompt and non-prompt yields is also studied. It is either doubled or neglected, and shows a minor impact of $1\%$ for all presented results.

4 Results

Results of $v_2$ measurements for prompt and non-prompt $J/\psi$ are shown as a function of $p_{\text {T}}$ in the range between 9 and 30 $\text {GeV}$ in Fig. 4, for three $p_{\text {T}}$ intervals. The centroid of each $p_{\text {T}}$ bin is determined by the average $p_{\text {T}}$ of the muon pairs in the corresponding bin and their values are 9.81, 12.17, and $17.6~\text {GeV}$. The horizontal error bars correspond to the bin width reflecting the kinematic range of the measurement. The vertical error bars are the fit errors associated with statistical uncertainties, and the shaded boxes are the systematic uncertainties. The data are consistent with a non-zero flow signal in the full kinematic range studied ($9<p_{\text {T}} <30~\text {GeV}$) for both prompt and non-prompt $J/\psi$. As shown in Fig. 3, for the lowest $p_{\text {T}}$ bin ($9< p_{\text {T}} < 11~\text {GeV}$), prompt $J/\psi$ $v_2$ deviates from 0 with a significance of $3\sigma$ and non-prompt $J/\psi$ with a significance of $1\sigma$. Prompt $J/\psi$ exhibits a decreasing trend with a maximum value of $v_2$ close to 0.09 that decreases by nearly a factor of two over the whole studied kinematic range. The results for non-prompt $J/\psi$ indicate a non-zero value with limited statistical significance. These $v_2$ values are consistent with being independent of $p_{\text {T}}$ and compatible within uncertainties with the $v_2$ values of prompt $J/\psi$, particularly at the highest $p_{\text {T}}$.

The rapidity dependence of $v_2$ is shown in Fig. 5 and the centrality dependence in Fig. 6 for both prompt and non-prompt $J/\psi$. Neither shows significant rapidity or centrality dependence. The prompt $J/\psi$ $v_2$ is larger than the non-prompt, in agreement with the larger values observed in the $p_{\text {T}}$ dependence integrated over rapidity and centrality. The measured value of $v_2$ for prompt $J/\psi$ stays approximately the same in central collisions as in non-central collisions within uncertainties, in agreement with the observation of Ref. [5]. This is similar to the case of non-prompt $J/\psi$ where no evident centrality dependence is observed within the uncertainties. This feature is in disagreement with the expected hydrodynamic behaviour for charm quarks and may manifest a transition at medium $p_{\text {T}}$ regime where there are thought to be different effects influencing $J/\psi$ production [6, 7, 9,10,11].

In Fig. 7 the available results for inclusive $J/\psi$ ($p_{\text {T}} <12~\text {GeV}$) from the ALICE experiment [4] and prompt $J/\psi$ (4 < $p_{\text {T}}$ < 30 $\text {GeV}$) from the CMS experiment [5] are compared with the results obtained in this analysis for prompt $J/\psi$ ($9< p_{\text {T}} < 30~\text {GeV}$) as a function of the $J/\psi$ transverse momentum. Despite different rapidity selections, the magnitudes of the elliptic flow coefficients are compatible with each other. Two features can be observed: first, the hydrodynamic peak is around 7 $\text {GeV}$, a value that is significantly higher than what is observed for charged particles [14,15,16,17,18] where the peak is around 3–4 $\text {GeV}$. This effect can be described qualitatively by thermalisation of charm quarks in the quark–gluon plasma medium with $J/\psi$ regeneration playing a dominant role in the flow formation [6, 7]. The second feature is that $v_2$ has a substantial magnitude at high $p_{\text {T}}$. This can be connected with the suppression of $J/\psi$ production due to mechanisms involving interactions with the medium such as absorption and melting [11] or energy loss [42, 43].

Fig. 4

Prompt (left) and non-prompt (right) $J/\psi$ $v_2$ as a function of transverse momentum for the rapidity interval $|y| < 2$ and centrality 0–60%. The statistical and systematic uncertainties are shown using vertical error bars and boxes respectively. The horizontal error bars represent the kinematic range of the measurement for each bin

Fig. 5

Prompt (left) and non-prompt (right) $J/\psi$ $v_2$ as a function of rapidity for transverse momentum in the range $9< p_{\text {T}} < 30~\text {GeV}$ and centrality 0–60%. The statistical and systematic uncertainties are shown using vertical error bars and boxes respectively. The horizontal error bars represent the kinematic range of the measurement for each bin

Fig. 6

Prompt (left) and non-prompt (right) $J/\psi$ $v_2$ as a function of average number of nucleons participating in the collision for transverse momentum in the range $9< p_{\text {T}} < 30~\text {GeV}$ and rapidity $|y| < 2$. The statistical and systematic uncertainties are shown using vertical error bars and boxes respectively. The centrality interval associated to a given value of $\langle N_{\text {part}} \rangle$ is written below each data point

Fig. 7

Results for $v_2$ as a function of the transverse momentum of prompt $J/\psi$ as measured by ATLAS in this analysis compared with inclusive $J/\psi$ with $p_{\text {T}} < 12~\text {GeV}$ as measured by ALICE at $5.02~\text {TeV}$ [4], and prompt $J/\psi$ with $p_{\text {T}}$ in the range $4< p_{\text {T}} < 30~\text {GeV}$ by CMS at $2.76~\text {TeV}$ [5]. The statistical and systematic uncertainties are shown using vertical error bars and boxes respectively

5 Summary

This paper presents measurements of the elliptic flow harmonic coefficients for $J/\psi$ particles in the dimuon decay channel in $0.42~\mathrm {nb}^{-1}$ of Pb+Pb collisions recorded at $\sqrt{s_{_\text {NN}}} =5.02~\text {TeV}$ with the ATLAS detector at the LHC. Results are presented for prompt and non-prompt $J/\psi$ as a function of transverse momentum, rapidity and centrality. The measurement is performed in the $J/\psi$ kinematic range $9<p_{\text {T}} <30~\text {GeV}$, $|y|<2$, and 0–60% centrality. The pseudo-proper decay time of the secondary vertex is used to separate the prompt and non-prompt components of $J/\psi$ production and both yields are analysed simultaneously to properly assess the correlation between the two contributions.

A significant flow signal is found for prompt $J/\psi$, which decreases with increasing $p_{\text {T}}$. With limited statistical significance, it is found that non-prompt $J/\psi$ $v_2$ is consistent with a flat behaviour over the studied $p_{\text {T}}$ range. At high $p_{\text {T}}$, the prompt and non-prompt $J/\psi$ $v_2$ values are compatible within the uncertainties. There is no evidence for a rapidity or centrality dependence for the prompt or non-prompt case. This suggests a similar underlying process describing the propagation of sufficiently high $p_{\text {T}}$ charm and bottom quarks through the medium. The idea is supported by the recent observation of $J/\psi$ yield suppression in Pb+Pb collisions by ATLAS, where a similar suppression pattern for prompt and non-prompt $J/\psi$ is observed at high $p_{\text {T}}$. Additionally, this measurement covers the high $p_{\text {T}}$ range of $J/\psi$ and is found to be in a good agreement with previous reports, despite the different beam energy and rapidity selections.