rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2016_CONF_2016_078

ATLAS ICHEP16 0-lepton SUSY search at 13 \text{TeV} with 13.2/fb
Experiment: ATLAS (LHC)
Status: VALIDATED
Authors:
  • Andy Buckley
No references listed
Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • BSM signal events.

ATLAS search for SUSY in 13 TeV $pp$ collisions at LHC Run 2, using 13.2/fb of integrated luminosity and events containing missing transverse momentum and no isolated high-energy leptons. Detailed info: https://atlas.web.cern.ch/Atlas/GROUPS/PHYSICS/CONFNOTES/ATLAS-CONF-2016-078/

Source code: ATLAS_2016_CONF_2016_078.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 0-lepton SUSY search, from 13/fb ICHEP'16 CONF note
 16  class ATLAS_2016_CONF_2016_078 : public Analysis {
 17  public:
 18
 19    /// Constructor
 20    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_CONF_2016_078);
 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 < 3.2);
 31      FastJets fj(calofs, FastJets::ANTIKT, 0.4);
 32      declare(fj, "TruthJets");
 33      declare(SmearedJets(fj, JET_SMEAR_ATLAS_RUN2, //JET_BTAG_ATLAS_RUN2_MV2C10
 34                          [](const Jet& j) {
 35                            if (j.abseta() > 2.5) return 0.;
 36                            return j.bTagged(Cuts::pT > 5*GeV) ? 0.77 : j.cTagged(Cuts::pT > 5*GeV) ? 1/6. : 1/134.;
 37                          }), "RecoJets");
 38
 39      MissingMomentum mm(calofs);
 40      declare(mm, "TruthMET");
 41      declare(SmearedMET(mm, MET_SMEAR_ATLAS_RUN2), "RecoMET");
 42
 43      PromptFinalState es(Cuts::abseta < 2.47 && Cuts::abspid == PID::ELECTRON, true, true);
 44      declare(es, "TruthElectrons");
 45      declare(SmearedParticles(es, ELECTRON_RECOEFF_ATLAS_RUN2, ELECTRON_SMEAR_ATLAS_RUN2), "RecoElectrons");
 46
 47      PromptFinalState mus(Cuts::abseta < 2.7 && Cuts::abspid == PID::MUON, true);
 48      declare(mus, "TruthMuons");
 49      declare(SmearedParticles(mus, MUON_EFF_ATLAS_RUN2, MUON_SMEAR_ATLAS_RUN2), "RecoMuons");
 50
 51
 52      // Book histograms/counters
 53      book(_h_2j_0800,"2j-0800");
 54      book(_h_2j_1200,"2j-1200");
 55      book(_h_2j_1600,"2j-1600");
 56      book(_h_2j_2000,"2j-2000");
 57      book(_h_3j_1200,"2j-2000");
 58      book(_h_4j_1000,"4j-1000");
 59      book(_h_4j_1400,"4j-1400");
 60      book(_h_4j_1800,"4j-1800");
 61      book(_h_4j_2200,"4j-2200");
 62      book(_h_4j_2600,"4j-2600");
 63      book(_h_5j_1400,"5j-1400");
 64      book(_h_6j_1800,"6j-1800");
 65      book(_h_6j_2200,"6j-2200");
 66
 67
 68      // Book cut-flows
 69      const vector<string> cuts23j = {"Pre-sel+MET+pT1+meff", "Njet", "Dphi_min(j123,MET)", "Dphi_min(j4+,MET)", "pT2", "eta_j12", "MET/sqrtHT", "m_eff(incl)"};
 70      _flows.addCutflow("2j-0800", cuts23j);
 71      _flows.addCutflow("2j-1200", cuts23j);
 72      _flows.addCutflow("2j-1600", cuts23j);
 73      _flows.addCutflow("2j-2000", cuts23j);
 74      _flows.addCutflow("3j-1200", cuts23j);
 75      const vector<string> cuts456j = {"Pre-sel+MET+pT1+meff", "Njet", "Dphi_min(j123,MET)", "Dphi_min(j4+,MET)", "pT4", "eta_j1234", "Aplanarity", "MET/m_eff(Nj)", "m_eff(incl)"};
 76      _flows.addCutflow("4j-1000", cuts456j);
 77      _flows.addCutflow("4j-1400", cuts456j);
 78      _flows.addCutflow("4j-1800", cuts456j);
 79      _flows.addCutflow("4j-2200", cuts456j);
 80      _flows.addCutflow("4j-2600", cuts456j);
 81      _flows.addCutflow("5j-1400", cuts456j);
 82      _flows.addCutflow("6j-1800", cuts456j);
 83      _flows.addCutflow("6j-2200", cuts456j);
 84
 85    }
 86
 87
 88    /// Perform the per-event analysis
 89    void analyze(const Event& event) {
 90
 91      _flows.fillinit();
 92
 93      // Same MET cut for all signal regions
 94      const Vector3 vmet = -apply<SmearedMET>(event, "RecoMET").vectorEt();
 95      const double met = vmet.mod();
 96      if (met < 250*GeV) vetoEvent;
 97
 98      // Get baseline electrons, muons, and jets
 99      Particles elecs = apply<ParticleFinder>(event, "RecoElectrons").particles(Cuts::pT > 10*GeV);
100      Particles muons = apply<ParticleFinder>(event, "RecoMuons").particles(Cuts::pT > 10*GeV);
101      Jets jets = apply<JetAlg>(event, "RecoJets").jetsByPt(Cuts::pT > 20*GeV && Cuts::abseta < 2.8); ///< @todo Pile-up subtraction
102
103      // Jet/electron/muons overlap removal and selection
104      // Remove electrons within dR = 0.2 of a b-tagged jet
105      for (const Jet& j : jets)
106        if (j.abseta() < 2.5 && j.pT() > 50*GeV && j.bTagged(Cuts::pT > 5*GeV))
107          ifilter_discard(elecs, deltaRLess(j, 0.2, RAPIDITY));
108      // Remove any |eta| < 2.8 jet within dR = 0.2 of a remaining electron
109      for (const Particle& e : elecs)
110        ifilter_discard(jets, deltaRLess(e, 0.2, RAPIDITY));
111      // Remove any electron with dR in [0.2, 0.4] of a remaining jet
112      for (const Jet& j : jets)
113        ifilter_discard(elecs, [&](const Particle& e) { return inRange(deltaR(e,j, RAPIDITY), 0.2, 0.4); });
114      // Remove any muon with dR close to a remaining jet, via a functional form
115      for (const Jet& j : jets)
116        ifilter_discard(muons, [&](const Particle& m) { return deltaR(m,j, RAPIDITY) < min(0.4, 0.04 + 10*GeV/m.pT()); });
117      // Remove any |eta| < 2.8 jet within dR = 0.2 of a remaining muon if track conditions are met
118      for (const Particle& m : muons)
119        /// @todo Add track efficiency random filtering
120        ifilter_discard(jets, [&](const Jet& j) {
121            if (deltaR(j,m, RAPIDITY) > 0.2) return false;
122            const Particles trks = j.particles(Cuts::abscharge > 0 && Cuts::pT > 0.5*GeV);
123            return trks.size() < 3 || (m.pT() > 2*j.pT() && m.pT() > 0.7*sum(trks, pT, 0.0));
124          });
125      // Loose electron selection
126      ifilter_select(elecs, ParticleEffFilter(ELECTRON_EFF_ATLAS_RUN2_LOOSE));
127
128      // Veto the event if there are any remaining baseline leptons
129      if (!elecs.empty()) vetoEvent;
130      if (!muons.empty()) vetoEvent;
131
132      // Passed presel & MET
133      _flows.fill(1);
134
135      // Get jets and their pTs
136      const Jets jets20 = jets;
137      const Jets jets50 = filterBy(jets, Cuts::pT > 50*GeV);
138      const size_t njets50 = jets50.size(), njets20 = jets20.size();
139      if (jets50.size() < 2) vetoEvent;
140      vector<double> jetpts20, jetpts50;
141      transform(jets20, jetpts20, pT);
142      transform(jets50, jetpts50, pT);
143
144      // Construct multi-jet observables
145      const double ht = sum(jetpts20, 0.0);
146      const double met_sqrtHT = met / sqrt(ht);
147      const double meff_incl = sum(jetpts50, met);
148      const double meff_4 = (njets50 >= 4) ? sum(head(jetpts50, 4), met) : -1;
149      const double meff_5 = (njets50 >= 5) ? sum(head(jetpts50, 5), met) : -1;
150      const double meff_6 = (njets50 >= 6) ? sum(head(jetpts50, 6), met) : -1;
151      const double met_meff_4 = met / meff_4;
152      const double met_meff_5 = met / meff_5;
153      const double met_meff_6 = met / meff_6;
154
155      // Jet |eta|s
156      vector<double> jetetas20; transform(jets20, jetetas20, abseta);
157      const double etamax_2 = (njets20 >= 2) ? max(head(jetetas20, 2)) : -1;
158      const double etamax_4 = (njets20 >= 4) ? max(head(jetetas20, 4)) : -1;
159      const double etamax_6 = (njets20 >= 6) ? max(head(jetetas20, 6)) : -1;
160
161      // Get dphis between MET and jets
162      vector<double> dphimets50; transform(jets50, dphimets50, deltaPhiWRT(vmet));
163      const vector<double> dphimets50_123 = head(dphimets50, 3);
164      const vector<double> dphimets50_more = tail(dphimets50, -3);
165      const double dphimin_123 = !dphimets50_123.empty() ? min(dphimets50_123) : -1;
166      const double dphimin_more = !dphimets50_more.empty() ? min(dphimets50_more) : -1;
167
168      // Jet aplanarity
169      Sphericity sph; sph.calc(jets50);
170      const double aplanarity = sph.aplanarity();
171
172
173      //////////////////
174
175
176      // 2 jet regions
177      if (dphimin_123 > 0.8 && dphimin_more > 0.4) {
178        if (jetpts50[1] > 200*GeV && etamax_2 < 0.8) { //< implicit pT[0] cut
179          if (met_sqrtHT > 14*sqrt(GeV) && meff_incl > 800*GeV) _h_2j_0800->fill();
180        }
181        if (jetpts50[1] > 250*GeV && etamax_2 < 1.2) { //< implicit pT[0] cut
182          if (met_sqrtHT > 16*sqrt(GeV) && meff_incl > 1200*GeV) _h_2j_1200->fill();
183          if (met_sqrtHT > 18*sqrt(GeV) && meff_incl > 1600*GeV) _h_2j_1600->fill();
184          if (met_sqrtHT > 20*sqrt(GeV) && meff_incl > 2000*GeV) _h_2j_2000->fill();
185        }
186      }
187
188      // 3 jet region
189      if (njets50 >= 3 && dphimin_123 > 0.4 && dphimin_more > 0.2) {
190        if (jetpts50[0] > 600*GeV && jetpts50[2] > 50*GeV) { //< implicit pT[1] cut
191          if (met_sqrtHT > 16*sqrt(GeV) && meff_incl > 1200*GeV) _h_3j_1200->fill();
192        }
193      }
194
195      // 4 jet regions (note implicit pT[1,2] cuts)
196      if (njets50 >= 4 && dphimin_123 > 0.4 && dphimin_more > 0.4 && jetpts50[0] > 200*GeV && aplanarity > 0.04) {
197        if (jetpts50[3] > 100*GeV && etamax_4 < 1.2 && met_meff_4 > 0.25*sqrt(GeV) && meff_incl > 1000*GeV) _h_4j_1000->fill();
198        if (jetpts50[3] > 100*GeV && etamax_4 < 2.0 && met_meff_4 > 0.25*sqrt(GeV) && meff_incl > 1400*GeV) _h_4j_1400->fill();
199        if (jetpts50[3] > 100*GeV && etamax_4 < 2.0 && met_meff_4 > 0.20*sqrt(GeV) && meff_incl > 1800*GeV) _h_4j_1800->fill();
200        if (jetpts50[3] > 150*GeV && etamax_4 < 2.0 && met_meff_4 > 0.20*sqrt(GeV) && meff_incl > 2200*GeV) _h_4j_2200->fill();
201        if (jetpts50[3] > 150*GeV &&                   met_meff_4 > 0.20*sqrt(GeV) && meff_incl > 2600*GeV) _h_4j_2600->fill();
202      }
203
204      // 5 jet region (note implicit pT[1,2,3] cuts)
205      if (njets50 >= 5 && dphimin_123 > 0.4 && dphimin_more > 0.2 && jetpts50[0] > 500*GeV) {
206        if (jetpts50[4] > 50*GeV && met_meff_5 > 0.3*sqrt(GeV) && meff_incl > 1400*GeV) _h_5j_1400->fill();
207      }
208
209      // 6 jet regions (note implicit pT[1,2,3,4] cuts)
210      if (njets50 >= 6 && dphimin_123 > 0.4 && dphimin_more > 0.2 && jetpts50[0] > 200*GeV && aplanarity > 0.08) {
211        if (jetpts50[5] >  50*GeV && etamax_6 < 2.0 && met_meff_6*sqrt(GeV) > 0.20 && meff_incl > 1800*GeV) _h_6j_1800->fill();
212        if (jetpts50[5] > 100*GeV &&                   met_meff_6*sqrt(GeV) > 0.15 && meff_incl > 2200*GeV) _h_6j_2200->fill();
213      }
214
215      // Cutflows
216      _flows["2j-0800"].filltail({true, dphimin_123 > 0.8, dphimin_more > 0.4, jetpts50[1] > 200*GeV, etamax_2 < 0.8, met_sqrtHT > 14*sqrt(GeV), meff_incl >  800*GeV});
217      _flows["2j-1200"].filltail({true, dphimin_123 > 0.8, dphimin_more > 0.4, jetpts50[1] > 250*GeV, etamax_2 < 1.2, met_sqrtHT > 16*sqrt(GeV), meff_incl > 1200*GeV});
218      _flows["2j-1600"].filltail({true, dphimin_123 > 0.8, dphimin_more > 0.4, jetpts50[1] > 250*GeV, etamax_2 < 1.2, met_sqrtHT > 18*sqrt(GeV), meff_incl > 1600*GeV});
219      _flows["2j-2000"].filltail({true, dphimin_123 > 0.8, dphimin_more > 0.4, jetpts50[1] > 250*GeV, etamax_2 < 1.2, met_sqrtHT > 20*sqrt(GeV), meff_incl > 2000*GeV});
220      _flows["3j-1200"].filltail({njets50 >= 3, dphimin_123 > 0.4, dphimin_more > 0.2, jetpts50[0] > 600*GeV && jetpts50[2] > 50*GeV, true, met_sqrtHT > 16*sqrt(GeV), meff_incl > 1200*GeV});
221      _flows["4j-1000"].filltail({njets50 >= 4, dphimin_123 > 0.4, dphimin_more > 0.4, jetpts50[0] > 200*GeV && jetpts50[3] > 100*GeV, etamax_4 < 1.2, aplanarity > 0.04, met_meff_4 > 0.25*sqrt(GeV), meff_incl > 1000*GeV});
222      _flows["4j-1400"].filltail({njets50 >= 4, dphimin_123 > 0.4, dphimin_more > 0.4, jetpts50[0] > 200*GeV && jetpts50[3] > 100*GeV, etamax_4 < 2.0, aplanarity > 0.04, met_meff_4 > 0.25*sqrt(GeV), meff_incl > 1400*GeV});
223      _flows["4j-1800"].filltail({njets50 >= 4, dphimin_123 > 0.4, dphimin_more > 0.4, jetpts50[0] > 200*GeV && jetpts50[3] > 100*GeV, etamax_4 < 2.0, aplanarity > 0.04, met_meff_4 > 0.20*sqrt(GeV), meff_incl > 1800*GeV});
224      _flows["4j-2200"].filltail({njets50 >= 4, dphimin_123 > 0.4, dphimin_more > 0.4, jetpts50[0] > 200*GeV && jetpts50[3] > 150*GeV, etamax_4 < 2.0, aplanarity > 0.04, met_meff_4 > 0.20*sqrt(GeV), meff_incl > 2200*GeV});
225      _flows["4j-2600"].filltail({njets50 >= 4, dphimin_123 > 0.4, dphimin_more > 0.4, jetpts50[0] > 200*GeV && jetpts50[3] > 150*GeV, true,           aplanarity > 0.04, met_meff_4 > 0.20*sqrt(GeV), meff_incl > 2600*GeV});
226      _flows["5j-1400"].filltail({njets50 >= 5, dphimin_123 > 0.4, dphimin_more > 0.2, jetpts50[0] > 500*GeV && jetpts50[4] > 50*GeV, true, true, met_meff_5 > 0.3*sqrt(GeV), meff_incl > 1400*GeV});
227      _flows["6j-1800"].filltail({njets50 >= 6, dphimin_123 > 0.4, dphimin_more > 0.2, jetpts50[0] > 200*GeV && jetpts50[5] >  50*GeV, etamax_6 < 2.0, aplanarity > 0.08, met_meff_6 > 0.20*sqrt(GeV), meff_incl > 1800*GeV});
228      _flows["6j-2200"].filltail({njets50 >= 6, dphimin_123 > 0.4, dphimin_more > 0.2, jetpts50[0] > 200*GeV && jetpts50[5] > 100*GeV, true,           aplanarity > 0.08, met_meff_6 > 0.15*sqrt(GeV), meff_incl > 2200*GeV});
229
230    }
231
232
233    /// Normalise counters after the run
234    void finalize() {
235
236      const double sf = 13.3*crossSection()/femtobarn/sumOfWeights();
237      scale(_h_2j_0800, sf); scale(_h_2j_1200, sf); scale(_h_2j_1600, sf);
238      scale(_h_2j_2000, sf); scale(_h_3j_1200, sf); scale(_h_4j_1000, sf);
239      scale(_h_4j_1400, sf); scale(_h_4j_1800, sf); scale(_h_4j_2200, sf);
240      scale(_h_4j_2600, sf); scale(_h_5j_1400, sf); scale(_h_6j_1800, sf);
241      scale(_h_6j_2200, sf);
242
243      _flows.scale(sf);
244      MSG_INFO("CUTFLOWS:\n\n" << _flows);
245
246    }
247
248    //@}
249
250
251  private:
252
253    /// @name Histograms
254    //@{
255    CounterPtr _h_2j_0800, _h_2j_1200, _h_2j_1600, _h_2j_2000, _h_3j_1200;
256    CounterPtr _h_4j_1000, _h_4j_1400, _h_4j_1800, _h_4j_2200, _h_4j_2600;
257    CounterPtr _h_5j_1400, _h_6j_1800, _h_6j_2200;
258    //@}
259
260    /// Cut-flows
261    Cutflows _flows;
262
263  };
264
265  // The hook for the plugin system
266  RIVET_DECLARE_PLUGIN(ATLAS_2016_CONF_2016_078);
267}