Suppression of Upsilon Production in d+Au and Au+Au Collisions at sqrt(s_NN) = 200 GeV

We report measurements of Upsilon meson production in p+p, d+Au, and Au+Au collisions using the STAR detector at RHIC. We compare the Upsilon yield to the measured cross section in p+p collisions in order to quantify any modifications of the yield in cold nuclear matter using d+Au data and in hot nuclear matter using Au+Au data separated into three centrality classes. Our p+p measurement is based on three times the statistics of our previous result. We obtain a nuclear modification factor for Upsilon(1S+2S+3S) in the rapidity range |y|<1 in d+Au collisions of R_dAu = 0.79 +/- 0.24 (stat.) +/- 0.03 (sys.) +/- 0.10 (pp sys.). A comparison with models including shadowing and initial state parton energy loss indicates the presence of additional cold-nuclear matter suppression. Similarly, in the top 10% most-central Au+Au collisions, we measure a nuclear modification factor of R_AA=0.49 +/- 0.1 (stat.) +/- 0.02 (sys.) +/- 0.06 (pp sys.), which is a larger suppression factor than that seen in cold nuclear matter. Our results are consistent with complete suppression of excited-state Upsilon mesons in Au+Au collisions. The additional suppression in Au+Au is consistent with the level expected in model calculations that include the presence of a hot, deconfined Quark-Gluon Plasma. However, understanding the suppression seen in d+Au is still needed before any definitive statements about the nature of the suppression in Au+Au can be made.

for Υ(1S+2S+3S) in the rapidity range |y| < 1 in d+Au collisions of R dAu = 0.79 ± 0.24(stat.) ± 0.03(sys.) ± 0.10(p+p sys.). A comparison with models including shadowing and initial state parton energy loss indicates the presence of additional cold-nuclear matter suppression. Similarly, in the top 10% most-central Au+Au collisions, we measure a nuclear modification factor of RAA = 0.49 ± 0.1(stat.) ± 0.02(sys.) ± 0.06(p+p sys.), which is a larger suppression factor than that seen in cold nuclear matter. Our results are consistent with complete suppression of excited-state Υ mesons in Au+Au collisions. The additional suppression in Au+Au is consistent with the level expected in model calculations that include the presence of a hot, deconfined Quark-Gluon Plasma. However, understanding the suppression seen in d+Au is still needed before any definitive statements about the nature of the suppression in Au+Au can be made.

