Measurement of the differential cross-sections of prompt and non-prompt production of $J/\psi$ and $\psi(2\mathrm{S})$ in $pp$ collisions at $\sqrt{s} = 7$ and $8$ TeV with the ATLAS detector

The production rates of prompt and non-prompt $J/\psi$ and $\psi(2\mathrm{S})$ mesons are measured using 2.1 $fb^{-1}$ and 11.4 $fb^{-1}$ of data collected with the ATLAS experiment at the LHC, in proton-proton collisions at $\sqrt{s}=7$ and 8 TeV respectively. Production cross-sections for both prompt and non-prompt production sources, ratios of $\psi(2\mathrm{S})$ to $J/\psi$ production, and fractions of non-prompt to inclusive production for $J/\psi$ and $\psi(2\mathrm{S})$ are measured double-differentially as a function of meson $p_{T}$ and rapidity. These measurements are made in a restricted fiducial volume and also corrected for geometrical acceptance after which they are compared to a variety of theoretical predictions.


Introduction
Measurements of heavy quark-antiquark bound states (quarkonia) production processes provide an insight into the nature of quantum chromodynamics (QCD) close to the boundary between the perturbative and non-perturbative regimes. More than forty years since the discovery of the J/ψ, the investigation of hidden heavy-flavour production in hadronic collisions still presents significant challenges to both theory and experiment.
In high-energy hadronic collisions, charmonium states can be produced either directly by short-lived QCD sources ("prompt" production), or by long-lived sources in the decay chains of beauty hadrons ("nonprompt" production). These can be separated experimentally using the distance between the protonproton primary interaction and the decay vertex of the quarkonium state. While Fixed-Order with Nextto-Leading-Log (FONLL) calculations [1,2], made within the framework of perturbative QCD, have been quite successful in describing non-prompt production of various quarkonium states, a satisfactory understanding of the prompt production mechanisms is still to be achieved.
Early attempts to describe the formation of charmonium [25][26][27][28][29][30][31][32] using leading-order perturbative QCD gave rise to a variety of models, none of which could explain the large production cross-sections measured at the Tevatron [3,13,[21][22][23]. Within the colour-singlet model (CSM) [33], next-to-next-to-leadingorder (NNLO) contributions to the hadronic production of S-wave quarkonia were calculated without introducing any new phenomenological parameters. However, technical difficulties have so far made it impossible to perform the full NNLO calculation, or to extend those calculations to the P-wave states. So it is not entirely surprising that the predictions of the model underestimate the experimental data for inclusive production of J/ψ and Υ states, where the feed-down is significant, but offer a better description for ψ(2S) production [18,34].
Non-relativistic QCD (NRQCD) calculations that include colour-octet (CO) contributions [35] introduce a number of phenomenological parameters -long-distance matrix elements (LDMEs) -which are determined from fits to the experimental data, and can hence describe the cross-sections and differential spectra satisfactorily [36]. However, the attempts to describe the polarization of S-wave quarkonium states using this approach have not been so successful [37], prompting a suggestion [38] that a more coherent approach is needed for the treatment of polarization within the QCD-motivated models of quarkonium production.
Neither the CSM nor the NRQCD model gives a satisfactory explanation for the measurement of prompt J/ψ production in association with the W [39] and Z [40] bosons: in both cases, the measured differential cross-section is larger than theoretical expectations [41][42][43][44]. It is therefore important to broaden the scope of comparisons between theory and experiment by providing a variety of experimental information about quarkonium production across a wider kinematic range. In this context, ATLAS has measured the inclusive differential cross-section of J/ψ production, with 2.3 pb −1 of integrated luminosity [18], at √ s = 7 TeV using the data collected in 2010, as well as the differential cross-sections of the production of χ c states (4.5 fb −1 ) [14], and of the ψ(2S) in its J/ψππ decay mode (2.1 fb −1 ) [9], at √ s = 7 TeV with data collected in 2011. The cross-section and polarization measurements from CDF [4], CMS [7,45,46], LHCb [8,10,12,[47][48][49] and ALICE [5,50,51], cover a considerable variety of charmonium production characteristics in a wide kinematic range (transverse momentum p T ≤ 100 GeV and rapidities |y| < 5), thus providing a wealth of information for a new generation of theoretical models. This paper presents a precise measurement of J/ψ and ψ(2S) production in the dimuon decay mode, both at √ s = 7 TeV and at √ s = 8 TeV. It is presented as a double-differential measurement in transverse momentum and rapidity of the quarkonium state, separated into prompt and non-prompt contributions, covering a range of transverse momenta 8 < p T ≤ 110 GeV and rapidities |y| < 2.0. The ratios of ψ(2S) to J/ψ cross-sections for prompt and non-prompt processes are also reported, as well as the non-prompt fractions of J/ψ and ψ(2S).

The ATLAS detector
The ATLAS experiment [52] is a general-purpose detector consisting of an inner tracker, a calorimeter and a muon spectrometer. The inner detector (ID) directly surrounds the interaction point; it consists of a silicon pixel detector, a semiconductor tracker and a transition radiation tracker, and is embedded in an axial 2 T magnetic field. The ID covers the pseudorapidity 1 range |η| = 2.5 and is enclosed by a calorimeter system containing electromagnetic and hadronic sections. The calorimeter is surrounded by a large muon spectrometer (MS) in a toroidal magnet system. The MS consists of monitored drift tubes and cathode strip chambers, designed to provide precise position measurements in the bending plane in the range |η| < 2.7. Momentum measurements in the muon spectrometer are based on track segments formed in at least two of the three precision chamber planes. The ATLAS trigger system [53] is separated into three levels: the hardware-based Level-1 trigger and the two-stage High Level Trigger (HLT), comprising the Level-2 trigger and Event Filter, which reduce the 20 MHz proton-proton collision rate to several-hundred Hz of events of interest for data recording to mass storage. At Level-1, the muon trigger searches for patterns of hits satisfying different transverse momentum thresholds with a coarse position resolution but a fast response time using resistive-plate chambers and thin-gap chambers in the ranges |η| < 1.05 and 1.05 < |η| < 2.4, respectively. Around these Level-1 hit patterns "Regions-of-Interest" (RoI) are defined that serve as seeds for the HLT muon reconstruction. The HLT uses dedicated algorithms to incorporate information from both the MS and the ID, achieving position and momentum resolution close to that provided by the offline muon reconstruction.

Candidate selection
The analysis is based on data recorded at the LHC in 2011 and 2012 during proton-proton collisions at centre-of-mass energies of 7 TeV and 8 TeV, respectively. This data sample corresponds to a total 1 ATLAS uses a right-handed coordinate system with its origin at the nominal interaction point (IP) in the centre of the detector and the z-axis along the beam pipe. The x-axis points from the IP to the centre of the LHC ring, and the y-axis points upward. Cylindrical coordinates (r, φ) are used in the transverse plane, φ being the azimuthal angle around the beam pipe.
The pseudorapidity η is defined in terms of the polar angle θ as η = − ln tan(θ/2) and the transverse momentum p T is defined as p T = p sin θ. The rapidity is defined as y = 0.5 ln (E + p z ) / (E − p z ) , where E and p z refer to energy and longitudinal momentum, respectively. The η-φ distance between two particles is defined as ∆R = (∆η) 2 + (∆φ) 2 .
Events were selected using a trigger requiring two oppositely charged muon candidates, each passing the requirement p T > 4 GeV. The muons are constrained to originate from a common vertex, which is fitted with the track parameter uncertainties taken into account. The fit is required to satisfy χ 2 < 20 for the one degree of freedom.
For 7 TeV data, the Level-1 trigger required only spatial coincidences in the MS [54]. For 8 TeV data, a 4 GeV muon p T threshold was also applied at Level-1, which reduced the trigger efficiency for low-p T muons.
The offline analysis requires events to have at least two muons, identified by the muon spectrometer and with matching tracks reconstructed in the ID [55]. Due to the ID acceptance, muon reconstruction is possible only for |η| < 2.5. The selected muons are further restricted to |η| < 2.3 to ensure high-quality tracking and triggering, and to reduce the contribution from misidentified muons. For the momenta of interest in this analysis (corresponding to muons with a transverse momentum of at most O(100) GeV), measurements of the muons are degraded by multiple scattering within the MS and so only the ID tracking information is considered. To ensure accurate ID measurements, each muon track must fulfil muon reconstruction and selection requirements [55]. The pairs of muon candidates satisfying these quality criteria are required to have opposite charges.
In order to allow an accurate correction for trigger inefficiencies, each reconstructed muon candidate is required to match a trigger-identified muon candidate within a cone of ∆R = (∆η) 2 + (∆φ) 2 = 0.01. Dimuon candidates are obtained from muon pairs, constrained to originate from a common vertex using ID track parameters and uncertainties, with a requirement of χ 2 < 20 of the vertex fit for the one degree of freedom. All dimuon candidates with an invariant mass within 2.6 < m(µµ) < 4.0 GeV and within the kinematic range p T (µµ) > 8 GeV, |y(µµ)| < 2.0 are retained for the analysis. If multiple candidates are found in an event (occurring in approximately 10 −6 of selected events), all candidates are retained. The properties of the dimuon system, such as invariant mass m(µµ), transverse momentum p T (µµ), and rapidity |y(µµ)| are determined from the result of the vertex fit.

Methodology
The measurements are performed in intervals of dimuon p T and absolute value of the rapidity (|y|). The term "prompt" refers to the J/ψ or ψ(2S) states -hereafter called ψ to refer to either -are produced from short-lived QCD decays, including feed-down from other charmonium states as long as they are also produced from short-lived sources. If the decay chain producing a ψ state includes long-lived particles such as b-hadrons, then such ψ mesons are labelled as "non-prompt". Using a simultaneous fit to the invariant mass of the dimuon and its "pseudo-proper decay time" (described below), prompt and nonprompt signal and background contributions can be extracted from the data.
The probability for the decay of a particle as a function of proper decay time t follows an exponential distribution, p(t) = 1/τ B · e −t/τ B where τ B is the mean lifetime of the particle. For each decay, the proper decay time can be calculated as t = Lm/p, where L is the distance between the particle production and decay vertices, p is the momentum of the particle, and m is its invariant mass. As the reconstruction of non-prompt ψ mesons, such as b-hadrons, does not fully describe the properties of the parent, the transverse momentum of the dimuon system and the reconstructed dimuon invariant mass are used to construct the "pseudo-proper decay time", τ = L xy m(µµ)/p T (µµ), where L xy ≡ L · p T (µµ)/p T (µµ) is the signed projection of the distance of the dimuon decay vertex from the primary vertex, L, onto its transverse momentum, p T (µµ). This is a good approximation of using the parent b-hadron information when the ψ and parent momenta are closely aligned, which is the case for the values of ψ transverse momenta considered here, and τ therefore can be used to distinguish statistically between the non-prompt and prompt processes (in which the latter are assumed to decay with vanishingly small lifetime). If the event contains multiple primary vertices [52], the primary vertex closest in z to the dimuon decay vertex is selected. The effect of selecting an incorrect vertex has been shown [56] to have a negligible impact on the extraction of prompt and non-prompt contributions. If any of the muons in the dimuon candidate contributes to the construction of the primary vertex, the corresponding tracks are removed and the vertex is refitted.

Double differential cross-section determination
The double differential dimuon prompt and non-prompt production cross-sections times branching ratio are measured separately for J/ψ and ψ(2S) mesons according to the equations: where Ldt is the integrated luminosity, ∆p T and ∆y are the interval sizes in terms of dimuon transverse momentum and rapidity, respectively, and N p(np) ψ is the number of observed prompt (non-prompt) ψ mesons in the slice under study, corrected for acceptance, trigger and reconstruction efficiencies. The intervals in ∆y combine the data from negative and positive rapidities.
The determination of the cross-sections proceeds in several steps. First, a weight is determined for each selected dimuon candidate equal to the inverse of the total efficiency for each candidate. The total weight, w tot , for each dimuon candidate includes three factors: the fraction of produced ψ → µ + µ − decays with both muons in the fiducial region p T (µ) > 4 GeV and |η(µ)| < 2.3 (defined as acceptance, A), the probability that a candidate within the acceptance satisfies the offline reconstruction selection (ǫ reco ), and the probability that a reconstructed event satisfies the trigger selection (ǫ trig ). The weight assigned to a given candidate when calculating the cross-sections is therefore given by: After the weight determination, an unbinned maximum-likelihood fit is performed to these weighted events in each (p T (µµ), |y(µµ)|) interval using the dimuon invariant mass, m(µµ), and pseudo-proper decay time, τ(µµ), observables. The fitted yields of J/ψ → µ + µ − and ψ(2S) → µ + µ − are determined separately for prompt and non-prompt processes. Finally, the differential cross-section times the ψ → µ + µ − branching fraction is calculated for each state by including the integrated luminosity and the p T and rapidity interval widths as shown in Eqs.

Non-prompt fraction
The non-prompt fraction f ψ b is defined as the number of non-prompt ψ (produced via the decay of a bhadron) divided by the number of inclusively produced ψ decaying to muon pairs after applying weighting corrections: where this fraction is determined separately for J/ψ and ψ(2S). Determining the fraction from this ratio is advantageous since acceptance and efficiencies largely cancel and the systematic uncertainty is reduced.

Ratio of ψ(2S) to J/ψ production
The ratio of ψ(2S) to J/ψ production, in their dimuon decay modes, is defined as: is the number of prompt (non-prompt) J/ψ or ψ(2S) mesons decaying into a muon pair in an interval of p T and y, corrected for selection efficiencies and acceptance.
For the ratio measurements, similarly to the non-prompt fraction, the acceptance and efficiency corrections largely cancel, thus allowing a more precise measurement. The theoretical uncertainties on such ratios are also smaller, as several dependencies, such as parton distribution functions and b-hadron production spectra, largely cancel in the ratio.

Acceptance
The kinematic acceptance A for a ψ → µ + µ − decay with p T and y is given by the probability that both muons pass the fiducial selection (p T (µ) > 4 GeV and |η(µ)| < 2.3). This is calculated using generator-level "accept-reject" simulations, based on the analytic formula described below. Detectorlevel corrections, such as bin migration effects due to detector resolution, are found to be small. They are applied to the results and are also considered as part of the systematic uncertainties.
The acceptance A depends on five independent variables (the two muon momenta are constrained by the m(µµ) mass condition), chosen as the p T , |y| and azimuthal angle φ of the ψ meson in the laboratory frame, and two angles characterizing the ψ → µ + µ − decay, θ ⋆ and φ ⋆ , described in detail in Ref. [57]. The angle θ ⋆ is the angle between the direction of the positive-muon momentum in the ψ rest frame and the momentum of the ψ in the laboratory frame, while φ ⋆ is defined as the angle between the dimuon production and decay planes in the laboratory frame. The ψ production plane is defined by the momentum of the ψ in the laboratory frame and the positive z-axis direction. The distributions in θ ⋆ and φ ⋆ differ for various possible spin-alignment scenarios of the dimuon system.
The spin-alignment of the ψ may vary depending on the production mechanism, which in turn affects the angular distribution of the dimuon decay. Predictions of various theoretical models are quite contradictory, while the recent experimental measurements [7] indicate that the angular dependence of J/ψ and ψ(2S) decays is consistent with being isotropic.
The coefficients λ θ , λ φ and λ θφ in are related to the spin-density matrix elements of the dimuon spin wave function.
Since the polarization of the ψ state may affect acceptance, seven extreme cases that lead to the largest possible variations of acceptance within the phase space of this measurement are identified. These cases, described in Table 1, are used to define a range in which the results may vary under any physically allowed spin-alignment assumptions. The same technique has also been used in other measurements [9,14,34]. This analysis adopts the isotropic distribution in both cos θ ⋆ and φ ⋆ as nominal, and the variation of the results for a number of extreme spin-alignment scenarios is studied and presented as sets of correction factors, detailed further in Appendix A.

Angular coefficients
For each of the two mass-points (corresponding to the J/ψ and ψ(2S) masses), two-dimensional maps are produced as a function of dimuon p T (µµ) and |y(µµ)| for the set of spin-alignment hypotheses. Each point on the map is determined from a uniform sampling over φ ⋆ and cos θ ⋆ , accepting those trials that pass the fiducial selections. To account for various spin-alignment scenarios, all trials are weighted according to Eq. 3. Acceptance maps are defined within the range 8 < p T (µµ) < 110 GeVand |y(µµ)| < 2.0, corresponding to the data considered in the analysis. The map is defined by 100 slices in |y(µµ)| and 4400 in p T (µµ), using 200k trials for each point, resulting in sufficiently high precision that the statistical uncertainty can be neglected. Due to the contributions of background, and the detector resolution of the signal, the acceptance for each candidate is determined from a linear interpolation of the two maps, which are generated for the J/ψ and ψ(2S) known masses, as a function of the reconstructed mass m(µµ). Figure 1 shows the acceptance, projected in p T for all the spin-alignment hypotheses for the J/ψ meson. The differences between the acceptance of the ψ(2S) and J/ψ meson, are independent of rapidity, except near |y| ≈ 2 at low p T . Similarly, the only dependence on p T is found below p T ≈ 9 GeV. The correction factors (as given in Appendix. A) vary most at low p T , ranging from −35% under longitudinal, to +100% for transverse-positive scenarios. At high p T , the range is between −14% for longitudinal, and +9% for transverse-positive scenarios. For the fraction and ratio measurements, the correction factor is determined from the appropriate ratio of the individual correction factors. [GeV]

Muon reconstruction and trigger efficiency determination
The technique for correcting the 7 TeV data for trigger and reconstruction inefficiencies is described in detail in Ref. [9,34]. For the 8 TeV data, a similar technique is used, however different efficiency maps are required for each set of data, and the 8 TeV corrections are detailed briefly below.
The single-muon reconstruction efficiency is determined from a tag-and-probe study in dimuon decays [40]. The efficiency map is calculated as a function of p T (µ) and q × η(µ), where q = ±1 is the electrical charge of the muon, expressed in units of e.
The trigger efficiency correction consists of two components. The first part represents the trigger efficiency for a single muon in intervals of p T (µ) and q × η(µ). For the dimuon system there is a second correction to account for reductions in efficiency due to closely spaced muons firing only a single RoI, vertex-quality cuts, and opposite-sign requirements. This correction is performed in three rapidity intervals: 0-1.0, 1.0-1.2 and 1.2-2.3. The correction is a function of ∆R(µµ) in the first two rapidity intervals and a function of ∆R(µµ) and |y(µµ)| in the last interval.
The combination of the two components (single-muon efficiency map and dimuon corrections) is illustrated in Figure 2 by plotting the average trigger-weight correction for the events in this analysis in terms of p T (µµ) and |y(µµ)|. The increased weight at low p T and |y| ≈ 1.25 is caused by the geometrical acceptance of the muon trigger system and the turn-on threshold behaviour of the muon trigger. At high p T the weight is increased due to the reduced opening angle between the two muons.

Fitting technique
To extract the corrected yields of prompt and non-prompt J/ψ and ψ(2S) mesons, two-dimensional weighted unbinned maximum-likelihood fits are performed on the dimuon invariant mass, m(µµ), and pseudo-proper decay time, τ(µµ), in intervals of p T (µµ) and |y(µµ)|. Each interval is fitted independently from all the others. In m(µµ), signal processes of ψ meson decays are statistically distinguished as narrow peaks convolved with the detector resolution, at their respective mass positions, on top of background continuum. In τ(µµ), decays originating with zero pseudo-proper decay time and those following an exponential decay distribution (both convolved with a detector resolution function) statistically distinguish prompt and non-prompt signal processes, respectively. Various sources of background processes include Drell-Yan processes, mis-reconstructed muon pairs from prompt and non-prompt sources, and semileptonic decays from separate b-hadrons.
The probability density function (PDF) for each fit is defined as a normalized sum, where each term represents a specific signal or background contribution, with a physically motivated mass and τ dependence. The PDF can be written in a compact form as where κ i represents the relative normalization of the i th term of the seven considered signal and background contributions (such that i κ i = 1), f i (m) is the mass-dependent term, and ⊗ represents the convolution of the τ-dependent function h i (τ) with the τ resolution term, R(τ). The latter is modelled by a double Gaussian distribution with both means fixed to zero and widths determined from the fit. have different decay constants, where E 5 (|τ|) is a double-sided exponential with the same decay constant on either side of τ = 0. The parameter ω represents the fractional contribution of the B and G mass signal functions, while the Dirac delta function, δ(τ), is used to represent the pseudo-proper decay time distribution of the prompt candidates. Table 2: Description of the fit model PDF in Eq. 4. Components of the probability density function used to extract the prompt (P) and non-prompt (NP) contributions for J/ψ and ψ(2S) signal and the P, NP, and incoherent or mis-reconstructed background (Bkg) contributions. i In order to make the fitting procedure more robust and to reduce the number of free parameters, a number of component terms share common parameters, which led to 22 free parameters per interval. In detail, the signal mass models are described by the sum of a Crystal Ball shape (B) and a Gaussian shape (G). For each of J/ψ and ψ(2S), the B and G share a common mean, and freely determined widths, with the ratio of the B and G widths common to J/ψ and ψ(2S). The B parameters α, and n, describing the transition point of the low-edge from a Gaussian to a power-law shape, and the shape of the tail, respectively, are fixed, and variations are considered as part of the fit model systematic uncertainties. The width of G for ψ(2S) is set to the width for J/ψ multiplied by a free parameter scaling term. The relative fraction of B and G is left floating, but common to J/ψ and ψ(2S).
The non-prompt signal decay shapes (E 1 ,E 2 ) are described by an exponential function (for positive τ only) convolved with a double Gaussian function, R(τ) describing the pseudo-proper decay time resolution for the non-prompt component, and the same Gaussian response functions to describe the prompt contributions. Each Gaussian resolution component has its mean fixed at τ = 0 and a free width. The decay constants of the J/ψ and ψ(2S) are separate free parameters in the fit.
The background contributions are described by a prompt and non-prompt component, as well as a doublesided exponential function convolved with a double Gaussian function describing mis-reconstructed or non-coherent muon pairs. The same resolution function as in signal is used to describe the background. For the non-resonant mass parameterizations, the non-prompt contribution is modelled by a first-order Chebyshev polynomial. The prompt mass contribution follows a flat distribution and the double-sided background uses an exponential function. Variations of this fit model are considered as systematic uncertainties.
The following quantities are extracted directly from the fit in each interval: the fraction of events that are signal (prompt or non-prompt J/ψ or ψ(2S)); the fraction of signal events that are prompt; the fraction of prompt signal that is ψ(2S); and the fraction of non-prompt signal that is ψ(2S). From these parameters, and the weighted sum of events, all measured values are calculated.
For 7 TeV data, 168 fits are performed across the range of 8 < p T < 100 GeV (8 < p T < 60 GeV) for J/ψ (ψ(2S)) and 0 < |y| < 2. For 8 TeV data, 172 fits are performed across the range of 8 < p T < 110 GeV and 0 < |y| < 2, excluding the area where p T is less than 10 GeV and simultaneously |y| is greater than 0.75. This region is excluded due to a steeply changing low trigger efficiency causing large systematic uncertainties in the measured cross-section. Figure 3 shows the fit results for one of the intervals considered in the analysis, projected onto the invariant mass and pseudo-proper decay time distributions, for 7 TeV data, weighted according to the acceptance and efficiency corrections. The fit projections are shown for the total prompt and total non-prompt contributions (shown as curves), and also for the individual contributions of the J/ψ and ψ(2S) prompt and non-prompt signal yields (shown as hashed areas of various types).
In Figure 4 the fit results are shown for one high-p T interval of 8 TeV data.

Bin migration corrections
To account for bin migration effects due to the detector resolution, which results in decays of ψ in one bin, being identified and accounted for in another, the numbers of acceptance-and efficiency-corrected dimuon decays extracted from the fits in each interval of p T (µµ) and rapidity are corrected for the differences between the true and reconstructed values of the dimuon p T . These corrections are derived from data by comparing analytic functions that are fitted to the p T (µµ) spectra of dimuon events with and without convolution by the experimental resolution in p T (µµ) (as determined from the fitted mass resolution and measured muon angular resolutions), as described in Ref. [34]. The correction factors applied to the fitted yields deviate from unity by no more than 1.5%, and for the majority of slices are smaller than 1%. The ratio measurement and non-prompt fractions are corrected by the corresponding ratios of bin migration correction factors. Using a similar technique, bin migration corrections as a function of |y| are found to differ from unity by negligible amounts.

Systematic uncertainties
The sources of systematic uncertainties that are applied to the ψ double differential cross-section measurements are from uncertainties in: the luminosity determination; muon and trigger efficiency corrections; inner detector tracking efficiencies; the fit model parametrization; and due to bin migration corrections. For the non-prompt fraction and ratio measurements the systematic uncertainties are assessed in the same manner as for the uncertainties on the cross-section, except that in these ratios some systematic uncertainties, such as the luminosity uncertainty, cancel out. The sources of systematic uncertainty evaluated for the prompt and non-prompt ψ cross-section measurements, along with the minimum, maximum and median values, are listed in Table 3. The largest contributions, which originate from the trigger and fit model uncertainties, are typically for the high p T intervals and are due to the limited statistics of the efficiency maps (for the trigger), and the data sample (for the fit model). Figures 5 and 6 show, for a representative interval, the impact of the considered uncertainties on the production cross-section, as well as the non-prompt fraction and ratios for 7 TeV data. The impact is very similar at 8 TeV.

Luminosity
The uncertainty on the integrated luminosity is 1.8% (2.8%) for the 7 TeV (8 TeV) data-taking period. The methodology used to determine these uncertainties is described in Ref. [59]. The luminosity uncertainty is only applied to the J/ψ and ψ(2S) cross-section results.

Muon reconstruction and trigger efficiencies
To determine the systematic uncertainty on the muon reconstruction and trigger efficiency maps, each of the maps is reproduced in 100 pseudo-experiments. The dominant uncertainty in each bin is statistical and hence any bin-to-bin correlations are neglected. For each pseudo-experiment a new map is created by varying independently each bin content according to a Gaussian distribution about its estimated value, determined from the original map. In each pseudo-experiment, the total weight is recalculated for each dimuon p T and |y| interval of the analysis. The RMS of the total weight pseudo-experiment distributions for each efficiency type is used as the systematic uncertainty, where any correlation effects between the muon and trigger efficiencies can be neglected.
The ID tracking efficiency is in excess of 99.5% [34], and an uncertainty of 1% is applied to account for the ID dimuon reconstruction inefficiency (0.5% per muon, added coherently). This uncertainty is applied to the differential cross-sections and is assumed to cancel in the fraction of non-prompt to inclusive production for J/ψ and ψ(2S) and in the ratios of ψ(2S) to J/ψ production.
For the trigger efficiency ǫ trig , in addition to the trigger efficiency map, there is an additional correction term that accounts for inefficiencies due to correlations between the two trigger muons, such as the dimuon opening angle. This correction is varied by its uncertainty, and the shift in the resultant total weight relative to its central value is added in quadrature to the uncertainty from the map. The choice of triggers is known [60] to introduce a small lifetime-dependent efficiency loss but it is determined to have a negligible effect on the prompt and non-prompt yields and no correction is applied in this analysis. Similarly, the muon reconstruction efficiency corrections of prompt and non-prompt signals are found to be consistent within the statistical uncertainties of the efficiency measurements, and no additional uncertainty is applied.

Fit model uncertainty
The uncertainty due to the fit procedure is determined by varying one component at a time in the fit model described in Section 4.6, creating a set of new fit models. For each new fit model, all measured quantities are recalculated, and in each p T and |y| interval the spread of variations around the central fit model is used as its systematic uncertainty. The variations of the fit model also account for possible uncertainties due to final-state radiation. The following variations to the central model fit are evaluated: • signal mass model -using double Gaussian models in place of the Crystal Ball plus Gaussian model; variation of the α and n parameters of the B model, which are originally fixed; • signal pseudo-proper decay time model -a double exponential function is used to describe the pseudo-proper decay time distribution for the ψ non-prompt signal; • background mass models -variations of the mass model using exponentials functions, or quadratic Chebyshev polynomials to describe the components of prompt, non-prompt and double-sided background terms; • background pseudo-proper decay time model -a single exponential function was considered for the non-prompt component; • pseudo-proper decay time resolution model -using a single Gaussian function in place of the double Gaussian function to model the lifetime resolution (also prompt lifetime model); and variation of the mixing terms for the two Gaussian components of this term.
Of the variations considered, it is typically the parametrizations of the signal mass model and pseudoproper decay time resolution model that dominate the contribution to the fit model uncertainty.

Bin migrations
As the corrections to the results due to bin migration effects are factors close to unity in all regions, the difference between the correction factor and unity is applied as the uncertainty.
The variation of the acceptance corrections with spin-alignment is treated separately, and scaling factors supplied in Appendix A.

Results
The J/ψ and ψ(2S) non-prompt and prompt production cross-sections are presented, corrected for acceptance and detector efficiencies while assuming isotropic decay, as described in Section 4.1. Also presented are the ratios of non-prompt production relative to the inclusive production for J/ψ and ψ(2S) mesons separately, described in Section 4.2, and the ratio of ψ(2S) to J/ψ production for prompt and nonprompt components separately, described in Section 4.3. Correction factors for various spin-alignment hypotheses for both 7 and 8 TeV data can be found in Tables 4-15 (in Appendix) and Tables 16-27 (in Appendix) respectively, in terms of p T and rapidity intervals.

Production cross-sections
Figures 7 and 8 show respectively the prompt and non-prompt differential cross-sections of J/ψ and ψ(2S) as functions of p T and |y|, together with the relevant theoretical predictions, which are described below.

Non-prompt production fractions
The results for the fractions of non-prompt production relative to the inclusive production of J/ψ and ψ(2S) are presented as a function of p T for slices of rapidity in Figure 9. In each rapidity slice, the nonprompt fraction is seen to increase as a function of p T and has no strong dependence on either rapidity or centre-of-mass energy.
Production ratios of ψ(2S) to J/ψ Figure 10 shows the ratios of ψ(2S) to J/ψ decaying to a muon pair in prompt and non-prompt processes, presented as a function of p T for slices of rapidity. The non-prompt ratio is shown to be relatively flat across the considered range of p T , for each slice of rapidity. For the prompt ratio, a slight increase as a function of p T is observed, with no strong dependence on rapidity or centre-of-mass energy. [nb GeV    [nb GeV [nb GeV  [nb GeV    [nb GeV [nb GeV      Figure 10: The ratio of ψ(2S) to J/ψ production times dimuon branching fraction for prompt (left) and non-prompt (right) processes as a function of p T (µµ) for each of the slices of rapidity. For each increasing rapidity slice, an additional factor of 0.1 is applied to the plotted points for visual clarity. The top (bottom) row shows the 7 TeV (8 TeV) results. The centre of each bin on the horizontal axis represents the mean of the weighted p T distribution. The horizontal error bars represent the range of p T for the bin, and the vertical error bar covers the statistical and systematic uncertainty.

Comparison with theory
For prompt production, as shown in Figure 11, the ratio of the NLO NRQCD theory calculations [61] to data, as a function of p T and in slices of rapidity, is provided for J/ψ and ψ(2S) at both the 7 and 8 TeV centre-of-mass energies. The theory predictions are based on the long-distance matrix elements (LDMEs) from Refs. [61,62], with uncertainties originating from the choice of scale, charm quark mass and LDMEs (see Refs. [61,62] for more details). Figure 11 shows fair agreement between the theoretical calculation and the data points for the whole p T range. The ratio of theory to data does not depend on rapidity.
For non-prompt ψ production, comparisons are made to FONLL theoretical predictions [1, 2], which describe the production of b-hadrons followed by their decay into ψ + X. Figure 12 shows the ratios of J/ψ and ψ(2S) FONLL predictions to data, as a function of p T and in slices of rapidity, for centre-of-mass energies of 7 and 8 TeV. For J/ψ, agreement is generally good, but the theory predicts slightly harder p T spectra than observed in the data. For ψ(2S), the shapes of data and theory appear to be in satisfactory agreement, but the theory predicts higher yields than in the data. There is no observed dependence on rapidity in the comparisons between theory and data for non-prompt J/ψ and ψ(2S) production.

Comparison of cross-sections 8 TeV with 7 TeV
It is interesting to compare the cross-section results between the two centre-of-mass energies, both for data and the theoretical predictions. Figure 13 shows the 8 TeV to 7 TeV cross-section ratios of prompt and non-prompt J/ψ and ψ(2S) for both data sets. For the theoretical ratios the uncertainties are neglected here, since the high correlation between them results in large cancellations.
Due to a finer granularity in p T for the 8 TeV data, a weighted average of the 8 TeV results is taken across equivalent intervals of the 7 TeV data to enable direct comparisons. Both data and theoretical predictions agree that the ratios become larger with increasing p T , however at the lower edge of the p T range the data tends to be slightly below theory.  Figure 12: The ratio of the FONLL theoretical predictions to data are presented for the differential non-prompt cross-section of J/ψ (left) and ψ(2S) (right) as a function of p T (µµ) for each rapidity slice. The top (bottom) row shows the 7 TeV (8 TeV) results. The error on the data is the relative error of each data point, while the error bars on the theory prediction are the relative error of each theory point.   Figure 13: The ratio of the 8 TeV and 7 TeV differential cross-sections are presented for prompt (top) and nonprompt (bottom) J/ψ (left) and ψ(2S) (right) for both data (red points with error bars) and theoretical predictions (green points). The theoretical predictions used are NRQCD for prompt and FONLL for non-prompt production. The uncertainty on the data ratio does not account for possible correlations between 7 and 8 TeV data, and no uncertainty is shown for the ratio of theory predictions.

Summary and conclusions
The prompt and non-prompt production cross-sections, the non-prompt production fraction of the J/ψ and ψ(2S) decaying into two muons, the ratio of prompt ψ(2S) to prompt J/ψ production, and the ratio of nonprompt ψ(2S) to non-prompt J/ψ production were measured in the rapidity range |y| < 2.0 for transverse momenta between 8 and 110 GeV. This measurement was carried out using 2.1 fb −1 (11.4 fb −1 ) of pp collision data at a centre-of-mass energy of 7 TeV (8 TeV) recorded by the ATLAS experiment at the LHC. It is the latest in a series of related measurements of the production of charmonium states made by ATLAS. In line with previous measurements, the central values were obtained assuming isotropic ψ → µµ decays. Correction factors for these cross-sections, computed for a number of extreme spin-alignment scenarios, are between −35% and +100% at the lowest transverse momenta studied, and between −14% and +9% at the highest transverse momenta, depending on the specific scenario.
The ATLAS measurements presented here extend the range of existing measurements to higher transverse momenta, and to a higher collision energy of √ s = 8 TeV, and, in overlapping phase-space regions, are consistent with previous measurements made by ATLAS and other LHC experiments. For the prompt production mechanism, the predictions from the NRQCD model, which includes colour-octet contributions with various matrix elements tuned to earlier collider data, are found to be in good agreement with the observed data points. For the non-prompt production, the fixed-order next-to-leading-logarithm calculations reproduce the data reasonably well, with a slight overestimation of the differential cross-sections at the highest transverse momenta reached in this analysis.

A. Spin-alignment correction factors
The measurement presented here assumes an unpolarized spin-alignment hypothesis for determining the correction factor. In principle, the polarization may be non-zero and may vary with p T . In order to correct these measurements when well-measured J/ψ and ψ(2S) polarizations are determined, a set of correction factors are provided in Tables 4-15 for the 7 TeV data, and in the Tables 16-27 for the 8 TeV data. These tables are created by altering the spin-alignment hypothesis for either the J/ψ or ψ(2S) meson and then determining the ratio of the mean sum-of-weights of the new hypotheses to the original flat hypothesis. The mean weight is calculated from all the events in each dimuon p T and rapidity analysis bin, selecting those dimuons within ±2σ of the ψ fitted mean mass position. The choice of spin-alignment hypothesis for each ψ meson has negligible effect on the results of the other ψ meson, and therefore these possible permutations are not considered. The definitions of each of the spin-alignment scenarios, which are given in the caption to the table, are defined in Table 1.      Table 9: Mean weight correction factor for J/ψ under the "off-(λ θ -λ φ )-plane negative" spin-alignment hypothesis for 7 TeV.            Table 20: Mean weight correction factor for J/ψ under the "off-(λ θ -λ φ )-plane positive" spin-alignment hypothesis for 8 TeV. Those intervals not measured in the analysis at low p T , high rapidity are also excluded here.  Table 21: Mean weight correction factor for J/ψ under the "off-(λ θ -λ φ )-plane negative" spin-alignment hypothesis for 8 TeV. Those intervals not measured in the analysis at low p T , high rapidity are also excluded here.      Table 26: Mean weight correction factor for ψ(2S) under the "off-(λ θ -λ φ )-plane positive" spin-alignment hypothesis for 8 TeV. Those intervals not measured in the analysis at low p T , high rapidity are also excluded here.  Table 27: Mean weight correction factor for ψ(2S) under the "off-(λ θ -λ φ )-plane negative" spin-alignment hypothesis for 8 TeV. Those intervals not measured in the analysis at low p T , high rapidity are also excluded here.