Rivet analyses referenceCMS_2011_I916908Traditional leading jet UE measurement at $\sqrt{s} = 0.9$ and 7 TeVExperiment: CMS (LHC) Inspire ID: 916908 Status: VALIDATED Authors:
Beam energies: (450.0, 450.0); (3500.0, 3500.0) GeV Run details:
A measurement of the underlying activity in scattering processes with a hard scale in the several-GeV region is performed in proton-proton collisions at Energies of 0.9 and 7 TeV, using data collected by the CMS experiment at the LHC. The production of charged particles with pseudorapidity |eta| < 2 and transverse momentum $p_\perp > 0.5$ GeV/$c$ is studied in the azimuthal region transverse to that of the leading set of charged particles forming a track-jet. Various comparisons are made between the two different energies and also beteen two sets of cuts on p_\perp for leading track jet p_\perp-leading $> 3$ GeV and pT-leading $> 20$ GeV. The activity is studied using 5 types of plots. Two profile plots for the multiplicity of charged particles and the scalar sum of p_\perp. and three distributions for the two previous quantities as well as p_\perp for all the particles in the transverse region. Source code: CMS_2011_I916908.cc 1// -*- C++ -*-
2
3#include "Rivet/Analysis.hh"
4#include "Rivet/Projections/FinalState.hh"
5#include "Rivet/Projections/ChargedFinalState.hh"
6#include "Rivet/Projections/FastJets.hh"
7
8namespace Rivet {
9
10
11 /// UE charged particles vs. leading jet
12 class CMS_2011_I916908 : public Analysis {
13 public:
14
15 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2011_I916908);
16
17
18 void init() {
19 const ChargedFinalState cfs(Cuts::abseta < 2.0 && Cuts::pT >= 500*MeV);
20 declare(cfs, "CFS");
21
22 const ChargedFinalState cfsforjet(Cuts::abseta < 2.5 && Cuts::pT >= 500*MeV);
23 const FastJets jetpro(cfsforjet, JetAlg::SISCONE, 0.5);
24 declare(jetpro, "Jets");
25
26 for (double eVal : allowedEnergies()) {
27 const string en = toString(int(eVal));
28 if (isCompatibleWithSqrtS(eVal)) _sqs = en;
29
30 if (en == "7000"s) {
31 book(_p[en+"Nch_vs_pT"], 1, 1, 1); // Nch vs. pT_max
32 book(_p[en+"Sum_vs_pT"], 2, 1, 1); // sum(pT) vs. pT_max
33 book(_h[en+"pT3_Nch"], 5, 1, 1); // transverse Nch, pT_max > 3GeV
34 book(_h[en+"pT3_Sum"], 6, 1, 1); // transverse sum(pT), pT_max > 3GeV
35 book(_h[en+"pT3_pT"], 7, 1, 1); // transverse pT, pT_max > 3GeV
36 book(_h[en+"pT20_Nch"], 8, 1, 1); // transverse Nch, pT_max > 20GeV
37 book(_h[en+"pT20_Sum"], 9, 1, 1); // transverse sum(pT), pT_max > 20GeV
38 book(_h[en+"pT20_pT"], 10, 1, 1); // transverse pT, pT_max > 20GeV
39 }
40 else if (en == "900"s) {
41 book(_p[en+"Nch_vs_pT"], 3, 1, 1); // Nch vs. pT_max
42 book(_p[en+"Sum_vs_pT"], 4, 1, 1); // sum(pT) vs. pT_max
43 book(_h[en+"pT3_Nch"], 11, 1, 1); // transverse Nch, pT_max > 3GeV
44 book(_h[en+"pT3_Sum"], 12, 1, 1); // transverse sum(pT), pT_max > 3GeV
45 book(_h[en+"pT3_pT"], 13, 1, 1); // transverse pT, pT_max > 3GeV
46 }
47 book(_c[en+"pT3"], "TMP/nch_tot_pT3"+en);
48 book(_c[en+"pT20"], "TMP/nch_tot_pT20"+en);
49 }
50 if (_sqs == "" && !merging()) {
51 throw BeamError("Invalid beam energy for " + name() + "\n");
52 }
53 }
54
55
56 /// Perform the per-event analysis
57 void analyze(const Event& event) {
58
59 // Find the lead jet, applying a restriction that the jets must be within |eta| < 2.
60 FourMomentum p_lead;
61 for (const Jet& j : apply<FastJets>(event, "Jets").jetsByPt(Cuts::pT > 1.0*GeV && Cuts::abseta < 2.0)) {
62 p_lead = j.momentum();
63 break;
64 }
65 if (p_lead.isZero()) vetoEvent;
66 const double philead = p_lead.phi();
67 const double pTlead = p_lead.pT();
68
69 Particles particles = apply<ChargedFinalState>(event, "CFS").particlesByPt();
70
71 int nTransverse = 0;
72 double ptSumTransverse = 0.;
73 for (const Particle& p : particles) {
74 double dphi = deltaPhi(philead, p.phi());
75 if (dphi>PI/3. && dphi<PI*2./3.) { // Transverse region
76 ++nTransverse;
77
78 const double pT = p.pT()/GeV;
79 ptSumTransverse += pT;
80
81 if (pTlead > 3*GeV) _h[_sqs+"pT3_pT"]->fill(pT/GeV);
82 if (_sqs == "7000"s && pTlead > 20*GeV) _h[_sqs+"pT20_pT"]->fill(pT/GeV);
83 }
84 }
85
86 const double area = 8./3. * PI;
87 _p[_sqs+"Nch_vs_pT"]->fill(pTlead/GeV, 1./area*nTransverse);
88 _p[_sqs+"Sum_vs_pT"]->fill(pTlead/GeV, 1./area*ptSumTransverse);
89 if (pTlead > 3.0*GeV) {
90 _h[_sqs+"pT3_Nch"]->fill(nTransverse);
91 _h[_sqs+"pT3_Sum"]->fill(ptSumTransverse/GeV);
92 _c[_sqs+"pT3"]->fill(nTransverse);
93 }
94 if (_sqs == "7000"s && pTlead > 20*GeV) {
95 _h[_sqs+"pT20_Nch"]->fill(nTransverse);
96 _h[_sqs+"pT20_Sum"]->fill(ptSumTransverse/GeV);
97 _c[_sqs+"pT20"]->fill(nTransverse);
98 }
99 }
100
101
102 /// Normalise histograms etc., after the run
103 void finalize() {
104
105 for (double eVal : allowedEnergies()) {
106 const string en = toString(int(eVal));
107 if (_h[en+"pT3_pT"]->sumW()) {
108 normalize(_h[en+"pT3_pT"], dbl(*_c[en+"pT3"]) / _h[en+"pT3_Nch"]->sumW());
109 }
110 if (en == "7000"s && _h[en+"pT20_pT"]->sumW()) {
111 normalize(_h[en+"pT20_pT"], dbl(*_c[en+"pT20"]) / _h[en+"pT20_Nch"]->sumW());
112 }
113 }
114 for (auto& item : _h) {
115 if (item.first.find("_pT") != string::npos) continue;
116 normalize(item.second);
117 }
118
119 }
120
121
122
123 private:
124
125 /// @{
126 map<string,CounterPtr> _c;
127 map<string,Profile1DPtr> _p;
128 map<string,Histo1DPtr> _h;
129
130 string _sqs = "";
131 /// @}
132
133 };
134
135
136
137 RIVET_DECLARE_ALIASED_PLUGIN(CMS_2011_I916908, CMS_2011_S9120041);
138
139}
|