INTRODUCTION
In the study of the properties of the Quark-Gluon Plasma (QGP) an extensive effort has been devoted to measuring quarkonium yields since these have been predicted to be sensitive to color deconfinement [1]. Studies have mainly focused on charmonium, but with the high collision energies available at the Relativistic Heavy Ion Collider (RHIC) and the Large Hadron Collider (LHC) we can now study bottomonium in hot nuclear matter with sufficient statistics. For a recent review of quarkonium in-medium, see e.g. Ref. [2], Sec. 5. One prediction is that excited quarkonium states are expected to dissociate at or above temperatures near that of the crossover to the deconfined QGP phase, T c ≈ 150 − 190 MeV [3][4][5]. The more tightly bound ground states are expected to dissociate at even higher temperatures. The details of the temperature dependence of the dissociation of the excited states and of the feed-down pattern of the excited states into the ground state lead to a sequential suppression pattern of the inclusive upsilon states with increasing temperature [6]. The binding energy of the Υ(2S) state (∼540 MeV) is about half that of the Υ(1S) state (∼1.1 GeV); the Υ(3S) is still more weakly bound at ∼200 MeV. Recent studies take into account not only the Debye screening effect on the heavy quark potential but also an imaginary part of the potential which modifies the widths of the various quarkonia states (e.g. [7][8][9][10]). In Ref. [8] it is estimated that the Υ(2S) state will melt at a temperature of T ≈ 250 MeV, whereas the ground state Υ(1S) will melt at temperatures near T ≈ 450 MeV.
We focus here on the measurement of bottomonium mesons in collisions at √ s N N =200 GeV. An observation of suppression in the bottomonium sector in hot nuclear matter is important for two reasons. First, it would be evidence for color deconfinement in the produced matter since the aforementioned effects are all ultimately based on studies of the high temperature phase of Quantum Chromodynamics (QCD) done on the lattice, where color is an active degree of freedom. Second, bottomonium suppression provides a way to estimate model-dependent bounds on the temperature with the bounds depending on the particular suppression pattern seen.
The cross section for bottomonium production is smaller than that of charmonium [11][12][13] making the experimental study of Υ production challenging. However, the theoretical interpretation of bottomonium suppression is less complicated than that of charmonium for several reasons. While charmonium production at RHIC and higher energies can be affected by the statistical recombination of charm quarks that are produced in different nucleon-nucleon collisions within the same nuclear interaction event, this effect is negligible for bottomonium due to the much smaller bb production cross section (σ bb is measured to be in the range 1.34 -1.84 µb [14] and calculated to be 1.87 +0.99 −0.67 µb [15], compared to σ cc ≈ 550 -1400 µb [16,17]). Another complication in the charmonium case is that even in a purely hadronic scenario, charmonium mesons can be suppressed due to their interaction with hadronic co-movers [18,19]. The cross section for inelastic collisions of Υ(1S) with hadrons is small [20]. Hence, absorption in the medium by the abundantly produced co-moving hadrons is predicted to be minimal. The cold-nuclear-matter (CNM) effects on Υ production, which are those seen in p+A collisions and can be due to shadowing of the parton distribution functions in the nucleus or energy-loss in the nucleus, can still be important. There is evidence of some Υ suppression in fixed target experiments at 800 GeV/c lab momentum from E772 [21]. However, the CNM suppression observed for Υ is smaller than that for J/ψ reported by NA50 [22]. For all these reasons, the Υ family is expected to be a cleaner and more direct probe of hot QCD, and of the corresponding color deconfinement effects.
In this letter we present measurements of Υ production in p+p, d+Au, and Au+Au collisions at √ s N N =200 GeV via the e + e − decay channel obtained by the STAR experiment. We extract invariant cross sections for all three collision systems studied. Using this p+p measurement as a baseline we obtain the nuclear modification factor (R dAu and R AA ) of the three states combined: Υ(1S+2S+3S). The ratio R AA is used to quantify deviations of the yields in d+Au and Au+Au compared to those expected from a superposition of elementary p+p collisions. The data were taken during 2008 (d+Au), 2009 (p+p) and 2010 (Au+Au) at RHIC, and correspond to integrated luminosities of 28.2 nb −1 , 20.0 pb −1 , and 1.08 nb −1 , respectively. All three datasets were taken with the same detector configuration. For this reason the data from our previous p+p result (2006) was not included in this analysis; the amount of material in the detector at that point was substantially larger than it was in the three datasets discussed here. We compare our data to model calculations of the cross section based on perturbative QCD (pQCD) [23], and to recent models of Υ production in d+Au and Au+Au collisions [24][25][26][27][28][29].

