rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2011_I879315

Measurement of the NSD charged particle multiplicity at $\sqrt{s} = 0.9$, 2.36, and 7 TeV
Experiment: CMS (LHC)
Inspire ID: 879315
Status: VALIDATED
Authors:
  • Romain Rougny
References: Beams: p+ p+
Beam energies: (450.0, 450.0); (1180.0, 1180.0); (3500.0, 3500.0) GeV
Run details:
  • Non-single-diffractive (NSD) events only. Should include double-diffractive (DD) events and non-diffractive (ND) events but NOT single-diffractive (SD) events. For example, in Pythia6 the SD processes to be turned off are 92 and 93 and in Pythia8 the SD processes are 103 and 104 (also called SoftQCD:singleDiffractive).

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.

Source code: CMS_2011_I879315.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/ChargedFinalState.hh"
  4#include "Rivet/Projections/Beam.hh"
  5#include "Rivet/Tools/HistoGroup.hh"
  6using namespace std;
  7
  8namespace Rivet {
  9
 10
 11  /// Measurement of the NSD charged particle multiplicity at 0.9, 2.36, and 7 TeV
 12  class CMS_2011_I879315 : public Analysis {
 13  public:
 14
 15    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2011_I879315);
 16
 17
 18    void init() {
 19      ChargedFinalState cfs(Cuts::abseta < 2.4);
 20      declare(cfs, "CFS");
 21
 22      for (double eVal : allowedEnergies()) {
 23        const string en = toString(int(eVal));
 24        if (isCompatibleWithSqrtS(eVal))  _sqs = en;
 25
 26        size_t offset = 0;
 27        if (en == "2360"s)  offset = 5;
 28        else if (en == "7000"s)  offset = 10;
 29        book(_g[en+"dNch_dn"], {0.0, 0.5, 1.0, 1.5, 2.0, 2.4});
 30        for (auto& b : _g[en+"dNch_dn"]->bins()) {
 31          book(b, b.index()+offset+1, 1, 1);
 32        }
 33
 34        if (en == "2360"s)  offset = 1;
 35        else if (en == "7000"s)  offset = 2;
 36        book(_h[en+"dNch_dn_pt500_eta24"], 20+offset, 1, 1);
 37        book(_p[en+"dmpt_dNch_eta24"],     23+offset, 1, 1);
 38      }
 39      if (_sqs == "" && !merging()) {
 40        throw BeamError("Invalid beam energy for " + name() + "\n");
 41      }
 42    }
 43
 44
 45    void analyze(const Event& event) {
 46
 47      // Get the charged particles
 48      const ChargedFinalState& charged = apply<ChargedFinalState>(event, "CFS");
 49
 50      // Resetting the multiplicity for the event to 0;
 51      size_t nBins = _g[_sqs+"dNch_dn"]->numBins();
 52      vector<int> _nch_in_Evt(nBins, 0);
 53      vector<int> _nch_in_Evt_pt500(nBins, 0);
 54      double sumpt = 0;
 55
 56      // Loop over particles in event
 57      for (const Particle& p : charged.particles()) {
 58        // Selecting only charged hadrons
 59        if (! PID::isHadron(p.pid())) continue;
 60
 61        const double pT = p.pT(); sumpt += pT;
 62        const double abseta = p.abseta();
 63        size_t ieta = nBins;
 64        while (ieta--) {
 65          if (abseta > _g[_sqs+"dNch_dn"]->bin(ieta+1).xMax()) break;
 66          ++_nch_in_Evt[ieta];
 67          if (pT > 0.5*GeV) ++_nch_in_Evt_pt500[ieta];
 68        }
 69      }
 70
 71      // Filling multiplicity-dependent histogramms
 72      for (auto& b : _g[_sqs+"dNch_dn"]->bins()) {
 73        b->fill(_nch_in_Evt[b.index()-1]);
 74      }
 75
 76      if (_nch_in_Evt[nBins-1] != 0) {
 77        _p[_sqs+"dmpt_dNch_eta24"]->fill(_nch_in_Evt[nBins-1], sumpt/GeV / _nch_in_Evt[nBins-1]);
 78      }
 79      _h[_sqs+"dNch_dn_pt500_eta24"]->fill(_nch_in_Evt_pt500[nBins-1]);
 80    }
 81
 82
 83    void finalize() {
 84      normalize(_g);
 85      normalize(_h);
 86    }
 87
 88
 89  private:
 90
 91    /// @{
 92    map<string, Histo1DPtr> _h;
 93    map<string, Profile1DPtr> _p;
 94    map<string, Histo1DGroupPtr> _g;
 95
 96    string _sqs = "";
 97    /// @}
 98
 99  };
100
101
102
103  RIVET_DECLARE_ALIASED_PLUGIN(CMS_2011_I879315, CMS_2011_S8884919);
104
105}