Rivet analyses referenceATLAS_2012_I1183818Pseudorapidity dependence of the total transverse energy at 7 TeVExperiment: ATLAS (LHC 7TeV) Inspire ID: 1183818 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
The measurement of the sum of the transverse energy of particles as a function of particle pseudorapidity, eta, in proton-proton collisions at a centre-of-mass energy, $\sqrt{s} = 7 \text{TeV}$ using the ATLAS detector at the Large Hadron Collider. The measurements are performed in the region $|\eta| < 4.8$ for two event classes: those requiring the presence of particles with a low transverse momentum and those requiring particles with a significant transverse momentum (dijet events where both jets have $E_T > 20$ GeV). In the second dataset measurements are made in the region transverse to the hard scatter. Source code: ATLAS_2012_I1183818.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
7
8/// @author Peter Wijeratne <paw@hep.ucl.ac.uk>
9/// @author Robindra Prabhu <prabhu@cern.ch>
10namespace Rivet {
11
12 // A very basic analysis sensitive to ET flow in minbias and dijet events
13 class ATLAS_2012_I1183818 : public Analysis {
14
15 public:
16
17 ATLAS_2012_I1183818()
18 : Analysis("ATLAS_2012_I1183818")
19 {}
20
21
22 public:
23
24 void init() {
25
26 const FinalState cnfs((Cuts::etaIn(-4.8, 4.8)));
27 const ChargedFinalState cfs((Cuts::etaIn(-2.5, 2.5) && Cuts::pT >= 250*MeV));
28 declare(cnfs, "FS");
29 declare(cfs, "CFS");
30
31 const FastJets jetsAntiKt4(cnfs, JetAlg::ANTIKT, 0.4);
32 declare(jetsAntiKt4, "AntiKt4Jets");
33
34 // ------- MINBIAS HISTOGRAMS --------
35 //
36 // MB event counter
37 book(m_chargedEvents, "m_chargedEvents_MB");
38
39 book(_h_ETflowEta ,1, 1, 1);
40 book(_h_SumETbin1 ,3, 1, 1);
41 book(_h_SumETbin2 ,4, 1, 1);
42 book(_h_SumETbin3 ,5, 1, 1);
43 book(_h_SumETbin4 ,6, 1, 1);
44 book(_h_SumETbin5 ,7, 1, 1);
45 book(_h_SumETbin6 ,8, 1, 1);
46
47 // ------- DIJET HISTOGRAMS --------
48 //
49 // Dijet event counter
50 book(m_events_dijets, "m_chargedEvents_dijets");
51
52 // sumET
53 book(_h_transETflowEta , 2, 1, 1);
54 book(_h_transSumETbin1 , 9, 1, 1);
55 book(_h_transSumETbin2 ,10, 1, 1);
56 book(_h_transSumETbin3 ,11, 1, 1);
57 book(_h_transSumETbin4 ,12, 1, 1);
58 book(_h_transSumETbin5 ,13, 1, 1);
59 book(_h_transSumETbin6 ,14, 1, 1);
60
61
62 }
63
64
65 void analyze(const Event& event) {
66
67 const FinalState& cfs = apply<FinalState>(event, "CFS");
68
69 bool isCharged = false;
70 if (cfs.size() >= 2) { // event selection: > 2 charged particles with pT > 250.MeV and |eta| < 2.5
71 isCharged = true;
72 m_chargedEvents->fill();
73 }
74
75 const FinalState& cnfs = apply<FinalState>(event, "FS");
76
77 Particles particles;
78 for( const Particle& p : cnfs.particles() ) {
79 // enforce truth selection representing detected particle sensitivity
80 double pp = p.p3().mod();
81 if (PID::charge3(p.pid()) != 0 && pp < 0.5*GeV) continue;
82 if (PID::charge3(p.pid()) == 0 && pp < 0.2*GeV) continue;
83
84 particles.push_back(p);
85 }
86
87
88 // get jets
89 const FastJets& jetsAntiKt4 = apply<FastJets>(event, "AntiKt4Jets");
90 const Jets& jets = jetsAntiKt4.jetsByPt(Cuts::pT > 20*GeV);
91
92 // initialise sumET variables
93 double sumETbin1 = 0;
94 double sumETbin2 = 0;
95 double sumETbin3 = 0;
96 double sumETbin4 = 0;
97 double sumETbin5 = 0;
98 double sumETbin6 = 0;
99
100 // if (passes event selection)
101 if (isCharged) {
102
103 for( const Particle& p : particles ) {
104
105 ///calculate variables
106 double ET = p.Et()/GeV;
107 double eta = p.abseta();
108
109 // fill histograms
110 _h_ETflowEta->fill(eta, ET);
111
112 if (eta < 0.8) sumETbin1 += ET;
113 else if (eta < 1.6) sumETbin2 += ET;
114 else if (eta < 2.4) sumETbin3 += ET;
115 else if (eta < 3.2) sumETbin4 += ET;
116 else if (eta < 4.0) sumETbin5 += ET;
117 else if (eta <= 4.8) sumETbin6 += ET;
118
119 } // end of for
120
121 _h_SumETbin1->fill(sumETbin1);
122 _h_SumETbin2->fill(sumETbin2);
123 _h_SumETbin3->fill(sumETbin3);
124 _h_SumETbin4->fill(sumETbin4);
125 _h_SumETbin5->fill(sumETbin5);
126 _h_SumETbin6->fill(sumETbin6);
127 }
128
129 // --- do dijet analysis ---
130
131 if ( jets.size() >= 2 && // require at least two jets
132 jets[0].Et() >= 20.*GeV && // require two leading jets to pass ET cuts
133 jets[1].Et() >= 20.*GeV &&
134 fabs(jets[0].eta()) < 2.5 && // require leading jets to be central
135 fabs(jets[1].eta()) < 2.5 &&
136 deltaPhi(jets[0], jets[1]) > 2.5 && // require back-to-back topology
137 jets[1].Et()/jets[0].Et() >= 0.5) { //require ET-balance
138
139 // found an event that satisfies dijet selection, now fill histograms...
140 // initialise dijet sumET variables
141 double trans_sumET_bin1 = 0.;
142 double trans_sumET_bin2 = 0.;
143 double trans_sumET_bin3 = 0.;
144 double trans_sumET_bin4 = 0.;
145 double trans_sumET_bin5 = 0.;
146 double trans_sumET_bin6 = 0.;
147
148 m_events_dijets->fill();
149
150 // loop over all particles and check their relation to leading jet
151 for( const Particle& particle : particles ) {
152
153 // calculate variables
154 double dPhi = deltaPhi( jets[0], particle.momentum() );
155 double ET = particle.Et()/GeV;
156 double eta = fabs(particle.eta());
157
158 // Transverse region
159 if ( dPhi > 1./3.*M_PI && dPhi < 2./3.*M_PI ) {
160 _h_transETflowEta->fill( eta, ET );
161 if (eta < 0.8) { trans_sumET_bin1 += ET; }
162 else if (eta < 1.6) { trans_sumET_bin2 += ET; }
163 else if (eta < 2.4) { trans_sumET_bin3 += ET; }
164 else if (eta < 3.2) { trans_sumET_bin4 += ET; }
165 else if (eta < 4.0) { trans_sumET_bin5 += ET; }
166 else if (eta <= 4.8) { trans_sumET_bin6 += ET; }
167 }
168
169 } // end loop over particles
170
171 _h_transSumETbin1->fill( trans_sumET_bin1);
172 _h_transSumETbin2->fill( trans_sumET_bin2);
173 _h_transSumETbin3->fill( trans_sumET_bin3);
174 _h_transSumETbin4->fill( trans_sumET_bin4);
175 _h_transSumETbin5->fill( trans_sumET_bin5);
176 _h_transSumETbin6->fill( trans_sumET_bin6);
177 } // end of dijet selection cuts
178
179 }
180
181
182 void finalize() {
183 /// several scale factors here:
184 /// 1. nEvents (m_chargedEvents)
185 /// 2. phase-space (2*M_PI)
186 /// 3. double binning due to symmetrisation (2)
187 scale( _h_ETflowEta, 1./m_chargedEvents->val()/(4.*M_PI) );
188 scale( _h_SumETbin1, 1./m_chargedEvents->val() );
189 scale( _h_SumETbin2, 1./m_chargedEvents->val() );
190 scale( _h_SumETbin3, 1./m_chargedEvents->val() );
191 scale( _h_SumETbin4, 1./m_chargedEvents->val() );
192 scale( _h_SumETbin5, 1./m_chargedEvents->val() );
193 scale( _h_SumETbin6, 1./m_chargedEvents->val() );
194
195 //Dijet analysis
196
197 // Dijet scale factors:
198 //1. number of events passing dijet selection
199 //2. phase-space: 1. / 2/3*M_PI
200 //3. double binning due to symmetrisation in |eta| plots : 1/2
201 scale( _h_transETflowEta, 1./m_events_dijets->val() * 1./(4./3.*M_PI) );
202 scale( _h_transSumETbin1, 1./m_events_dijets->val() );
203 scale( _h_transSumETbin2, 1./m_events_dijets->val() );
204 scale( _h_transSumETbin3, 1./m_events_dijets->val() );
205 scale( _h_transSumETbin4, 1./m_events_dijets->val() );
206 scale( _h_transSumETbin5, 1./m_events_dijets->val() );
207 scale( _h_transSumETbin6, 1./m_events_dijets->val() );
208 }
209
210 private:
211
212 // Event counts
213 CounterPtr m_chargedEvents;
214 CounterPtr m_events_dijets;
215
216 // Minbias-analysis: variable + histograms
217 Histo1DPtr _h_ETflowEta;
218 Histo1DPtr _h_SumETbin1;
219 Histo1DPtr _h_SumETbin2;
220 Histo1DPtr _h_SumETbin3;
221 Histo1DPtr _h_SumETbin4;
222 Histo1DPtr _h_SumETbin5;
223 Histo1DPtr _h_SumETbin6;
224
225 // Transverse region
226 Histo1DPtr _h_transETflowEta;
227 Histo1DPtr _h_transSumETbin1;
228 Histo1DPtr _h_transSumETbin2;
229 Histo1DPtr _h_transSumETbin3;
230 Histo1DPtr _h_transSumETbin4;
231 Histo1DPtr _h_transSumETbin5;
232 Histo1DPtr _h_transSumETbin6;
233
234 };
235
236
237 RIVET_DECLARE_PLUGIN(ATLAS_2012_I1183818);
238
239}
|