Rivet analyses referenceALICE_2019_I1672860Rho meson production in pp and Pb-Pb at 2.76 TeVExperiment: ALICE (LHC) Inspire ID: 1672860 Status: VALIDATED Authors:
Beam energies: (1380.0, 1380.0); (287040.0, 287040.0) GeV
The production of the $\rho(770)0$ meson has been measured at midrapidity ($|y|<0.5$) in pp and centrality differential Pb-Pb collisions at $\sqrt{s_\text{NN}}= 2.76$ TeV with the ALICE detector at the Large Hadron Collider. The particles have been reconstructed in the $\rho(770)0\to\pi^{++}\pi^{-}$ decay channel in the transverse-momentum ($p_\text{T}$) range 0.511 GeV/c. A centrality-dependent suppression of the ratio of the integrated yields $2\rho(770)0/(\pi^{++}\pi^{-})$ is observed. The ratio decreases by $\sim40\%$ from pp to central Pb-Pb collisions. A study of the $p_\mathrm{T}$-differential $2\rho(770)0/(\pi^{++}\pi^{-})$ ratio reveals that the suppression occurs at low transverse momenta, $p_\mathrm{T}<2$ GeV/c. At higher momentum, particle ratios measured in heavy-ion and pp collisions are consistent. The observed suppression is very similar to that previously measured for the $K^\ast(892)0/K$ ratio and is consistent with EPOS3 predictions that may imply that rescattering in the hadronic phase is a dominant mechanism for the observed suppression. Source code: ALICE_2019_I1672860.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/Beam.hh"
4#include "Rivet/Projections/FinalState.hh"
5#include "Rivet/Projections/UnstableParticles.hh"
6#include "Rivet/Tools/Cuts.hh"
7#include "Rivet/Projections/SingleValueProjection.hh"
8#include "Rivet/Tools/AliceCommon.hh"
9#include "Rivet/Projections/AliceCommon.hh"
10#include "Rivet/Projections/HepMCHeavyIon.hh"
11
12namespace Rivet {
13
14
15 /// @brief Rho meson production in pp and Pb-Pb at 2.76 TeV
16 class ALICE_2019_I1672860 : public Analysis {
17
18 public:
19
20 /// Constructor
21 RIVET_DEFAULT_ANALYSIS_CTOR(ALICE_2019_I1672860);
22
23 void init() {
24
25 // Find out the beam type
26 const ParticlePair& beam = beams();
27 if (beam.first.pid() == PID::PROTON && beam.second.pid() == PID::PROTON) isHI = false;
28 else if (beam.first.pid() == PID::LEAD && beam.second.pid() == PID::LEAD) {
29 isHI = true;
30 }
31 else {
32 MSG_WARNING("No beam found. You are likely in REENTRANT status.");
33 isHI = true;
34 }
35
36 if (isHI) {
37 declare(HepMCHeavyIon(), "HepMC");
38
39 declareCentrality(ALICE::V0MMultiplicity(), "ALICE_2015_CENT_PBPB", "V0M", "V0M");
40
41 _centrality_regions.clear();
42 _centrality_regions = {{0., 20.}, {20., 40.}, {40., 60.}, {60., 80.}};
43 }
44
45
46
47 // Charged, primary particles with |eta| < 0.5
48 declare(ALICE::PrimaryParticles(Cuts::absrap < 0.5 && Cuts::abscharge > 0), "APRIM");
49
50 // Resonances
51 declare(UnstableParticles(Cuts::absrap<0.5), "RSN");
52
53 // Booking histograms
54 // rho pt spectrum in pp (Table 6 of HEPData and Fig. 5 in Article)
55 book(_hist_rho_PP, 6, 1, 1);
56 book(_counterSOW_PP, "/TMP/counterSOW_PP");
57
58 // NOTE. pion pt spectrum is needed for the ratio and
59 // therefore it requires the same binning as the rho pt spectrum
60 std::string name_pion_PP_Fig5 = mkAxisCode(6, 1, 1) + "-pion-pp";
61 book(_hist_pion_PP, name_pion_PP_Fig5, refData(6, 1, 1));
62
63 // ratio to pion yield (Table 13, Fig. 9)
64 book(_rho_pion_ratio_PP, 13, 1, 1);
65
66
67 book(_counter_temp, "/TMP/counter.temp"); // counter for PbPb
68 book(_counterNcoll_temp, "/TMP/counter.ncoll.temp"); // Ncoll counter for PbPb
69
70
71
72 //----------------------------------------------------------------------------------
73 // Loop over all histograms
74 for (size_t ihist = 0; ihist < NHISTOS; ++ihist) {
75
76 const string nameCounterPbPb = "/TMP/counter.pbpb." + std::to_string(ihist);
77 book(_counterSOW[PBPB][ihist], nameCounterPbPb); // Sum of weights counter for PbPb
78
79 const string nameCounterNcoll = "/TMP/counter.ncoll." + std::to_string(ihist);
80 book(_counterNcoll[ihist], nameCounterNcoll); // Ncoll counter for PbPb
81
82 // rho pt spectra in PbPb (Tables 7-10 in HEPData and Fig. 6 in Article)
83 book(_hist_rho[PBPB][ihist], ihist+7, 1, 1);
84
85 // NOTE. pion pt spectra are needed for the ratios and
86 // therefore require the same binning as the rho pt spectra
87 const string name_pion_PbPb = "/TMP/"+mkAxisCode(ihist+7,1,1) + "-pion";
88 book(_hist_pion[PBPB][ihist], name_pion_PbPb, refData(ihist+7, 1, 1));
89
90 // NOTE. Only two out of four ratios (0-20% and 60-80%) are presented in the article
91 // HEPData (Tables 14, 15) and Article (Fig. 10)
92 const string name_rho_pion_ratio_PbPb = "/TMP/"+mkAxisCode(ihist+7,1,1) + "-rho_pion_ratio";
93 if (ihist==0) book(_rho_pion_ratio[PBPB][ihist], 14, 1, 1);
94 else if (ihist==3) book(_rho_pion_ratio[PBPB][ihist], 15, 1, 1);
95 else book(_rho_pion_ratio[PBPB][ihist], name_rho_pion_ratio_PbPb, refData(ihist+7, 1, 1).xEdges());
96
97 // Next lines are for RAA
98 // Initialize pp objects. In principle, only one pp histogram would be
99 // needed since centrality does not make any difference here. However,
100 // in some cases in this analysis the binning differ from each other,
101 // so this is easy-to-implement way to account for that.
102 const string nameCounterpp = "/TMP/counter.pp." + std::to_string(ihist);
103 book(_counterSOW[PP][ihist], nameCounterpp); // Sum of weights counter for pp
104 const string namePP = "/TMP/"+mkAxisCode(ihist+7,1,1) + "-pp";
105 book(_hist_rho[PP][ihist], namePP, refData(ihist+7, 1, 1));
106 const string name_pion_PP = "/TMP/"+mkAxisCode(ihist+7,1,1) + "-pion-pp";
107 book(_hist_pion[PP][ihist], name_pion_PP, refData(ihist+7, 1, 1));
108 const string name_rho_pion_ratio_PP = "/TMP/"+mkAxisCode(ihist+7,1,1) + "-rho_pion_ratio-pp";
109 book(_rho_pion_ratio[PP][ihist], name_rho_pion_ratio_PP, refData(ihist+7, 1, 1).xEdges());
110 // RAA are in Tables 16-19 in HEPData and Fig. 11 in Article
111 book(_hist_RAA[ihist], ihist+16, 1, 1);
112 } // end loop over histograms
113
114 // integrated yields, ratios and mean pt vs. multiplicity
115 // Table 11, Fig. 8, left
116 book(_hist_integrated_rho_pion_ratio, 11, 1, 1);
117 book(_hist_integrated_yield_rho, "/TMP/integrated_yield_rho", refData( 11, 1, 1));
118 book(_hist_integrated_yield_pion, "/TMP/integrated_yield_pion", refData( 11, 1, 1));
119 book(_hist_mean_pt_rho, 12, 1, 1);
120
121 } // end init
122
123 // Perform the per-event analysis
124 void analyze(const Event& event) {
125
126 if (int_edges.empty()) {
127 int_edges.push_back(_hist_integrated_yield_pion->bin(1).xMid());
128 int_edges.push_back(_hist_integrated_yield_pion->bin(3).xMid());
129 int_edges.push_back(_hist_integrated_yield_pion->bin(5).xMid());
130 int_edges.push_back(_hist_integrated_yield_pion->bin(7).xMid());
131 int_edges.push_back(_hist_integrated_yield_pion->bin(9).xMid());
132 }
133
134 // Charged, primary particles in eta range of |eta| < 0.5
135 Particles chargedParticles = apply<ALICE::PrimaryParticles>(event,"APRIM").particlesByPt();
136
137 // Resonances
138 const UnstableParticles &rsn = apply<UnstableParticles>(event, "RSN");
139
140 if (isHI) {
141
142 const HepMCHeavyIon &hi = apply<HepMCHeavyIon>(event, "HepMC");
143 if (!hi.ok()) {
144 MSG_WARNING("HEPMC Heavy ion container needed for this analysis, but not "
145 "found for this event. Skipping.");
146 vetoEvent;
147 }
148
149 _counter_temp->fill(hi.Ncoll());
150 _counterNcoll_temp->fill(hi.Ncoll());
151
152
153 // Prepare centrality projection and value
154 const CentralityProjection& centrProj = apply<CentralityProjection>(event, "V0M");
155 double centr = centrProj();
156 // Veto event for too large centralities since those are not used
157 // in the analysis at all
158 if ((centr < 0.) || (centr > 80.)) vetoEvent;
159
160 for (size_t ihist = 0; ihist < NHISTOS; ++ihist) {
161
162 const double low_edge_rho = _hist_rho[PBPB][ihist]->xMin();
163 const double high_edge_rho = _hist_rho[PBPB][ihist]->xMax();
164 const double low_edge_pion = _hist_pion[PBPB][ihist]->xMin();
165 const double high_edge_pion = _hist_pion[PBPB][ihist]->xMax();
166
167
168 if (inRange(centr, _centrality_regions[ihist].first, _centrality_regions[ihist].second)) {
169
170 _counterSOW[PBPB][ihist]->fill();
171 _counterNcoll[ihist]->fill(hi.Ncoll());
172
173 for (const Particle &p : rsn.particles()) {
174 if (p.abspid() != PID::RHO0) continue;
175 _hist_integrated_yield_rho->fill(int_edges[4-ihist]);
176 const double pT = p.pT()/GeV;
177 _hist_mean_pt_rho->fill(int_edges[4-ihist], pT);
178 if (pT > low_edge_rho && pT < high_edge_rho) {
179 _hist_rho[PBPB][ihist]->fill(pT);
180 } // condition on pT
181 } // end loop over resonances
182
183 //----------------------------------------------------------------------------------
184 for (const Particle& p : chargedParticles) {
185 if (p.abspid() != PID::PIPLUS) continue;
186 _hist_integrated_yield_pion->fill(int_edges[4-ihist]);
187 const double pT = p.pT()/GeV;
188 if (pT > low_edge_pion && pT < high_edge_pion) {
189 _hist_pion[PBPB][ihist]->fill(pT);
190 }
191 } // end loop over charged primary particles
192 } // centrality
193 } // histo loop
194 } // end PbPb event
195 else { // PP event
196
197 _counterSOW_PP->fill();
198
199 double low_edge_rho = _hist_rho_PP->xMin();
200 double high_edge_rho = _hist_rho_PP->xMax();
201
202 for (const Particle &p : rsn.particles()) {
203 if (p.abspid() != PID::RHO0) continue;
204 _hist_integrated_yield_rho->fill(int_edges[0]); // fill first bin for pp
205 const double pT = p.pT()/GeV;
206 _hist_mean_pt_rho->fill(int_edges[0], pT);
207 if (pT > low_edge_rho && pT < high_edge_rho) {
208 _hist_rho_PP->fill(pT);
209 } // condition on pT
210 } // end loop over resonances
211
212 double low_edge_pion = _hist_pion_PP->xMin();
213 double high_edge_pion = _hist_pion_PP->xMax();
214
215 for (const Particle& p : chargedParticles) {
216 if (p.abspid() != PID::PIPLUS) continue;
217 _hist_integrated_yield_pion->fill(int_edges[0]); // fill first bin for pp
218 const double pT = p.pT()/GeV;
219 if (pT > low_edge_pion && pT < high_edge_pion) {
220 _hist_pion_PP->fill(pT);
221 }
222 } // end loop over charged primary particles
223
224 // next histograms are needed only for the RAA
225 // NOTE. for completeness and consistency checks also histo for pions are provided
226 // although they are not needed (see above for rho/pion ratio in pp).
227 for (size_t ihist = 0; ihist < NHISTOS; ++ihist) {
228
229 _counterSOW[PP][ihist]->fill();
230
231 low_edge_rho = _hist_rho[PP][ihist]->xMin();
232 high_edge_rho = _hist_rho[PP][ihist]->xMax();
233
234 for (const Particle &p : rsn.particles()) {
235 if (p.abspid() != PID::RHO0) continue;
236 double pT = p.pT()/GeV;
237 if (pT > low_edge_rho && pT < high_edge_rho) {
238 _hist_rho[PP][ihist]->fill(pT);
239 } // condition on pT
240 } // end loop over resonances
241
242 low_edge_pion = _hist_pion[PP][ihist]->xMin();
243 high_edge_pion = _hist_pion[PP][ihist]->xMax();
244
245 for (const Particle& p : chargedParticles) {
246 if (p.abspid() != PID::PIPLUS) continue;
247 const double pT = p.pT()/GeV;
248 if (pT > low_edge_pion && pT < high_edge_pion) {
249 _hist_pion[PP][ihist]->fill(pT);
250 }
251 } // end loop over charged primary particles
252 } // loop over histos
253 } // end pp event
254 } // end analyze
255
256
257
258
259 // Normalise histograms etc., after the run
260 void finalize() {
261
262 if (_counterSOW_PP->sumW() > 0.) {
263 scale(_hist_rho_PP, 1. / _counterSOW_PP->sumW());
264 scale(_hist_pion_PP, 1. / _counterSOW_PP->sumW());
265
266 if (_hist_rho_PP->numEntries() > 0 && _hist_pion_PP->numEntries() > 0) {
267 divide(_hist_rho_PP, _hist_pion_PP, _rho_pion_ratio_PP);
268 scale(_rho_pion_ratio_PP, 2.);
269 }
270
271 }
272
273 // Scaling of the histograms with their individual weights.
274 for (size_t itype = 0; itype < EVENT_TYPES; ++itype ) {
275 for (size_t ihist = 0; ihist < NHISTOS; ++ihist) {
276
277 if (_counterSOW[itype][ihist]->sumW() > 0.) {
278 scale(_hist_rho[itype][ihist], 1./ _counterSOW[itype][ihist]->sumW());
279 scale(_hist_pion[itype][ihist], 1./ _counterSOW[itype][ihist]->sumW());
280
281 if (_hist_rho[itype][ihist]->numEntries() > 0 && _hist_pion[itype][ihist]->numEntries() > 0) {
282 divide(_hist_rho[itype][ihist], _hist_pion[itype][ihist], _rho_pion_ratio[itype][ihist]);
283 scale(_rho_pion_ratio[itype][ihist], 2.);
284 }
285 }
286
287 } // end histo loop
288 } // PP and PBPB
289
290
291 // Postprocessing for RAA
292 for (size_t ihist = 0; ihist < NHISTOS; ++ihist) {
293 // If there are entries in histograms for both beam types
294 if (_hist_rho[PP][ihist]->numEntries() > 0 && _hist_rho[PBPB][ihist]->numEntries() > 0) {
295 // Initialize and fill R_AA histograms
296 divide(_hist_rho[PBPB][ihist], _hist_rho[PP][ihist], _hist_RAA[ihist]);
297
298 // Scale by Ncoll. Unfortunately some generators does not provide
299 // Ncoll value (eg. JEWEL), so the following scaling will be done
300 // only if there are entries in the counters
301 double ncoll = _counterNcoll[ihist]->sumW();
302
303 //PHYSICAL REVIEW C 88, 044909 (2013)
304 if(ihist==0) ncoll=1210.85;
305 else if(ihist==1) ncoll=438.;
306 else if(ihist==2) ncoll=127.7;
307 else ncoll=26.7;
308
309 double sow = _counterSOW[PBPB][ihist]->sumW();
310 if (ncoll > 1e-6 && sow > 1e-6) {
311 scale(_hist_RAA[ihist], 1. / ncoll);
312 }
313
314 }
315 } // loop over histos
316
317
318 // Postprocessing for integrated yield vs. multiplicity
319 if (_hist_integrated_yield_rho->numEntries() > 0. && _hist_integrated_yield_pion->numEntries() > 0.) {
320 divide( _hist_integrated_yield_rho, _hist_integrated_yield_pion, _hist_integrated_rho_pion_ratio);
321 }
322
323
324
325 } // end finalize
326 //=================================================================================
327
328
329 private:
330
331 bool isHI;
332 static const int NHISTOS = 4;
333 static const int EVENT_TYPES = 2;
334 static const int PP = 0;
335 static const int PBPB = 1;
336
337
338 // pt spectrum for rho in pp@2.76 TeV
339 // Table 6 in HEPData and Fig. 5 in Article
340 // NOTE. The histo for pions is need for the rho/pion ratio
341 // Table 9 in HEPData and Fig. 9 in Article
342 CounterPtr _counterSOW_PP;
343 Histo1DPtr _hist_rho_PP;
344 Histo1DPtr _hist_pion_PP;
345 Estimate1DPtr _rho_pion_ratio_PP;
346
347 // pt spectra for rho in PbPb@2.76 TeV in 0-20%, 20-40%, 40-60%, 60-80%
348 // Tables 7-10 in HEPData and Fig. 6 in Article
349 // NOTE1. For EVENT_TYPES=PP, histos are pt spectra for rho in pp@2.76 TeV
350 // with the same binning as the pt spectra in PbPb.
351 // They are needed for RAA
352 // Tables 16-19 in HEPData and Fig. 11 in Article
353 // NOTE2. histos for pions are needed for the rho/pion ratio
354 // Tables 14 and 15 in HEPData and Fig. 10 in Article
355 Histo1DPtr _hist_rho[EVENT_TYPES][NHISTOS];
356 Histo1DPtr _hist_pion[EVENT_TYPES][NHISTOS];
357 CounterPtr _counterSOW[EVENT_TYPES][NHISTOS];
358 CounterPtr _counterNcoll[NHISTOS];
359
360 CounterPtr _counter_temp;
361 CounterPtr _counterNcoll_temp;
362
363 Estimate1DPtr _rho_pion_ratio[EVENT_TYPES][NHISTOS];
364 Estimate1DPtr _hist_RAA[NHISTOS];
365
366 // integrated yields vs. multiplicity and mean pT
367 // Table 11-12, Fig. 8, left and right
368 Histo1DPtr _hist_integrated_yield_rho;
369 Histo1DPtr _hist_integrated_yield_pion;
370 Estimate1DPtr _hist_integrated_rho_pion_ratio;
371
372 Profile1DPtr _hist_mean_pt_rho;
373
374
375 std::vector<std::pair<double, double>> _centrality_regions;
376 vector<double> int_edges;
377
378
379 };
380
381
382 RIVET_DECLARE_PLUGIN(ALICE_2019_I1672860);
383
384}
|