Rivet analyses referenceCMS_2011_S8884919Measurement of the NSD charged particle multiplicity at $\sqrt{s} = 0.9$, 2.36, and 7 TeVExperiment: CMS (LHC) Inspire ID: 879315 Status: VALIDATED Authors:
Beam energies: (450.0, 450.0); (1180.0, 1180.0); (3500.0, 3500.0) GeV Run details:
Measurements of primary charged hadron multiplicity distributions are presented for non-single-diffractive events in proton-proton collisions at centre-of-mass energies of $\sqrt{s} = 0.9$, 2.36, and 7 TeV, in five pseudorapidity ranges from $|\eta| < 0.5$ to $|\eta| < 2.4$. The data were collected with the minimum-bias trigger of the CMS experiment during the LHC commissioning runs in 2009 and the 7 TeV run in 2010. The average transverse momentum as a function of the multiplicity is also presented. The measurement of higher-order moments of the multiplicity distribution confirms the violation of Koba-Nielsen-Olesen scaling that has been observed at lower energies. Beam energy must be specified (in GeV) as analysis option "ENERGY" when rivet-merging samples. Source code: CMS_2011_S8884919.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/ChargedFinalState.hh"
4#include "Rivet/Projections/Beam.hh"
5using namespace std;
6
7namespace Rivet {
8
9
10 /// Measurement of the NSD charged particle multiplicity at 0.9, 2.36, and 7 TeV
11 class CMS_2011_S8884919 : public Analysis {
12 public:
13
14 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2011_S8884919);
15
16
17 void init() {
18 ChargedFinalState cfs((Cuts::etaIn(-2.4, 2.4)));
19 declare(cfs, "CFS");
20
21 // eta bins
22 _etabins.push_back(0.5);
23 _etabins.push_back(1.0);
24 _etabins.push_back(1.5);
25 _etabins.push_back(2.0);
26 _etabins.push_back(2.4) ;
27
28 if (isCompatibleWithSqrtS(900)) {
29 for (size_t ietabin=0; ietabin < _etabins.size(); ietabin++) {
30 _h_dNch_dn.push_back( Histo1DPtr() );
31 book( _h_dNch_dn.back(), 2 + ietabin, 1, 1);
32 }
33 book(_h_dNch_dn_pt500_eta24 ,20, 1, 1);
34 book(_h_dmpt_dNch_eta24 ,23, 1, 1);
35 }
36
37 if (isCompatibleWithSqrtS(2360)) {
38 for (size_t ietabin=0; ietabin < _etabins.size(); ietabin++) {
39 _h_dNch_dn.push_back( Histo1DPtr() );
40 book(_h_dNch_dn.back(), 7 + ietabin, 1, 1);
41 }
42 book(_h_dNch_dn_pt500_eta24 ,21, 1, 1);
43 book(_h_dmpt_dNch_eta24 ,24, 1, 1);
44 }
45
46 if (isCompatibleWithSqrtS(7000)) {
47 for (size_t ietabin=0; ietabin < _etabins.size(); ietabin++) {
48 _h_dNch_dn.push_back( Histo1DPtr() );
49 book(_h_dNch_dn.back(), 12 + ietabin, 1, 1);
50 }
51 book(_h_dNch_dn_pt500_eta24 ,22, 1, 1);
52 book(_h_dmpt_dNch_eta24 ,25, 1, 1);
53 }
54 }
55
56
57 void analyze(const Event& event) {
58
59 // Get the charged particles
60 const ChargedFinalState& charged = apply<ChargedFinalState>(event, "CFS");
61
62 // Resetting the multiplicity for the event to 0;
63 vector<int> _nch_in_Evt;
64 vector<int> _nch_in_Evt_pt500;
65 _nch_in_Evt.assign(_etabins.size(), 0);
66 _nch_in_Evt_pt500.assign(_etabins.size(), 0);
67 double sumpt = 0;
68
69 // Loop over particles in event
70 for (const Particle& p : charged.particles()) {
71 // Selecting only charged hadrons
72 if (! PID::isHadron(p.pid())) continue;
73
74 double pT = p.pT();
75 double eta = p.eta();
76 sumpt += pT;
77 for (size_t ietabin = _etabins.size(); ietabin > 0; --ietabin) {
78 if (fabs(eta) > _etabins[ietabin-1]) break;
79 ++_nch_in_Evt[ietabin-1];
80 if (pT > 0.5/GeV) ++_nch_in_Evt_pt500[ietabin-1];
81 }
82 }
83
84 // Filling multiplicity-dependent histogramms
85 for (size_t ietabin = 0; ietabin < _etabins.size(); ietabin++) {
86 _h_dNch_dn[ietabin]->fill(_nch_in_Evt[ietabin]);
87 }
88
89 // Do only if eta bins are the needed ones
90 if (_etabins[4] == 2.4 && _etabins[0] == 0.5) {
91 if (_nch_in_Evt[4] != 0) {
92 _h_dmpt_dNch_eta24->fill(_nch_in_Evt[4], sumpt/GeV / _nch_in_Evt[4]);
93 }
94 _h_dNch_dn_pt500_eta24->fill(_nch_in_Evt_pt500[4]);
95 } else {
96 MSG_WARNING("You changed the number of eta bins, but forgot to propagate it everywhere !!");
97 }
98 }
99
100
101 void finalize() {
102 for (size_t ietabin = 0; ietabin < _etabins.size(); ietabin++){
103 normalize(_h_dNch_dn[ietabin]);
104 }
105 normalize(_h_dNch_dn_pt500_eta24);
106 }
107
108
109 private:
110
111 /// @{
112 vector<Histo1DPtr> _h_dNch_dn;
113 Histo1DPtr _h_dNch_dn_pt500_eta24;
114 Profile1DPtr _h_dmpt_dNch_eta24;
115 /// @}
116
117 vector<double> _etabins;
118
119 };
120
121
122
123 RIVET_DECLARE_ALIASED_PLUGIN(CMS_2011_S8884919, CMS_2011_I879315);
124
125}
|