Rivet analyses referenceATLAS_2016_CONF_2016_054ATLAS 2016 1-lepton SUSY search at 13 \text{TeV}, from 14.8/fb CONF noteExperiment: ATLAS (LHC) Status: UNVALIDATED Authors:
Beams: p+ p+ Beam energies: (6500.0, 6500.0) GeV Run details:
A search for squarks and gluinos in final states with an isolated electron or muon, multiple jets and large missing transverse momentum using proton--proton collision data at a centre-of-mass energy of $\sqrt{s} = 13 \text{TeV}$. The dataset corresponds to an integrated luminosity of 14.8/fb, recorded in 2015 and 2016 by the ATLAS experiment at the Large Hadron Collider. Source code: ATLAS_2016_CONF_2016_054.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/PromptFinalState.hh"
5#include "Rivet/Projections/FastJets.hh"
6#include "Rivet/Projections/Sphericity.hh"
7#include "Rivet/Projections/SmearedParticles.hh"
8#include "Rivet/Projections/SmearedJets.hh"
9#include "Rivet/Projections/SmearedMET.hh"
10
11namespace Rivet {
12
13
14 /// @brief ATLAS 2016 1-lepton SUSY search, from 14.8/fb CONF note
15 class ATLAS_2016_CONF_2016_054 : public Analysis {
16 public:
17
18 /// Constructor
19 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_CONF_2016_054);
20
21
22 /// @name Analysis methods
23 /// @{
24
25 /// Book histograms and initialise projections before the run
26 void init() {
27
28 // Initialise and register projections
29 FinalState calofs(Cuts::abseta < 4.9);
30 FastJets fj(calofs, JetAlg::ANTIKT, 0.4);
31 declare(fj, "TruthJets");
32 declare(SmearedJets(fj, JET_SMEAR_ATLAS_RUN2, [](const Jet& j) {
33 if (j.abseta() > 2.5) return 0.;
34 return j.bTagged(Cuts::pT > 5*GeV) ? 0.77 :
35 j.cTagged(Cuts::pT > 5*GeV) ? 1/6.2 : 1/134.; }), "Jets");
36
37 MissingMomentum mm(calofs);
38 declare(mm, "TruthMET");
39 declare(SmearedMET(mm, MET_SMEAR_ATLAS_RUN2), "MET");
40
41 FinalState es(Cuts::abseta < 2.47 && Cuts::pT > 7*GeV && Cuts::abspid == PID::ELECTRON);
42 declare(es, "TruthElectrons");
43 declare(SmearedParticles(es, ELECTRON_RECOEFF_ATLAS_RUN2, ELECTRON_SMEAR_ATLAS_RUN2), "Electrons");
44
45 FinalState mus(Cuts::abseta < 2.5 && Cuts::pT > 6*GeV && Cuts::abspid == PID::MUON);
46 declare(mus, "TruthMuons");
47 declare(SmearedParticles(mus, MUON_EFF_ATLAS_RUN2, MUON_SMEAR_ATLAS_RUN2), "Muons");
48
49
50 // Book histograms/counters
51 book(_h_gg2j,"GG-2j");
52 book(_h_gg6j0,"GG-6j-0bulk");
53 book(_h_gg6j1,"GG-6j-1highmass");
54 book(_h_gg4j0,"GG-4j-0lowx");
55 book(_h_gg4j1,"GG-4j-1lowxbveto");
56 book(_h_gg4j2,"GG-4j-2highx");
57 book(_h_ss4j0,"SS-4j-0x12");
58 book(_h_ss4j1,"SS-4j-1lowx");
59 book(_h_ss5j0,"SS-5j-0x12");
60 book(_h_ss5j1,"SS-5j-1highx");
61
62 }
63
64
65 /// Perform the per-event analysis
66 void analyze(const Event& event) {
67
68 // Get baseline electrons, muons, and jets
69 Particles elecs = apply<ParticleFinder>(event, "Electrons").particles();
70 Particles muons = apply<ParticleFinder>(event, "Muons").particles();
71 Jets jets = apply<JetFinder>(event, "Jets").jetsByPt(Cuts::pT > 20*GeV && Cuts::abseta < 4.5);
72
73 // Jet/electron/muons overlap removal and selection
74 // Remove any jet within dR = 0.2 of an electron
75 for (const Particle& e : elecs)
76 idiscard(jets, deltaRLess(e, 0.2, RAPIDITY));
77 // Remove any electron within dR = 0.01 of a muon
78 for (const Particle& m : muons)
79 idiscard(elecs, deltaRLess(m, 0.01, RAPIDITY));
80 // Assemble b-jets collection, and remove muons within dR = 0.2 of a b-tagged jet
81 Jets bjets;
82 for (const Jet& j : jets) {
83 if (j.abseta() < 2.5 && j.pT() > 30*GeV && j.bTagged(Cuts::pT > 5*GeV)) {
84 bjets += j;
85 idiscard(muons, deltaRLess(j, 0.2, RAPIDITY));
86 }
87 }
88 // Remove any jet within dR = 0.2 of a muon if track conditions are met
89 for (const Particle& m : muons)
90 idiscard(jets, [&](const Jet& j){
91 if (deltaR(j,m) > 0.2) return false;
92 /// @todo Add track efficiency random filtering
93 const Particles trks = j.particles(Cuts::abscharge > 0 && Cuts::pT > 0.5*GeV);
94 return trks.size() < 4 || m.pT()/j.pT() > 0.7;
95 });
96 // Remove any muon within dR = 0.2 of a remaining jet if the same track conditions are *not* met
97 /// @todo There must be nicer way to do complementary removal...
98 for (const Jet& j : jets) {
99 /// @todo Add track efficiency random filtering
100 const size_t ntrks = j.particles(Cuts::abscharge > 0 && Cuts::pT > 0.5*GeV).size();
101 idiscard(muons, [&](const Particle& m){
102 if (deltaR(j,m) > 0.2) return false;
103 return ntrks > 3 && m.pT()/j.pT() < 0.7;
104 });
105 }
106 // Remove any muon with dR close to a remaining jet, via a functional form
107 for (const Jet& j : jets)
108 idiscard(muons, [&](const Particle& m) { return deltaR(m,j, RAPIDITY) < min(0.4, 0.04 + 10*GeV/m.pT()); });
109
110
111 // Signal jet selection
112 const Jets sigjets = select(jets, Cuts::pT > 30*GeV && Cuts::abseta < 2.8);
113 const Jets sigbjets = bjets;
114
115 // "Gradient-loose" signal lepton selection
116 const ParticleEffFilter grad_loose_filter([](const Particle& e) { return e.pT() > 60*GeV ? 0.98 : 0.95; });
117 Particles sigelecs = select(elecs, grad_loose_filter);
118 Particles sigmuons = select(muons, grad_loose_filter);
119 // Tight electron selection (NB. assuming independent eff to gradient-loose... hmm)
120 iselect(sigelecs, ParticleEffFilter(ELECTRON_EFF_ATLAS_RUN2_TIGHT));
121
122
123 // MET calculation (NB. done generically, with smearing, rather than via explicit physics objects)
124 const Vector3 vmet = -apply<SmearedMET>(event, "MET").vectorEt();
125 const double etmiss = vmet.mod();
126
127
128 //////////////////
129
130
131 // Event selection cuts
132 if (sigelecs.size() + sigmuons.size() != 1) vetoEvent;
133 const Particle siglepton = sigelecs.empty() ? sigmuons.front() : sigelecs.front();
134
135 // mT and m_eff
136 const double mT = sqrt(2*siglepton.pT()*etmiss*(1-cos(deltaPhi(siglepton,vmet))));
137 const double meff = siglepton.pT() + sum(sigjets, pT, 0.0) + etmiss;
138
139 // Aplanarities
140 Sphericity sph;
141 vector<FourMomentum> vecs;
142 transform(sigjets, vecs, mom);
143 sph.calc(vecs);
144 const double jet_aplanarity = sph.aplanarity();
145 vecs += siglepton.mom();
146 sph.calc(vecs);
147 const double lepton_aplanarity = sph.aplanarity();
148
149
150 //////////////////
151
152
153 // Fill counters
154 // GG
155 if (siglepton.pT() < 35*GeV && sigjets.size() >= 2 &&
156 sigjets[0].pT() > 200*GeV && sigjets[1].pT() > 30*GeV &&
157 mT > 100*GeV && etmiss > 460*GeV && etmiss/meff > 0.35) _h_gg2j->fill();
158 if (siglepton.pT() > 35*GeV && sigjets.size() >= 6 &&
159 sigjets[0].pT() > 125*GeV && sigjets[5].pT() > 30*GeV &&
160 mT > 225*GeV && etmiss > 250*GeV && meff > 1000*GeV && etmiss/meff > 0.2 &&
161 jet_aplanarity > 0.04) _h_gg6j0->fill();
162 if (siglepton.pT() > 35*GeV && sigjets.size() >= 6 &&
163 sigjets[0].pT() > 125*GeV && sigjets[5].pT() > 30*GeV &&
164 mT > 225*GeV && etmiss > 250*GeV && meff > 2000*GeV && etmiss/meff > 0.1 &&
165 jet_aplanarity > 0.04) _h_gg6j1->fill();
166 if (sigjets.size() >= 4 && sigjets[3].pT() > 100*GeV &&
167 mT > 125*GeV && etmiss > 250*GeV && meff > 2000*GeV && jet_aplanarity > 0.06) _h_gg4j0->fill();
168 if (sigjets.size() >= 4 && sigjets[3].pT() > 100*GeV && sigbjets.empty() &&
169 mT > 125*GeV && etmiss > 250*GeV && meff > 2000*GeV && jet_aplanarity > 0.03) _h_gg4j1->fill();
170 if (siglepton.pT() > 35*GeV &&
171 sigjets.size() >= 4 && sigjets[0].pT() > 400*GeV && inRange(sigjets[3].pT(), 30*GeV, 100*GeV) &&
172 mT > 475*GeV && etmiss > 250*GeV && meff > 1600*GeV && etmiss/meff > 0.3) _h_gg4j2->fill();
173 // SS
174 if (siglepton.pT() > 35*GeV && sigjets.size() >= 4 && sigjets[3].pT() > 50*GeV &&
175 mT > 175*GeV && etmiss > 300*GeV && meff > 1200*GeV && lepton_aplanarity > 0.08) _h_ss4j0->fill();
176 if (siglepton.pT() > 35*GeV && sigjets.size() >= 5 && sigjets[4].pT() > 50*GeV && sigbjets.empty() &&
177 mT > 175*GeV && etmiss > 300*GeV && etmiss/meff > 0.2) _h_ss5j0->fill();
178 if (siglepton.pT() > 35*GeV && sigjets.size() >= 4 && sigjets[0].pT() > 250*GeV && sigjets[3].pT() > 30*GeV &&
179 inRange(mT, 150*GeV, 400*GeV) && etmiss > 250*GeV && lepton_aplanarity > 0.03) _h_ss4j1->fill();
180 if (siglepton.pT() > 35*GeV && sigjets.size() >= 5 && sigjets[4].pT() > 30*GeV &&
181 mT > 400*GeV && etmiss > 400*GeV && lepton_aplanarity > 0.03) _h_ss5j1->fill();
182
183 }
184
185
186 /// Normalise counters after the run
187 void finalize() {
188
189 const double sf = 14.8*crossSection()/femtobarn/sumOfWeights();
190 scale(_h_gg2j, sf); scale(_h_gg6j0, sf); scale(_h_gg6j1, sf);
191 scale(_h_gg4j0, sf); scale(_h_gg4j1, sf); scale(_h_gg4j2, sf);
192 scale(_h_ss4j0, sf); scale(_h_ss4j1, sf); scale(_h_ss5j0, sf);
193 scale(_h_ss5j1, sf);
194
195 }
196
197 /// @}
198
199
200 private:
201
202 /// @name Histograms
203 /// @{
204 CounterPtr _h_gg2j, _h_gg6j0, _h_gg6j1, _h_gg4j0, _h_gg4j1, _h_gg4j2;
205 CounterPtr _h_ss4j0, _h_ss4j1, _h_ss5j0,_h_ss5j1;
206 /// @}
207
208
209 };
210
211
212
213 RIVET_DECLARE_PLUGIN(ATLAS_2016_CONF_2016_054);
214
215
216}
|