EXPERIMENTAL METHODS
The main detectors used are the STAR Time Projection Chamber (TPC) [30] for tracking and the STAR Barrel Electro-Magnetic Calorimeter (BEMC) [31] for triggering. Both the TPC and BEMC are used for particle identification. The starting point is the STAR Υ trigger whose main components are a fast hardware Level-0 (L0) trigger, which fires when a tower in the BEMC has energy E L0−BEM C ≥ 4.2 GeV, and a software Level-2 (L2) trigger, which requires the presence of two high-energy clusters in the BEMC (>4.5 GeV and >3.0 GeV). The cluster pair must also have an opening angle greater than 90 • and an invariant mass above 5 GeV/c 2 (6.5 GeV/c 2 ) in p+p (d+Au). Note that energy measured at the triggering level is partially calibrated leading to small but random biases. Hence, triggering thresholds are not precise in energy. The Υ trigger is required to be in coincidence with the STAR minimum bias trigger. For p+p collisions the minimum bias trigger is based on the STAR Beam-Beam Counters, while for d+Au and Au+Au it is based on the STAR Zero-Degree Calorimeters (ZDC) and the Vertex-Position Detectors (VPD). The L0-L2 combination was used for the d+Au data in 2008 and for p+p data in 2009. In the Au+Au 2010 run, an upgrade to the STAR data acquisition system allowed the processing of all the L0 triggers above the E L0−BEM C = 4.2 GeV threshold, thus removing the need for a Level-2 trigger.
Some of the key components common to all these analyses are the tracking, matching between TPC tracks and BEMC L0 and L2 clusters, and electron identification techniques. The main differences between the three datasets are summarized as follows. For Au+Au collisions we use the charged particle multiplicity measured in the TPC in order to determine the centrality of the collision. Using a Glauber model simulation, the multiplicity classes in the collision are used to estimate the average number of participants (N part ) and number of binary collisions (N coll ). The trigger, tracking, and electron identification efficiencies in the Au+Au case were studied as a function of centrality (see Tab. I). The presence of the underlying Au+Au event background increases the energy measured in the calorimeter towers and results in a slight increase in the trigger efficiency with increasing N part (more central collisions). Similarly, the increase in the track density in the TPC results in a decrease in the tracking efficiency which is especially noticeable at high N part . We used the specific ionization of the tracks in the TPC gas (dE/dx) for electron identification. In addition, the projection of the track onto the location of the BEMC shower maximum position was required to match the measured BEMC cluster position. Once a track was matched to a calorimeter cluster, the ratio of the energy of the cluster to the TPC momentum (E/p) was also used for electron identification. The combined acceptance times efficiency for detecting an Υ at midrapidity (|y| < 0.5) in Au+Au taking into account all aforementioned effects was found to vary from ∼ 12% in peripheral collisions to ∼ 9% in central collisions.
The cuts used in these analyses were chosen such that the tracking and electron identification efficiencies would be similar across the three datasets, allowing the systematic uncertainties to approximately cancel in the measurement of R AA . For further detail, the techniques used in these Υ measurements were described extensively for our previous p+p measurement [13] based on a 7.9 pb −1 dataset. All evaluated sources of systematic uncertainty are summarized in Tab. II. An important effect in addition to those discussed in [13] is the change in tracking and mass resolution with increasing detector occupancy. Simulated Υ events were embedded in real data and their reconstructed line shapes were studied as a function of collision system and detector occupancy. In p+p collisions, we find a mass resolution of 1.3% for reconstructed Υ(1S). Due to additional TPC alignment errors for the d+Au and Au+Au data the mass resolutions of the Υ(1S) increased to 2.7% in d+Au and peripheral Au+Au collisions and 2.9% in central collisions. This decreased mass resolution was accounted for in the binary scaling estimates of Υ(1S+2S+3S) yields (see gray bands in Figs. 1 and 4). Systematic uncertainties in those scaling estimates (line shapes) are included in the errors in Table II.

Source
Relative uncertainty Luminosity, Vernier scan (p+p) For all results we quote, the Υ data are integrated over all transverse momenta.
RESULTS AND DISCUSSION Figure 1 shows the invariant mass distributions of electron pairs for p+p (top) and d+Au (bottom) in the kinematic region |y Υ | < 0.5. Unlike-sign pairs are shown as red filled circles and like-sign pairs as hollow blue circles.
The data are fit with a parameterization consisting of the sum of various contributions to the electron-pair invariant-mass spectrum. The fit is performed simultaneously with the like-sign and unlike-sign spectra using a maximum-likelihood method. The lines in Fig. 1 show the yield from the combinatorial background (dashed blue line), the result of adding the physics background from Drell-Yan and bb pairs (dot-dashed green line), and finally the inclusion of the Υ contribution (solid red line). The shape of the Drell-Yan continuum is obtained via a next-to-leading order (NLO) pQCD calculation from Vogt [32]. PYTHIA 8 was used to calculate the shape of the bb contribution [33]. We model each of the Υ states with a Crystal Ball function [34], which incorporates detector resolution and losses from bremsstrahlung in the detector material.
The fit is done to the unlike and like-sign data simultaneously. The fit to the combinatorial background component extracted from the like-sign data is shared by the functional form used to parameterize the unlike-sign data. In the usual like-sign subtraction procedure some information would be lost. In contrast, by performing a simultaneous fit to both the like-sign and unlike-sign signals we optimize the statistical power of our data. The L2 trigger condition has the effect of cutting off the lower invariant masses. This cut-off shape is parameterized in the fits using an error function.
We integrate the unlike-sign invariant mass distribution in the region 8.8 − 11 GeV/c 2 and subtract from the data the fit to the combinatorial, Drell-Yan, and bb background components in order to obtain the yield of Υ(1S+2S+3S). After accounting for efficiencies and sampled luminosity, we calculate a production cross section in p+p collisions of: B ee × dσ/dy| |y|<0.5 = 64 ± 10(stat. + fit) +14 −12 pb. Our previous result of 114 ± 38 +23

