Rivet analyses referenceATLAS_2016_I1426695Track-based minimum bias at 8 TeVExperiment: ATLAS (LHC) Inspire ID: 1426695 Status: VALIDATED Authors:
Beam energies: (4000.0, 4000.0) GeV Run details:
This paper presents measurements of distributions of charged particles which are produced in proton-proton collisions at a centre-of-mass energy of $\sqrt{s} = 8$ TeV and recorded by the ATLAS detector at the LHC. A special dataset recorded in 2012 with a small number of interactions per beam crossing (below 0.004) and corresponding to an integrated luminosity of 160 $\mu\text{b}^{-1}$ was used. A minimum-bias trigger was utilised to select a data sample of more than 9 million collision events. The multiplicity, pseudorapidity, and transverse momentum distributions of charged particles are shown in different regions of kinematics and charged-particle multiplicity, including measurements of final states at high multiplicity. The results are corrected for detector effects and are compared to the predictions of various Monte Carlo event generator models which simulate the full hadronic final state. Source code: ATLAS_2016_I1426695.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/ChargedFinalState.hh"
4
5namespace Rivet {
6
7
8 /// @brief Track-based minimum bias at 8 TeV
9 class ATLAS_2016_I1426695 : public Analysis {
10 public:
11
12 //phase space regions
13 enum regionID {
14 k_pt100_nch2 = 0,
15 k_pt500_nch1,
16 k_pt500_nch6,
17 k_pt500_nch20,
18 k_pt500_nch50,
19 kNregions
20 };
21
22
23 /// Constructor
24 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_I1426695);
25
26 /// @name Analysis methods
27 /// @{
28
29 /// Book histograms and initialise projections before the run
30 void init() {
31
32 for (int iR=0; iR < kNregions; ++iR) {
33 book(_sumW[iR], "_sumW" + to_str(iR)) ;
34 }
35
36 // Initialise and register projections
37 declare(ChargedFinalState(Cuts::abseta < 2.5 && Cuts::pT > 100*MeV), "CFS_100");
38 declare(ChargedFinalState(Cuts::abseta < 2.5 && Cuts::pT > 500*MeV), "CFS_500");
39
40 // Book histograms
41 for (int iR=0; iR < kNregions; ++iR) {
42 if (iR == k_pt100_nch2 || iR == k_pt500_nch1) {
43 book(_hist_nch [iR], 2 + iR, 1, 1);
44 book(_hist_ptnch[iR], 14 + iR, 1, 1);
45 }
46 book(_hist_pt [iR], 4 + iR, 1, 1);
47 book(_hist_eta[iR], 9 + iR, 1, 1);
48 }
49 }
50
51 void fillPtEtaNch(const Particles particles, int nMin, int iRegion) {
52
53 //skip if event fails multiplicity cut
54 int nch =particles.size();
55 if (nch < nMin) return;
56
57 _sumW[iRegion]->fill();
58
59 // Fill nch
60 if (iRegion == k_pt100_nch2 || iRegion == k_pt500_nch1) {
61 _hist_nch[iRegion]->fill(nch);
62 }
63
64 for (const Particle &p : particles) {
65 // Loop over particles, fill pT, eta and ptnch
66 const double pt = p.pT()/GeV;
67 const double eta = p.eta();
68
69 _hist_pt [iRegion]->fill(pt , 1.0/pt);
70 _hist_eta[iRegion]->fill(eta);
71
72 if (iRegion == k_pt100_nch2 || iRegion == k_pt500_nch1) {
73 _hist_ptnch[iRegion]->fill(nch, pt);
74 }
75 } //end loop over particles
76 }
77
78 /// Perform the per-event analysis
79 void analyze(const Event& event) {
80
81 // Get charged particles, omitting some strange heavies
82 const Cut& pcut = (
83 (Cuts::abspid!=PID::SIGMAMINUS) && (Cuts::abspid!=PID::SIGMAPLUS) &&
84 (Cuts::abspid!=PID::XIMINUS) && (Cuts::abspid!=PID::OMEGAMINUS));
85 const Particles& p_100 = apply<ChargedFinalState>(event, "CFS_100").particles(pcut);
86 const Particles& p_500 = apply<ChargedFinalState>(event, "CFS_500").particles(pcut);
87
88 fillPtEtaNch(p_100, 2, 0);
89 fillPtEtaNch(p_500, 1, 1);
90 fillPtEtaNch(p_500, 6, 2);
91 fillPtEtaNch(p_500, 20, 3);
92 fillPtEtaNch(p_500, 50, 4);
93 }
94
95
96 void finalize() {
97
98 for (int iR = 0; iR < kNregions; ++iR) {
99 if (_sumW[iR]->val() > 0) {
100 if (iR == k_pt100_nch2 || iR == k_pt500_nch1) {
101 scale(_hist_nch[iR], 1.0/ *_sumW[iR]);
102 }
103 scale(_hist_pt [iR], 1.0/ dbl(*_sumW[iR])/TWOPI/5.);
104 scale(_hist_eta[iR], 1.0/ *_sumW[iR]);
105 }
106 }
107 }
108
109 /// @}
110
111
112 private:
113
114 CounterPtr _sumW[kNregions];
115
116 /// @name Histograms
117 Histo1DPtr _hist_nch [kNregions];
118 Histo1DPtr _hist_pt [kNregions];
119 Histo1DPtr _hist_eta [kNregions];
120 Profile1DPtr _hist_ptnch [kNregions];
121
122 };
123
124
125 RIVET_DECLARE_PLUGIN(ATLAS_2016_I1426695);
126
127}
|