Rivet analyses referenceCMS_2015_PAS_FSQ_15_007Underlying Event Measurements with Leading Particles and Jets in pp collisions at 13 TeVExperiment: CMS (LHC) Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
A measurement of the underlying event (UE) activity is performed in proton-proton collisions at the centre-of-mass energy of 13 TeV with the CMS experiment at the LHC. The measurement is performed using leading charged-particles as well as lead- ing charged-particle jets as reference objects. A leading charged-particle or charged- particle jet is required to be produced in the central pseudorapidity region ( |eta| < 2) jet and with transverse momentum pT greater than 0.5 (pT greater than 1) GeV for leading charged-particle (charge particle jet). The UE activity is measured in terms of the average multiplicity and scalar sum of pT of the charged-particles, in the azimuthal region transverse to the direction of the leading reference object. Source code: CMS_2015_PAS_FSQ_15_007.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/ChargedFinalState.hh"
5#include "Rivet/Projections/FastJets.hh"
6
7namespace Rivet {
8
9
10 // UE charged particles vs. leading jet
11 class CMS_2015_PAS_FSQ_15_007 : public Analysis {
12 public:
13 /// Constructor
14 CMS_2015_PAS_FSQ_15_007() : Analysis("CMS_2015_PAS_FSQ_15_007") {}
15
16 void init() {
17 const ChargedFinalState cfs(Cuts::abseta < 2.0 && Cuts::pt > 500 * MeV);
18 declare(cfs, "CFS");
19
20 const ChargedFinalState cfsforjet(Cuts::abseta < 2.5 && Cuts::pt > 500 * MeV);
21 const FastJets jetpro(cfsforjet, JetAlg::SISCONE, 0.5);
22 declare(jetpro, "Jets");
23
24 book(_h_PtSum_vs_leadTrackPt_transMin, 1, 1, 1);
25 book(_h_PtSum_vs_leadTrackPt_transMax, 2, 1, 1);
26 book(_h_PtSum_vs_leadTrackPt_transDiff, 3, 1, 1);
27 book(_h_PtSum_vs_leadTrackPt_transAvg, 4, 1, 1);
28
29 book(_h_Nch_vs_leadTrackPt_transMin, 5, 1, 1);
30 book(_h_Nch_vs_leadTrackPt_transMax, 6, 1, 1);
31 book(_h_Nch_vs_leadTrackPt_transDiff, 7, 1, 1);
32 book(_h_Nch_vs_leadTrackPt_transAvg, 8, 1, 1);
33
34 book(_h_PtSum_vs_leadJetPt_transMin, 9, 1, 1);
35 book(_h_PtSum_vs_leadJetPt_transMax, 10, 1, 1);
36 book(_h_PtSum_vs_leadJetPt_transDiff, 11, 1, 1);
37 book(_h_PtSum_vs_leadJetPt_transAvg, 12, 1, 1);
38
39 book(_h_Nch_vs_leadJetPt_transMin, 13, 1, 1);
40 book(_h_Nch_vs_leadJetPt_transMax, 14, 1, 1);
41 book(_h_Nch_vs_leadJetPt_transDiff, 15, 1, 1);
42 book(_h_Nch_vs_leadJetPt_transAvg, 16, 1, 1);
43 }
44
45 /// Perform the per-event analysis
46 void analyze(const Event& event) {
47
48 // Find the lead jet, applying a restriction that the jets must be within |eta| < 2.
49 FourMomentum p_leadjet;
50 Jets js = apply<FastJets>(event, "Jets").jetsByPt(Cuts::pT > 1*GeV && Cuts::abseta < 2.0);
51 for (const Jet& j : js) {
52 p_leadjet = j.momentum();
53 break;
54 }
55
56 // Find the lead track, also within |eta| < 2
57 FourMomentum p_leadtrack;
58 const Particles ts = apply<ChargedFinalState>(event, "CFS").particlesByPt(Cuts::abseta < 2.0 && Cuts::pT > 0.5*GeV);
59 for (const Particle& t : ts) {
60 p_leadtrack = t.momentum();
61 break;
62 }
63
64 if (p_leadjet.isZero() && p_leadtrack.isZero())
65 vetoEvent;
66 const double phileadjet = p_leadjet.phi();
67 const double pTleadjet = p_leadjet.pT();
68
69 const double phileadtrack = p_leadtrack.phi();
70 const double pTleadtrack = p_leadtrack.pT();
71
72 Particles particles = apply<ChargedFinalState>(event, "CFS").particlesByPt();
73
74 // double ptSumTransverse_leadjet = 0.;
75 int nTransverse1_leadjet = 0;
76 double ptSumTransverse1_leadjet = 0.;
77 int nTransverse2_leadjet = 0;
78 double ptSumTransverse2_leadjet = 0.;
79 int nTransverseMin_leadjet = 0;
80 double ptSumTransverseMin_leadjet = 0.;
81 int nTransverseMax_leadjet = 0;
82 double ptSumTransverseMax_leadjet = 0.;
83
84 // double ptSumTransverse_leadtrack = 0.;
85 int nTransverse1_leadtrack = 0;
86 double ptSumTransverse1_leadtrack = 0.;
87 int nTransverse2_leadtrack = 0;
88 double ptSumTransverse2_leadtrack = 0.;
89 int nTransverseMin_leadtrack = 0;
90 double ptSumTransverseMin_leadtrack = 0.;
91 int nTransverseMax_leadtrack = 0;
92 double ptSumTransverseMax_leadtrack = 0.;
93 // double ptSumTowards_leadtrack = 0.;
94 // double ptSumAway_leadtrack = 0.;
95
96 for (const Particle& p : particles) {
97 const double pT = p.pT() / GeV;
98
99 if (!p_leadjet.isZero()) {
100 double dphi_leadjet = p.phi() - phileadjet;
101 while (dphi_leadjet > PI) {
102 dphi_leadjet = dphi_leadjet - 2.0 * PI;
103 }
104 while (dphi_leadjet < -PI) {
105 dphi_leadjet = dphi_leadjet + 2. * PI;
106 }
107
108 if (dphi_leadjet > PI / 3. && dphi_leadjet < PI * 2. / 3.) { // Transverse1 region
109 // ptSumTransverse_leadjet += pT;
110 nTransverse1_leadjet++;
111 ptSumTransverse1_leadjet += pT;
112 }
113
114 if (dphi_leadjet < -PI / 3. && dphi_leadjet > -PI * 2. / 3.) { // Transverse2 region
115 // ptSumTransverse_leadjet += pT;
116 nTransverse2_leadjet++;
117 ptSumTransverse2_leadjet += pT;
118 }
119
120 } //jet found
121
122 if (!p_leadtrack.isZero()) {
123 double dphi_leadtrack = p.phi() - phileadtrack;
124 while (dphi_leadtrack > PI) {
125 dphi_leadtrack = dphi_leadtrack - 2.0 * PI;
126 }
127 while (dphi_leadtrack < -PI) {
128 dphi_leadtrack = dphi_leadtrack + 2. * PI;
129 }
130
131 if (dphi_leadtrack > PI / 3. && dphi_leadtrack < PI * 2. / 3.) { // Transverse1 region
132 // ptSumTransverse_leadtrack += pT;
133 nTransverse1_leadtrack++;
134 ptSumTransverse1_leadtrack += pT;
135 }
136
137 if (dphi_leadtrack < -PI / 3. && dphi_leadtrack > -PI * 2. / 3.) { // Transverse2 region
138 // ptSumTransverse_leadtrack += pT;
139 nTransverse2_leadtrack++;
140 ptSumTransverse2_leadtrack += pT;
141 }
142
143 } //< track found
144
145 } //< loop over particles
146
147 const double fullarea = 8. / 3. * PI;
148 const double halfarea = 4. / 3. * PI;
149
150 if (!p_leadjet.isZero()) {
151 if (nTransverse2_leadjet > nTransverse1_leadjet) {
152 nTransverseMax_leadjet = nTransverse2_leadjet;
153 nTransverseMin_leadjet = nTransverse1_leadjet;
154 } else {
155 nTransverseMax_leadjet = nTransverse1_leadjet;
156 nTransverseMin_leadjet = nTransverse2_leadjet;
157 }
158
159 if (ptSumTransverse2_leadjet > ptSumTransverse1_leadjet) {
160 ptSumTransverseMax_leadjet = ptSumTransverse2_leadjet;
161 ptSumTransverseMin_leadjet = ptSumTransverse1_leadjet;
162 }
163
164 else {
165 ptSumTransverseMax_leadjet = ptSumTransverse1_leadjet;
166 ptSumTransverseMin_leadjet = ptSumTransverse2_leadjet;
167 }
168
169 _h_Nch_vs_leadJetPt_transDiff->fill(pTleadjet / GeV,
170 1. / halfarea * (nTransverseMax_leadjet - nTransverseMin_leadjet));
171 _h_PtSum_vs_leadJetPt_transDiff->fill(
172 pTleadjet / GeV, 1. / halfarea * (ptSumTransverseMax_leadjet - ptSumTransverseMin_leadjet));
173 _h_Nch_vs_leadJetPt_transAvg->fill(pTleadjet / GeV,
174 1. / fullarea * (nTransverseMax_leadjet + nTransverseMin_leadjet));
175 _h_PtSum_vs_leadJetPt_transAvg->fill(pTleadjet / GeV,
176 1. / fullarea * (ptSumTransverseMax_leadjet + ptSumTransverseMin_leadjet));
177 _h_Nch_vs_leadJetPt_transMax->fill(pTleadjet / GeV, 1. / halfarea * nTransverseMax_leadjet);
178 _h_PtSum_vs_leadJetPt_transMax->fill(pTleadjet / GeV, 1. / halfarea * ptSumTransverseMax_leadjet);
179 _h_Nch_vs_leadJetPt_transMin->fill(pTleadjet / GeV, 1. / halfarea * nTransverseMin_leadjet);
180 _h_PtSum_vs_leadJetPt_transMin->fill(pTleadjet / GeV, 1. / halfarea * ptSumTransverseMin_leadjet);
181
182 } //for leading jet
183
184 if (!p_leadtrack.isZero()) {
185 if (nTransverse2_leadtrack > nTransverse1_leadtrack) {
186 nTransverseMax_leadtrack = nTransverse2_leadtrack;
187 nTransverseMin_leadtrack = nTransverse1_leadtrack;
188 }
189
190 else {
191 nTransverseMax_leadtrack = nTransverse1_leadtrack;
192 nTransverseMin_leadtrack = nTransverse2_leadtrack;
193 }
194
195 if (ptSumTransverse2_leadtrack > ptSumTransverse1_leadtrack) {
196 ptSumTransverseMax_leadtrack = ptSumTransverse2_leadtrack;
197 ptSumTransverseMin_leadtrack = ptSumTransverse1_leadtrack;
198 }
199
200 else {
201 ptSumTransverseMax_leadtrack = ptSumTransverse1_leadtrack;
202 ptSumTransverseMin_leadtrack = ptSumTransverse2_leadtrack;
203 }
204
205 _h_Nch_vs_leadTrackPt_transDiff->fill(pTleadtrack / GeV,
206 1. / halfarea * (nTransverseMax_leadtrack - nTransverseMin_leadtrack));
207 _h_PtSum_vs_leadTrackPt_transDiff->fill(
208 pTleadtrack / GeV, 1. / halfarea * (ptSumTransverseMax_leadtrack - ptSumTransverseMin_leadtrack));
209 _h_Nch_vs_leadTrackPt_transAvg->fill(pTleadtrack / GeV,
210 1. / fullarea * (nTransverseMax_leadtrack + nTransverseMin_leadtrack));
211 _h_PtSum_vs_leadTrackPt_transAvg->fill(
212 pTleadtrack / GeV, 1. / fullarea * (ptSumTransverseMax_leadtrack + ptSumTransverseMin_leadtrack));
213
214 _h_Nch_vs_leadTrackPt_transMax->fill(pTleadtrack / GeV, 1. / halfarea * nTransverseMax_leadtrack);
215 _h_PtSum_vs_leadTrackPt_transMax->fill(pTleadtrack / GeV, 1. / halfarea * ptSumTransverseMax_leadtrack);
216 _h_Nch_vs_leadTrackPt_transMin->fill(pTleadtrack / GeV, 1. / halfarea * nTransverseMin_leadtrack);
217 _h_PtSum_vs_leadTrackPt_transMin->fill(pTleadtrack / GeV, 1. / halfarea * ptSumTransverseMin_leadtrack);
218
219 } //for leading track
220 }
221
222 /// Normalise histograms etc., after the run
223 void finalize() {}
224
225 private:
226 Profile1DPtr _h_Nch_vs_leadJetPt_transMax;
227 Profile1DPtr _h_PtSum_vs_leadJetPt_transMax;
228 Profile1DPtr _h_Nch_vs_leadJetPt_transMin;
229 Profile1DPtr _h_PtSum_vs_leadJetPt_transMin;
230 Profile1DPtr _h_Nch_vs_leadJetPt_transDiff;
231 Profile1DPtr _h_PtSum_vs_leadJetPt_transDiff;
232 Profile1DPtr _h_Nch_vs_leadJetPt_transAvg;
233 Profile1DPtr _h_PtSum_vs_leadJetPt_transAvg;
234
235 Profile1DPtr _h_Nch_vs_leadTrackPt_transMax;
236 Profile1DPtr _h_PtSum_vs_leadTrackPt_transMax;
237 Profile1DPtr _h_Nch_vs_leadTrackPt_transMin;
238 Profile1DPtr _h_PtSum_vs_leadTrackPt_transMin;
239 Profile1DPtr _h_Nch_vs_leadTrackPt_transDiff;
240 Profile1DPtr _h_PtSum_vs_leadTrackPt_transDiff;
241 Profile1DPtr _h_Nch_vs_leadTrackPt_transAvg;
242 Profile1DPtr _h_PtSum_vs_leadTrackPt_transAvg;
243 };
244
245 RIVET_DECLARE_PLUGIN(CMS_2015_PAS_FSQ_15_007);
246} // namespace Rivet
|