Rivet analyses referenceATLAS_2012_I1125575Studies of the underlying event at 7 TeV with track-jetsExperiment: ATLAS (LHC) Inspire ID: 1125575 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
Distributions sensitive to the underlying event are studied in events containing one or more charged particles. Jets are reconstructed using the anti-$k_t$ algorithm with radius parameter $R$ varying between 0.2 and 1.0. Distributions of the charged-particle multiplicity, the scalar sum of the transverse momentum of charged particles, and the average charged-particle pT are measured as functions of $p_\perp^\text{jet}$ in regions transverse to and opposite the leading jet for $4 \text{GeV} < p_\perp^\text{jet} < 100 \text{GeV}$. In addition, the $R$-dependence of the mean values of these observables is studied. Source code: ATLAS_2012_I1125575.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/ChargedFinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Tools/BinnedHistogram.hh"
6
7namespace Rivet {
8
9
10 /// ATLAS charged particle jet underlying event and jet radius dependence
11 class ATLAS_2012_I1125575 : public Analysis {
12 public:
13
14 /// @name Constructors etc.
15 //@{
16
17 /// Constructor
18 ATLAS_2012_I1125575()
19 : Analysis("ATLAS_2012_I1125575")
20 { }
21
22 //@}
23
24
25 /// @name Analysis methods
26 //@{
27
28 /// Book histograms and initialise projections before the run
29 void init() {
30
31 const ChargedFinalState jet_input((Cuts::etaIn(-2.5, 2.5) && Cuts::pT >= 0.5*GeV));
32 declare(jet_input, "JET_INPUT");
33
34 const ChargedFinalState track_input((Cuts::etaIn(-1.5, 1.5) && Cuts::pT >= 0.5*GeV));
35 declare(track_input, "TRACK_INPUT");
36
37 const FastJets jets02(jet_input, FastJets::ANTIKT, 0.2);
38 declare(jets02, "JETS_02");
39
40 const FastJets jets04(jet_input, FastJets::ANTIKT, 0.4);
41 declare(jets04, "JETS_04");
42
43 const FastJets jets06(jet_input, FastJets::ANTIKT, 0.6);
44 declare(jets06, "JETS_06");
45
46 const FastJets jets08(jet_input, FastJets::ANTIKT, 0.8);
47 declare(jets08, "JETS_08");
48
49 const FastJets jets10(jet_input, FastJets::ANTIKT, 1.0);
50 declare(jets10, "JETS_10");
51
52 // Mean number of tracks
53 initializeProfiles(_h_meanNch, 1);
54
55 // Mean of the average track pT in each region
56 initializeProfiles(_h_meanPtAvg, 2);
57
58 // Mean of the scalar sum of track pT in each region
59 initializeProfiles(_h_meanPtSum, 3);
60
61 // Distribution of Nch, in bins of leading track-jet pT
62 initializeHistograms(_h_Nch, 4);
63
64 // Distribution of average track-jet pT, in bins of leading track-jet pT
65 initializeHistograms(_h_PtAvg, 5);
66
67 // Distribution of sum of track-jet pT, in bins of leading track-jet pT
68 initializeHistograms(_h_PtSum, 6);
69
70 for (int i = 0; i < 5; ++i)
71 book(_nEvents[i], "nEvents_"+to_str(i));
72 }
73
74
75 void initializeProfiles(Profile1DPtr plots[5][2], int distribution) {
76 for (int i = 0; i < 5; ++i) {
77 for (int j = 0; j < 2; ++j) {
78 book(plots[i][j] ,distribution, i+1, j+1);
79 }
80 }
81 }
82
83
84 void initializeHistograms(BinnedHistogram plots[5][2], int distribution) {
85 Scatter2D refscatter = refData(1, 1, 1);
86 for (int i = 0; i < 5; ++i) {
87 for (int y = 0; y < 2; ++y) {
88 for (size_t j = 0; j < refscatter.numPoints(); ++j) {
89 int histogram_number = ((j+1)*2)-((y+1)%2);
90 double low_edge = refscatter.point(j).xMin();
91 double high_edge = refscatter.point(j).xMax();
92 Histo1DPtr tmp;
93 plots[i][y].add(low_edge, high_edge, book(tmp, distribution, i+1, histogram_number));
94 }
95 }
96 }
97 }
98
99
100 /// Perform the per-event analysis
101 void analyze(const Event& event) {
102 vector<Jets*> all_jets;
103 Jets jets_02 = apply<FastJets>(event, "JETS_02").jetsByPt(Cuts::pT > 4*GeV && Cuts::abseta < 1.5);
104 all_jets.push_back(&jets_02);
105 Jets jets_04 = apply<FastJets>(event, "JETS_04").jetsByPt(Cuts::pT > 4*GeV && Cuts::abseta < 1.5);
106 all_jets.push_back(&jets_04);
107 Jets jets_06 = apply<FastJets>(event, "JETS_06").jetsByPt(Cuts::pT > 4*GeV && Cuts::abseta < 1.5);
108 all_jets.push_back(&jets_06);
109 Jets jets_08 = apply<FastJets>(event, "JETS_08").jetsByPt(Cuts::pT > 4*GeV && Cuts::abseta < 1.5);
110 all_jets.push_back(&jets_08);
111 Jets jets_10 = apply<FastJets>(event, "JETS_10").jetsByPt(Cuts::pT > 4*GeV && Cuts::abseta < 1.5);
112 all_jets.push_back(&jets_10);
113
114 // Count the number of tracks in the away and transverse regions, for each set of jets
115 double n_ch[5][2] = { {0,0}, {0,0}, {0,0}, {0,0}, {0,0} };
116 // Also add up the sum pT
117 double sumpt[5][2] = { {0,0}, {0,0}, {0,0}, {0,0}, {0,0} };
118 // ptmean = sumpt / n_ch
119 double ptavg[5][2] = { {0,0}, {0,0}, {0,0}, {0,0}, {0,0} };
120 // lead jet pT defines which bin we want to fill
121 double lead_jet_pts[5] = {0.0};
122
123 // Loop over each of the jet radii:
124 for (int i = 0; i < 5; ++i) {
125 if (all_jets[i]->size() < 1) continue;
126
127 // Find the lead jet pT
128 lead_jet_pts[i] = all_jets[i]->at(0).pT();
129
130 // Loop over each of the charged particles
131 const Particles& tracks = apply<ChargedFinalState>(event, "TRACK_INPUT").particlesByPt();
132 for(const Particle& t : tracks) {
133
134 // Get the delta-phi between the track and the leading jet
135 double dphi = deltaPhi(all_jets[i]->at(0), t);
136
137 // Find out which region this puts it in.
138 // 0 = away region, 1 = transverse region, 2 = toward region
139 int region = region_index(dphi);
140
141 // If the track is in the toward region, ignore it.
142 if (region == 2) continue;
143
144 // Otherwise, increment the relevant counters
145 ++n_ch[i][region];
146 sumpt[i][region] += t.pT();
147
148 }
149 // Calculate the pT_avg for the away and transverse regions.
150 // (And make sure we don't try to divide by zero.)
151 ptavg[i][0] = (n_ch[i][0] == 0 ? 0.0 : sumpt[i][0] / n_ch[i][0]);
152 ptavg[i][1] = (n_ch[i][1] == 0 ? 0.0 : sumpt[i][1] / n_ch[i][1]);
153
154 _nEvents[i]->fill();
155 }
156
157 fillProfiles(_h_meanNch, n_ch, lead_jet_pts, 1.0 / (2*PI));
158 fillProfiles(_h_meanPtAvg, ptavg, lead_jet_pts, 1.0);
159 fillProfiles(_h_meanPtSum, sumpt, lead_jet_pts, 1.0 / (2*PI));
160
161 fillHistograms(_h_Nch, n_ch, lead_jet_pts);
162 fillHistograms(_h_PtAvg, ptavg, lead_jet_pts);
163 fillHistograms(_h_PtSum, sumpt, lead_jet_pts);
164 }
165
166
167 void fillProfiles(Profile1DPtr plots[5][2], double var[5][2], double lead_pt[5], double scale) {
168 for (int i=0; i<5; ++i) {
169 double pt = lead_pt[i];
170 for (int j=0; j<2; ++j) {
171 double v = var[i][j];
172 plots[i][j]->fill(pt, v*scale);
173 }
174 }
175 }
176
177
178 void fillHistograms(BinnedHistogram plots[5][2], double var[5][2], double lead_pt[5]) {
179 for (int i=0; i<5; ++i) {
180 double pt = lead_pt[i];
181 for (int j=0; j<2; ++j) {
182 double v = var[i][j];
183 plots[i][j].fill(pt, v);
184 }
185 }
186 }
187
188
189 int region_index(double dphi) {
190 assert(inRange(dphi, 0.0, PI, CLOSED, CLOSED));
191 if (dphi < PI/3.0) return 2;
192 if (dphi < 2*PI/3.0) return 1;
193 return 0;
194 }
195
196
197 /// Normalise histograms etc., after the run
198 void finalize() {
199 finalizeHistograms(_h_Nch);
200 finalizeHistograms(_h_PtAvg);
201 finalizeHistograms(_h_PtSum);
202 }
203
204
205 void finalizeHistograms(BinnedHistogram plots[5][2]) {
206 for (int i = 0; i < 5; ++i) {
207 for (int j = 0; j < 2; ++j) {
208 vector<Histo1DPtr> histos = plots[i][j].histos();
209 for(Histo1DPtr h : histos) {
210 scale(h, 1.0/ *_nEvents[i]);
211 }
212 }
213 }
214 }
215
216 //@}
217
218
219 private:
220
221 // Data members like post-cuts event weight counters go here
222 CounterPtr _nEvents[5];
223
224 Profile1DPtr _h_meanNch[5][2];
225 Profile1DPtr _h_meanPtAvg[5][2];
226 Profile1DPtr _h_meanPtSum[5][2];
227
228 BinnedHistogram _h_Nch[5][2];
229 BinnedHistogram _h_PtAvg[5][2];
230 BinnedHistogram _h_PtSum[5][2];
231
232 };
233
234
235
236 // The hook for the plugin system
237 RIVET_DECLARE_PLUGIN(ATLAS_2012_I1125575);
238
239}
|