rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2016_CONF_2016_054

ATLAS 2016 1-lepton SUSY search at 13 \text{TeV}, from 14.8/fb CONF note
Experiment: ATLAS (LHC)
Status: UNVALIDATED
Authors:
  • Andy Buckley
No references listed
Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • BSM signal events

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