−24
pb [13] is consistent with our new measurement. The greater sampled luminosity and decreased detector material in 2009 led to improved statistics and lower systematic uncertainties in the present measurement. In Fig. 1(b), the gray band shows the expected signal from the p+p data scaled by the number of binary collisions. Due to differences in detector occupancy and detector calibrations the width of the Υ signal differs between collision systems and centralities. As discussed in the previous section, a misalignment in the TPC in the d+Au and Au+Au datasets led to a broadening of the Υ line shapes compared to the p+p dataset. This can be seen by examining the line shapes for the Υ states in Fig. 1(a) (p+p) and Fig. 1(b) (d+Au). The average detector occupancy is comparable between the two systems, however the d+Au dataset has a noticeably broader line shape due to the aforementioned differences in calibration. The effects of the broadening of the line shapes are taken into account in systematic uncertainties (Tab. II). The comparison of the gray band with the d+Au data in panel (b) indicates a suppression of Υ production with respect to binary-collision scaling.
A similar procedure is followed for the region 0.5 < |y Υ | < 1 in p+p collisions. We combine the results to obtain the differential cross section: B ee × dσ/dy| |y|<1 = 58 ± 12(stat. + fit) +12 −11 pb. In d+Au collisions, we analyze the yields separately in the regions −1 < y Υ < −0.5 and 0.5 < y Υ < 1 because the d+Au system is not symmetric about y = 0. Hence, averaging between forward and reverse rapidities is not warranted as it is in p+p. Throughout this paper, the positive rapidity region is the deuteron-going direction, and the negative rapidity region is the Au-going direction. Integrating over our measured range (|y Υ | < 1), the cross section in d+Au collisions is found to be B ee × dσ/dy| |y|<1 = 19 ± 3(stat. + fit) ± 3(syst.) nb.
We extract the Υ(1S) yield directly by integrating over a narrower mass window (8.8-9.8 GeV/c 2 ). This mass window was chosen due to its high acceptance rate for Υ(1S) and its high rejection rate for the excited states. To account for sensitivity to the shape of the Υ signal, we varied the parameters of the line shape obtained from simulations and data-driven methods discussed previously by their measured uncertainties and varied the excited states from unsuppressed to completely suppressed. We then recalculated both efficiency and pu-rity (see Tab. II, Υ(1S) purity). Those variations were taken into account as additional systematics when quoting Υ(1S) results. Figure 2(a) shows the extracted Υ(1S+2S+3S) cross section for p+p and d+Au as a function of rapidity. The p+p measurements are shown as blue stars and the d+Au measurements as red circles. The p+p result in the region 0.5 < |y Υ | < 1.0 is displayed as a star at y = 0.75 and also as a hollow star at y = −0.75 to illustrate that the latter is not an independent measurement. The data from PHENIX at forward rapidity for p+p (filled blue diamonds) and d+Au (hollow red diamonds) are also shown [35].
The cross sections in p+p are compared to an NLO pQCD calculation of Υ production in the Color Evaporation Model (CEM) [23], which is consistent with our data within the statistical and systematic uncertainties. The same calculation is performed for d+Au including shadowing effects [24]. The EPS09 set of nuclear Parton Distribution Functions (nPDF) [36] were used. The model is in agreement with our data except for the midrapidity point which is lower than the prediction. To study this observation for d+Au further, we make a closer comparison to models and to previous measurements of Υ production in p+A collisions.
To focus on expected shadowing effects, we obtain the nuclear modification factor R dAu as a function of rapidity. The nuclear modification factor is defined in nucleusnucleus collisions as where the first factor accounts for the difference in inelastic cross section in p+p to d+Au or Au+Au collisions. The second factor accounts for the average number of nucleon collisions in a d+Au or Au+Au collision as calculated by a Glauber model. The third factor accounts for the measured Υ production in p+p, d+Au or Au+Au. We used the following total inelastic cross sections: σ pp = 42 mb, σ dAu = 2.2 b, and σ AuAu = 6 b.
Our results for R dAu are shown in Fig. 2(b) and summarized in Tab. III. Our data (red stars) are compared to CEM calculations with the uncertainty from the EPS09 nPDF shown as the shaded region. Note that this prediction for R dAu , which includes modification of the nuclear PDFs but does not include absorption, implies a modification factor of R dAu ≈ 1.1. A calculation in Ref. [25] explored various nPDFs (EKS98, EPS08, and nDSg) and also gave R dAu values above 1 with enhancements in the range of 5-20%. The models are in agreement with the data except in the y ∼ 0 region. An additional effect which can suppress the Υ yield is initial-state parton energy loss. A calculation by Arleo and Peigné [26] incorporating this effect is shown as the dashed line. The  calculation for a combination of energy loss and shadowing using EPS09 is shown as the dashed-dotted line.
The energy-loss model is also in agreement with the data except for the mid-rapidity point. The model from [26] does not include absorption from interactions with spectator nucleons. However, those effects only play a role in the rapidity region y 1.2, where the Υ mesons would be closer to the frame of the Au spectators. Therefore, the suppression at mid-rapidity is indicative of effects beyond shadowing, initial-state parton energy loss, or absorption by spectator nucleons.
We compare our measurements with results from E772 at √ s N N = 40 GeV, where suppression of the Υ states in p+A was observed. This is illustrated in Fig. 3(a), which shows the ratio of the cross section in d+Au collisions for STAR (p+A for E772) to that of p+p collisions normalized by the mass number A. E772 plotted a ratio of extracted cross sections normalized by the data where the proton beam hit a liquid deuterium target (A = 2). Assuming that the cross section scales as σ pA = A α σ pp , and using their p+d result as the baseline, the solid line shows that the ratio should scale as (A/2) α−1 . Our measurement in d+Au for the Υ(1S) state (red star) is consistent with the fit to the E772 data, shown as hollow blue circles for Υ(1S) and hollow green squares for Υ(2S+3S). Our results cover the rapidity range |y| < 1 whereas the E772 measurements were in the forward region 0 < y < 1.05. To better compare our rapidity coverage, we plot the α value as a function of Feynman-x (x F ) in Fig. 3(b). The larger suppression we observe at mid-rapidity is also consistent with the larger suppression seen in E772 for x F ∼ 0.
We next turn to the measurements in Au+Au collisions. The Au+Au invariant mass spectrum is fit in 3 centrality bins: 30-60% ( Fig. 4(a), 10-30% ( Fig. 4(b), and 0-10% (Fig. 4(c). As in Fig. 1 we show the fits including, in succession, combinatorial background (dashed blue line), the physics background from Drell-Yan and bb pairs (dot-dashed green line), and the Υ contribution (solid red line). The absence of the L2 trigger in the Au+Au dataset removes the cut-off effect. One can therefore see the background (modeled as the sum of two exponentials), dominated by the combinatorial component, rising at lower invariant mass. Measured cross sections are summarized in Tab This suppression is quantified in Fig. 5, which displays the nuclear modification factor, R AA , plotted as a function of N part with the 0-10% most-central collisions corresponding to N part = 326 ± 4. Figure 5(a) shows the data for all three states in the rapidity range |y| < 1, while Fig. 5(b) is for the narrower |y| < 0.5 range. . Note that the hollow star at y = −0.75 is a reflection of the filled one at y = 0.75 since these are not independent measurements. Results obtained by PHENIX are shown as filled diamonds. Systematic errors are shown as boxes around the data. The shaded bands are from nextto-leading-order pQCD color evaporation model calculations. The d+Au prediction uses the EPS09 nPDF which includes shadowing [24]. (b) R dAu vs. y for STAR (red stars) and PHENIX (green diamonds) results. The band on the right shows the overall normalization uncertainty for the STAR results due to systematic uncertainties in the p+p measurement. The shaded band shows the prediction for R dAu from EPS09 and its uncertainty. The dashed curve shows suppression due to initial-state parton energy loss and the dot-dashed curve shows the same model with EPS09 incorporated [26].   IV: Υ production cross sections in Au+Au collisions. The first uncertainty listed is the combination of the statistical and fit uncertainties and the second is the systematic uncertainty.
In the narrower rapidity range (Fig. 5(b)), we see an indication of a lower R dAu as discussed earlier. Our data and the E772 data show a larger suppression at y ∼ 0 or x F ∼ 0 than that expected from shadowing. The level of suppression we observe for |y| < 0.5 stays approximately constant from d+Au up to central Au+Au collisions. This suggests that suppression in d+Au in this kinematic range needs to be understood before interpreting the suppression in Au+Au.
For the 0-10% most-central collisions we find R AA (1S) = 0.66 ± 0.13 (Au+Au stat.) ± 0.10 (p+p stat.) +0.02 −0.05 (Au+Au syst.) ± 0.08 (p+p syst.). Similar suppression is found by CMS in Pb+Pb collisions (R AA (1S)≈ 0.45 at similar N part ) [37][38][39]. We observe the nuclear modification factor for the Υ(1S) as a function of N part to be consistent with unity in d+Au through mid-central Au+Au collisions (see Fig. 5c). In the most central Au+Au collisions, we see an indication of suppression of the Υ(1S) at the 2.7σ level. In the context of suppression of the excited states, if the feed-down fraction remains ∼ 49% as measured at higher energies and high-p ⊥ it is possible that an R AA (1S) as low as 0.51 could be due solely to suppression of the excited states [43].
One can relate the R AA of the combined states to that of the ground state via the equation R AA (1S+2S+3S) = R AA (1S)×(1 + N AA (2S + 3S)/N AA (1S))/(1 + N pp (2S + 3S)/N pp (1S)). The ratio of the excited states to the ground state can be obtained from measurements by CMS and Fermilab experiments [40,41] and alternatively from combining theoretical calculations [23] with measured branching ratios from the PDG [42]. In the case where N AA (2S + 3S) = 0, R AA (1S+2S+3S) ≈ R AA (1S)×0.7. This is consistent with our observed R AA values, and can also be inferred by examining the mass range 10-11 GeV/c 2 in Fig. 4, where no significant 2S or 3S signals are seen.
By applying the methods described in [44], we can cal- culate an upper limit on the R AA of the combined 2S and 3S states. Using the fit to Drell-Yan and bb (dashed, green curve) as the background, we find an upper limit of 29 signal counts with 95% confidence in the mass range 10-11 GeV/c 2 for 0-60% centrality collisions. To transform this upper limit into an upper limit on R AA (2S+3S), we assumed that the purity of excited states in this mass range is the same as in the p+p case. While the excited states are likely more suppressed than the ground state in the Au+Au case, using the p+p purity gives us an upper limit in the Au+Au purity which can be used to calculate an upper limit on the R AA . The 2S+3S cross section in p+p was extracted from the full cross section, assuming the purity can be obtained based on the PDG branching ratios [42] and the relative production cross sections of the three states. In the centrality range of 0-60%, we thus obtain a 95%-confidence upper limit of R AA (2S+3S) < 0.32 (see Fig. 6).
Our data are also compared to model calculations incorporating hot-nuclear-matter effects for Au+Au [27][28][29]. These aim to incorporate lattice-QCD results pertinent to screening and broadening of bottomonium and to model the dynamical propagation of the Υ meson in the colored medium. Both models are in agreement with the level of suppression seen in Au+Au. The model proposed by Emerick, Zhao, and Rapp (EZR), Ref. [28], includes possible CNM effects, modeled as an absorption cross section of up to 3 mb which can account for a value of R AA as low as 0.7. In this model the additional suppression to bring R AA down to ≈ 0.5 is due to hot-nuclear-matter effects. The calculation by Liu et al. [29] in Fig. 5(c) is for the inclusive Υ(1S) R AA , using the internal energy as the heavy-quark potential and an initial temperature of the fireball of T = 340 MeV, which given the input from lattice QCD results, is not hot enough to melt the directly produced Υ(1S). Hence, the suppression is mostly driven by the dissociation of the excited states (both the S-states and the P-states). The initial temperature used in the EZR model is 330 MeV (with a formation time of 0.6 fm/c). The temperatures of the QGP needed in Strickland's model, Ref. [27], are in the range 428 -442 MeV. However, it should be noted that neither the Strickland model, nor the calculation from Liu et al. include any CNM effects.
Considering two possible sources of suppression, CNM and QGP effects, we used a Monte Carlo pseudoexperiment to compare our results to different possible sources of suppression. We investigated four possible scenarios: (1) No suppression compared to p+p; (2) Suppression due to CNM effects only; (3) QGP suppression only; (4) Suppression from both CNM and QGP effects. We simulated Υ production in p+p, d+Au, and Au+Au collisions via a Poisson generator. CNM effects were included via the suppression parametrization used by E772 [21] and presented in Fig. 3(a). We used the predictions from the Strickland model [27] to estimate suppression from QGP effects. For scenario (4), the expected suppression is simply taken to be the product of the suppression from scenario (2) and scenario (3). For this pseudoexperiment we assumed a flat prior based on the allowed R AA given in Strickland-Bazow [27], depicted as the band for this calculation in Fig. 5, stemming from the choice of 1 < 4πη/S < 3.
A summary of the pseudoexperiment results is shown in Fig. 7. Panel (a) shows our result for R dAu in the range |y| < 1.0 compared to scenarios (1) and (2), shown as the solid line and dotted histogram, respectively. The 'nocold-suppression-scenarios' (1 and 3) are excluded while the CNM effects from E772 parameterization are consistent with our observation. Panel (b) shows R AA for the most-central Au+Au bin in the range |y| < 1.0. By comparing the results of the pseudoexperiments with our measurements, we are able to exclude scenario (1) at a ∼ 5σ confidence level. Finally, we see that hypothesis (4) (dot-dashed curve), including both hot and cold nuclear effects, is consistent with our measurements when both the d+Au and Au+Au results are taken into account.
We repeated this procedure for the rapidity range |y| < 0.5. The results are shown in Figs. 7 (c) and (d). In the mid-rapidity range we find a larger amount of suppression in d+Au than what we observe in the range |y| < 1.0. Furthermore, R dAu is comparable to R AA in 0-10% for this rapidity range. This could indicate that suppression of bottomonium already occurs in d+Au collisions. However, given the uncertainties in our current results, no particular model of Υ suppression in d+Au is favored. Hence, further investigation of cold-nuclearmatter effects on Υ production is highly warranted. The suppression effects seen in d+Au, which are not explained . We show our results for two systems and two rapidity ranges: (a) d+Au |y| < 1.0, (b) Au+Au |y| < 1.0, (c) d+Au |y| < 0.5, (d) Au+Au |y| < 0.5. Our data is shown as a red vertical line with systematics shown by the pink box. The QGP effects are modeled in [27].
by the models discussed here, still need to be understood before the Au+Au results can be fully interpreted.

CONCLUSIONS
In conclusion we studied Υ(1S+2S+3S) production in p+p, d+Au, and Au+Au collisions at √ s=200 GeV. We measured the cross section in p+p collisions to be B ee × dσ/dy| |y|<1 = 61 ± 8(stat. + fit) +13 −12 (syst.) pb and find it to be consistent within errors with NLO calculations. The cross section in d+Au collisions is found to be B ee × dσ/dy| |y|<1 = 19 ± 3(stat. + fit) ± 3(syst.) nb. We obtain a nuclear modification factor in this rapidity region (|y| < 1) of R dAu (1S + 2S + 3S) = 0.79 ± 0.22 (d+Au stat.)±0.10 (p+p stat.)±0.03 (d+Au syst.)± 0.09(p+p syst.). Models of Υ production in cold nuclear matter, which include shadowing and initial-state partonic energy loss, are consistent with the cross-sections we observe in our d+Au data. Higher statistics d+Au data are required to further investigate the 3σ deviation we observe at |y| < 0.5. We measured the Υ(1S+2S+3S) nuclear modification factor in Au+Au collisions at √ s N N =200 GeV as a function of centrality. In the range |y| < 1 and in 0-10% most-central collisions we find R AA (1S + 2S + 3S) = 0.49 ± 0.13 (Au+Au stat.) ± 0.07 (p+p stat.)±0.02 (Au+Au syst.)±0.06 (p+p syst.), indicating additional Υ suppression in hot nuclear matter compared to cold nuclear matter. In 0-60% centrality we find a 95%-confidence upper limit on the nuclear modification of the excited states of R AA (2S+3S)< 0.32. Calculations of the centrality dependence of Υ R AA using models based on lattice QCD calculations of bottomonium melting in a hot medium are found to be consistent with our data. Therefore, the suppression seen in central Au+Au collisions is indicative of the presence of deconfined matter in heavy-ion collisions. It would be desirable to have a higher statistics d+Au dataset in order to strengthen the conclusions regarding cold-nuclear modifications to Υ production before a stronger connection between parton deconfinement, Debye screening, and the observed Υ suppression in Au+Au can be made.