Rivet analyses referenceATLAS_2012_I1117704High jet multiplicity squark and gluino searchExperiment: ATLAS (LHC) Inspire ID: 1117704 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
Search for SUSY using events with 6 or more jets in association with missing transverse momentum produced in proton-proton collisions at a centre-of-mass energy of 7 TeV. The data sample has a total integrated luminosity of 4.7 fb$^{-1}$. Distributions in the W and top control regions are not produced, while in addition to the plots from the paper the count of events in the different signal regions is included. Source code: ATLAS_2012_I1117704.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/ChargedFinalState.hh"
5#include "Rivet/Projections/VisibleFinalState.hh"
6#include "Rivet/Projections/VetoedFinalState.hh"
7#include "Rivet/Projections/IdentifiedFinalState.hh"
8#include "Rivet/Projections/FastJets.hh"
9#include "Rivet/Tools/RivetMT2.hh"
10
11namespace Rivet {
12
13
14 class ATLAS_2012_I1117704 : public Analysis {
15 public:
16
17 /// @name Constructors etc.
18 //@{
19
20 /// Constructor
21 ATLAS_2012_I1117704()
22 : Analysis("ATLAS_2012_I1117704")
23 { }
24
25 //@}
26
27
28 public:
29
30 /// @name Analysis methods
31 //@{
32
33 /// Book histograms and initialise projections before the run
34 void init() {
35
36 // projection to find the electrons
37 IdentifiedFinalState elecs(Cuts::abseta < 2.47 && Cuts::pT > 20*GeV);
38 elecs.acceptIdPair(PID::ELECTRON);
39 declare(elecs, "elecs");
40
41 // projection to find the muons
42 IdentifiedFinalState muons(Cuts::abseta < 2.4 && Cuts::pT > 10*GeV);
43 muons.acceptIdPair(PID::MUON);
44 declare(muons, "muons");
45
46 // for pTmiss
47 declare(VisibleFinalState(Cuts::abseta < 4.9), "vfs");
48
49 VetoedFinalState vfs;
50 vfs.addVetoPairId(PID::MUON);
51
52 /// Jet finder
53 declare(FastJets(vfs, FastJets::ANTIKT, 0.4), "AntiKtJets04");
54
55 // all tracks (to do deltaR with leptons)
56 declare(ChargedFinalState(Cuts::abseta < 3),"cfs");
57
58 /// Book histograms
59 book(_etmiss_HT_7j55 ,"etmiss_HT_7j55", 8, 0., 16.);
60 book(_etmiss_HT_8j55 ,"etmiss_HT_8j55", 8, 0., 16.);
61 book(_etmiss_HT_9j55 ,"etmiss_HT_9j55", 8, 0., 16.);
62 book(_etmiss_HT_6j80 ,"etmiss_HT_6j80", 8, 0., 16.);
63 book(_etmiss_HT_7j80 ,"etmiss_HT_7j80", 8, 0., 16.);
64 book(_etmiss_HT_8j80 ,"etmiss_HT_8j80", 8, 0., 16.);
65
66 book(_hist_njet55 ,"hist_njet55", 11, 2.5, 13.5);
67 book(_hist_njet80 ,"hist_njet80", 11, 2.5, 13.5);
68
69 book(_count_7j55 ,"count_7j55", 1, 0., 1.);
70 book(_count_8j55 ,"count_8j55", 1, 0., 1.);
71 book(_count_9j55 ,"count_9j55", 1, 0., 1.);
72 book(_count_6j80 ,"count_6j80", 1, 0., 1.);
73 book(_count_7j80 ,"count_7j80", 1, 0., 1.);
74 book(_count_8j80 ,"count_8j80", 1, 0., 1.);
75
76 }
77
78
79 /// Perform the per-event analysis
80 void analyze(const Event& event) {
81 const double weight = 1.0;
82
83 // get the jet candidates
84 Jets cand_jets;
85 for (const Jet& jet :
86 apply<FastJets>(event, "AntiKtJets04").jetsByPt(20.0*GeV) ) {
87 if ( fabs( jet.eta() ) < 2.8 ) {
88 cand_jets.push_back(jet);
89 }
90 }
91
92 // candidate muons
93 Particles cand_mu;
94 Particles chg_tracks =
95 apply<ChargedFinalState>(event, "cfs").particles();
96 for ( const Particle& mu :
97 apply<IdentifiedFinalState>(event, "muons").particlesByPt() ) {
98 double pTinCone = -mu.pT();
99 for ( const Particle& track : chg_tracks ) {
100 if ( deltaR(mu.momentum(),track.momentum()) <= 0.2 )
101 pTinCone += track.pT();
102 }
103 if ( pTinCone < 1.8*GeV )
104 cand_mu.push_back(mu);
105 }
106
107 // candidate electrons
108 Particles cand_e =
109 apply<IdentifiedFinalState>(event, "elecs").particlesByPt();
110
111 // resolve jet/lepton ambiguity
112 Jets recon_jets;
113 for ( const Jet& jet : cand_jets ) {
114 // candidates after |eta| < 2.8
115 if ( fabs( jet.eta() ) >= 2.8 ) continue;
116 bool away_from_e = true;
117 for ( const Particle& e : cand_e ) {
118 if ( deltaR(e.momentum(),jet.momentum()) <= 0.2 ) {
119 away_from_e = false;
120 break;
121 }
122 }
123 if ( away_from_e ) recon_jets.push_back( jet );
124 }
125
126 // only keep electrons more than R=0.4 from jets
127 Particles recon_e;
128 for ( const Particle& e : cand_e ) {
129 bool away = true;
130 for ( const Jet& jet : recon_jets ) {
131 if ( deltaR(e.momentum(),jet.momentum()) < 0.4 ) {
132 away = false;
133 break;
134 }
135 }
136 if ( away )
137 recon_e.push_back( e );
138 }
139
140 // only keep muons more than R=0.4 from jets
141 Particles recon_mu;
142 for ( const Particle& mu : cand_mu ) {
143 bool away = true;
144 for ( const Jet& jet : recon_jets ) {
145 if ( deltaR(mu.momentum(),jet.momentum()) < 0.4 ) {
146 away = false;
147 break;
148 }
149 }
150 if ( away )
151 recon_mu.push_back( mu );
152 }
153
154 // pTmiss
155 Particles vfs_particles =
156 apply<VisibleFinalState>(event, "vfs").particles();
157 FourMomentum pTmiss;
158 for ( const Particle& p : vfs_particles ) {
159 pTmiss -= p.momentum();
160 }
161 double eTmiss = pTmiss.pT();
162
163 // now only use recon_jets, recon_mu, recon_e
164
165 // reject events with electrons and muons
166 if ( ! ( recon_mu.empty() && recon_e.empty() ) ) {
167 MSG_DEBUG("Charged leptons left after selection");
168 vetoEvent;
169 }
170
171 // calculate H_T
172 double HT=0;
173 for ( const Jet& jet : recon_jets ) {
174 if ( jet.pT() > 40 * GeV )
175 HT += jet.pT() ;
176 }
177
178 // number of jets
179 unsigned int njet55=0, njet80=0;
180 for (unsigned int ix=0;ix<recon_jets.size();++ix) {
181 if(recon_jets[ix].pT()>80.*GeV) ++njet80;
182 if(recon_jets[ix].pT()>55.*GeV) ++njet55;
183 }
184
185 if(njet55==0) vetoEvent;
186
187 double ratio = eTmiss/sqrt(HT);
188
189 if(ratio>4.) {
190 _hist_njet55->fill(njet55,weight);
191 _hist_njet80->fill(njet80,weight);
192 // 7j55
193 if(njet55>=7)
194 _count_7j55->fill( 0.5, weight);
195 // 8j55
196 if(njet55>=8)
197 _count_8j55->fill( 0.5, weight) ;
198 // 8j55
199 if(njet55>=9)
200 _count_9j55->fill( 0.5, weight) ;
201 // 6j80
202 if(njet80>=6)
203 _count_6j80->fill( 0.5, weight) ;
204 // 7j80
205 if(njet80>=7)
206 _count_7j80->fill( 0.5, weight) ;
207 // 8j80
208 if(njet80>=8)
209 _count_8j80->fill( 0.5, weight) ;
210 }
211
212 if(njet55>=7)
213 _etmiss_HT_7j55->fill( ratio, weight);
214 // 8j55
215 if(njet55>=8)
216 _etmiss_HT_8j55->fill( ratio, weight) ;
217 // 8j55
218 if(njet55>=9)
219 _etmiss_HT_9j55->fill( ratio, weight) ;
220 // 6j80
221 if(njet80>=6)
222 _etmiss_HT_6j80->fill( ratio, weight) ;
223 // 7j80
224 if(njet80>=7)
225 _etmiss_HT_7j80->fill( ratio, weight) ;
226 // 8j80
227 if(njet80>=8)
228 _etmiss_HT_8j80->fill( ratio, weight) ;
229
230 }
231
232 //@}
233
234 void finalize() {
235 double norm = crossSection()/femtobarn*4.7/sumOfWeights();
236
237 scale(_etmiss_HT_7j55,2.*norm);
238 scale(_etmiss_HT_8j55,2.*norm);
239 scale(_etmiss_HT_9j55,2.*norm);
240 scale(_etmiss_HT_6j80,2.*norm);
241 scale(_etmiss_HT_7j80,2.*norm);
242 scale(_etmiss_HT_8j80,2.*norm);
243
244 scale(_hist_njet55,norm);
245 scale(_hist_njet80,norm);
246
247 scale(_count_7j55,norm);
248 scale(_count_8j55,norm);
249 scale(_count_9j55,norm);
250 scale(_count_6j80,norm);
251 scale(_count_7j80,norm);
252 scale(_count_8j80,norm);
253 }
254
255 private:
256
257 /// @name Histograms
258 //@{
259 Histo1DPtr _etmiss_HT_7j55;
260 Histo1DPtr _etmiss_HT_8j55;
261 Histo1DPtr _etmiss_HT_9j55;
262 Histo1DPtr _etmiss_HT_6j80;
263 Histo1DPtr _etmiss_HT_7j80;
264 Histo1DPtr _etmiss_HT_8j80;
265
266 Histo1DPtr _hist_njet55;
267 Histo1DPtr _hist_njet80;
268
269 Histo1DPtr _count_7j55;
270 Histo1DPtr _count_8j55;
271 Histo1DPtr _count_9j55;
272 Histo1DPtr _count_6j80;
273 Histo1DPtr _count_7j80;
274 Histo1DPtr _count_8j80;
275 //@}
276
277 };
278
279 // The hook for the plugin system
280 RIVET_DECLARE_PLUGIN(ATLAS_2012_I1117704);
281
282}
|