Rivet analyses referenceATLAS_2011_CONF_2011_098b-jets search for supersymmetry with 0-leptonsExperiment: ATLAS (LHC) Status: OBSOLETE Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
Search for supersymmmetric particles by ATLAS at 7 TeV in events with b-jets, large missing energy, and no leptons. Event counts in four signal regions (1 b-jet, $m_eff>500$\,GeV; 1 b-jet, $m_eff>700$\,GeV; 2 b-jets, $m_eff>500$\,GeV; 2 b-jets, $m_eff>700$\,GeV) are implemented as one-bin histograms. Histograms for missing transverse energy, effective mass, and pT of the leading jet are implemented for the 1 b-tag and 2 b-tag signal regions. Source code: ATLAS_2011_CONF_2011_098.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/IdentifiedFinalState.hh"
7#include "Rivet/Projections/FastJets.hh"
8#include "Rivet/Tools/Random.hh"
9
10namespace Rivet {
11
12
13 class ATLAS_2011_CONF_2011_098 : public Analysis {
14 public:
15
16 /// Constructor
17 ATLAS_2011_CONF_2011_098()
18 : Analysis("ATLAS_2011_CONF_2011_098"),
19 //debug variables
20 threeJA(0), threeJB(0), threeJC(0), threeJD(0), bj(0), jets(0), zerolept(0), eTmisscut(0)
21 { }
22
23
24 /// @name Analysis methods
25 /// @{
26
27 /// Book histograms and initialise projections before the run
28 void init() {
29
30 // projection to find the electrons
31 IdentifiedFinalState elecs(Cuts::abseta < 2.47 && Cuts::pT > 20*GeV);
32 elecs.acceptIdPair(PID::ELECTRON);
33 declare(elecs, "elecs");
34
35
36 // projection to find the muons
37 IdentifiedFinalState muons(Cuts::abseta < 2.4 && Cuts::pT > 10*GeV);
38 muons.acceptIdPair(PID::MUON);
39 declare(muons, "muons");
40
41 /// Jet finder
42 declare(FastJets(FinalState(), JetAlg::ANTIKT, 0.4), "AntiKtJets04");
43
44
45 // all tracks (to do deltaR with leptons)
46 declare(ChargedFinalState(Cuts::abseta < 3.0),"cfs");
47
48 // for pTmiss
49 declare(VisibleFinalState(Cuts::abseta < 4.9),"vfs");
50
51
52 /// Book histograms
53 book(_count_threeJA ,"count_threeJA", 1, 0., 1.);
54 book(_count_threeJB ,"count_threeJB", 1, 0., 1.);
55 book(_count_threeJC ,"count_threeJC", 1, 0., 1.);
56 book(_count_threeJD ,"count_threeJD", 1, 0., 1.);
57 book(_hist_meff_1bjet ,"meff_1bjet", 32, 0., 1600.);
58 book(_hist_eTmiss_1bjet ,"eTmiss_1bjet", 6, 0., 600.);
59 book(_hist_pTj_1bjet ,"pTjet_1bjet", 20, 0., 800.);
60 book(_hist_meff_2bjet ,"meff_2bjet", 32, 0., 1600.);
61 book(_hist_eTmiss_2bjet ,"eTmiss_2bjet", 6, 0., 600.);
62 book(_hist_pTj_2bjet ,"pTjet_2bjet", 20, 0., 800.);
63
64
65 }
66
67
68 /// Perform the per-event analysis
69 void analyze(const Event& event) {
70
71
72 // Temp: calorimeter module failure with 10% acceptance loss;
73 // region unknown ==> randomly choose 10% of events to be vetoed
74
75 if ( rand01() < 0.1 ) vetoEvent;
76
77 Jets tmp_cand_jets = apply<FastJets>(event, "AntiKtJets04").jetsByPt(Cuts::pT>20*GeV && Cuts::abseta < 2.8);
78
79 Particles cand_e =
80 apply<IdentifiedFinalState>(event, "elecs").particlesByPt();
81 Particles cand_mu =
82 apply<IdentifiedFinalState>(event, "muons").particlesByPt();
83 Particles chg_tracks =
84 apply<ChargedFinalState>(event, "cfs").particles();
85
86//cerr << "cand_e.size(): " << cand_e.size() << " cand_mu.size(): " << cand_mu.size() << '\n';
87
88
89 Jets cand_jets;
90 for ( const Jet& jet : tmp_cand_jets ) {
91 if ( fabs( jet.eta() ) >= 2.8 )
92 cand_jets.push_back( jet );
93 else {
94 bool away_from_e = true;
95 for ( const Particle & e : cand_e ) {
96 if ( deltaR(e.momentum(),jet.momentum()) <= 0.2 ) {
97 away_from_e = false;
98 break;
99 }
100 }
101 if ( away_from_e )
102 cand_jets.push_back( jet );
103 }
104 }
105
106 Particles cand_lept;
107
108 bool isolated_e;
109 for ( const Particle & e : cand_e ) {
110 isolated_e = true;
111 for ( const Jet& jet : cand_jets ) {
112 if ( deltaR(e.momentum(),jet.momentum()) < 0.4 )
113 isolated_e = false;
114 }
115 if ( isolated_e == true )
116 cand_lept.push_back( e );
117 }
118
119
120 bool isolated_mu;
121 for ( const Particle & mu : cand_mu ) {
122 isolated_mu = true;
123 for ( const Jet& jet : cand_jets ) {
124 if ( deltaR(mu.momentum(),jet.momentum()) < 0.4 )
125 isolated_mu = false;
126 }
127 if ( isolated_mu == true)
128 cand_lept.push_back( mu );
129 }
130
131
132 // pTmiss
133 Particles vfs_particles
134 = apply<VisibleFinalState>(event, "vfs").particles();
135 FourMomentum pTmiss;
136 for ( const Particle & p : vfs_particles ) {
137 pTmiss -= p.momentum();
138 }
139 double eTmiss = pTmiss.pT();
140
141
142 // bjets
143 Jets bjets,recon_jets;
144 for (const Jet& j : cand_jets) {
145 if(fabs( j.eta() ) <= 2.8) {
146 recon_jets.push_back(j);
147 if ( fabs( j.eta() ) <= 2.5 && j.perp()>50. &&
148 j.bTagged() && rand01() < 0.5 )
149 bjets.push_back(j);
150 }
151 }
152
153 if (bjets.empty()) {
154 MSG_DEBUG("No b-jet axes in acceptance");
155 vetoEvent;
156 }
157
158 ++bj;
159
160
161
162 // Jets event selection
163 if ( recon_jets.size() < 3 )
164 vetoEvent;
165 if ( recon_jets[0].pT() <= 130*GeV )
166 vetoEvent;
167 if ( recon_jets[1].pT() <= 50*GeV ||
168 recon_jets[2].pT() <= 50*GeV )
169 vetoEvent;
170 ++jets;
171
172 // eTmiss cut
173 if ( eTmiss <= 130*GeV )
174 vetoEvent;
175
176 ++eTmisscut;
177
178 // 0-lepton requirement
179 if ( !cand_lept.empty() )
180 vetoEvent;
181 ++zerolept;
182
183 // m_eff cut
184 double m_eff = eTmiss
185 + recon_jets[0].pT()
186 + recon_jets[1].pT()
187 + recon_jets[2].pT();
188
189 if ( eTmiss / m_eff <= 0.25 )
190 vetoEvent;
191
192
193 // min_dPhi
194 double min_dPhi = 999.999;
195 for ( int i = 0; i < 3; ++i ) {
196 double dPhi = deltaPhi( pTmiss.phi(), recon_jets[i].phi() );
197 min_dPhi = min( min_dPhi, dPhi );
198 }
199
200 if ( min_dPhi <= 0.4 )
201 vetoEvent;
202
203
204
205 // ==================== FILL ====================
206
207
208 // 1 bjet
209 if ( bjets.size() >= 1 ) {
210
211 _hist_meff_1bjet->fill(m_eff);
212 _hist_eTmiss_1bjet->fill(eTmiss);
213 _hist_pTj_1bjet->fill(recon_jets[0].pT());
214
215 // 3JA region
216 if ( m_eff > 200*GeV ) {
217 ++threeJA;
218 _count_threeJA->fill(0.5);
219 }
220
221 // 3JB region
222 if ( m_eff > 700*GeV ) {
223 ++threeJB;
224 _count_threeJB->fill(0.5);
225 }
226 }
227
228 // 2 bjets
229 if ( bjets.size() >= 2 ) {
230
231 _hist_meff_2bjet->fill(m_eff);
232 _hist_eTmiss_2bjet->fill(eTmiss);
233 _hist_pTj_2bjet->fill(recon_jets[0].pT());
234
235 // 3JC region
236 if ( m_eff > 500*GeV ) {
237 ++threeJC;
238 _count_threeJC->fill(0.5);
239 }
240
241 // 3JD region
242 if ( m_eff > 700*GeV ) {
243 ++threeJD;
244 _count_threeJD->fill(0.5);
245 }
246 }
247
248
249
250
251 }
252
253 /// @}
254
255
256 void finalize() {
257 scale( _hist_meff_1bjet, 50. * 830. * crossSection()/picobarn/sumOfWeights() );
258 scale( _hist_eTmiss_1bjet, 100. * 830. * crossSection()/picobarn/sumOfWeights() );
259 scale( _hist_pTj_1bjet, 40. * 830. * crossSection()/picobarn/sumOfWeights() );
260 scale( _hist_meff_2bjet, 50. * 830. * crossSection()/picobarn/sumOfWeights() );
261 scale( _hist_eTmiss_2bjet, 100. * 830. * crossSection()/picobarn/sumOfWeights() );
262 scale( _hist_pTj_2bjet, 40. * 830. * crossSection()/picobarn/sumOfWeights() );
263
264// cerr<< '\n'<<'\n'
265// << "Saw "
266// << bj << " events aft bjets cut, "
267// << jets << " events aft jet cuts, "
268// << eTmisscut << " events aft eTmiss cut, "
269// << zerolept << " events after 0-lept cut. "
270// << '\n'
271// << threeJA << " 3JA events, "
272// << threeJB << " 3JB events, "
273// << threeJC << " 3JC events, "
274// << threeJD << " 3JD events. "
275// << '\n'
276// ;
277
278 }
279
280
281 private:
282
283 /// @name Histograms
284 /// @{
285 Histo1DPtr _count_threeJA;
286 Histo1DPtr _count_threeJB;
287 Histo1DPtr _count_threeJC;
288 Histo1DPtr _count_threeJD;
289 Histo1DPtr _hist_meff_1bjet;
290 Histo1DPtr _hist_eTmiss_1bjet;
291 Histo1DPtr _hist_pTj_1bjet;
292 Histo1DPtr _hist_meff_2bjet;
293 Histo1DPtr _hist_eTmiss_2bjet;
294 Histo1DPtr _hist_pTj_2bjet;
295
296 /// @}
297
298
299// debug variables
300int threeJA;
301int threeJB;
302int threeJC;
303int threeJD;
304int bj;
305int jets;
306int zerolept;
307int eTmisscut;
308
309 };
310
311
312
313 RIVET_DECLARE_PLUGIN(ATLAS_2011_CONF_2011_098);
314
315}
|