1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5#include "Rivet/Tools/Cutflow.hh"
  6
  7namespace Rivet {
  8
  9
 10  /// Search for R-parity-violating SUSY in multi-jet final states at 13 TeV
 11  class ATLAS_2018_I1667046 : public Analysis {
 12  public:
 13
 14    /// Constructor
 15    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2018_I1667046);
 16
 17
 18    /// @name Analysis methods
 19    /// @{
 20
 21    /// Book histograms and initialise projections before the run
 22    void init() {
 23
 24      // Projections
 25      const FinalState fs (Cuts::abseta < 4.9);
 26      declare("SmallRJ", FastJets(fs, JetAlg::ANTIKT, 0.4));
 27      declare("LargeRJ", FastJets(fs, JetAlg::ANTIKT, 1.0));
 28
 29      // Book histograms
 30      book(_h_sigmaM, "sigmaM",   50, 200, 2000);
 31      book(_h_modeta, "ModEta12", 42,   0,  4.2);
 32
 33      // Cutflows
 34      book(_flows, {"CutFlow1", "CutFlow2"},
 35                   { {"NJet >= 4 ", "Delta12 < 1.4", "PJet1 > 400 GeV", "M SumJ > 1.0 ",
 36                      "NbJet > 0", "M SumJ > 1.0  & NbJet > 0"},
 37                     {"NJet >= 4 ", "Delta12 < 1.4", "NJet >= 5 ", "M SumJ > 0.8 ",
 38                      "NbJet > 0", "M SumJ > 0.8  & NbJet > 0"} });
 39    }
 40
 41
 42    /// Perform the per-event analysis
 43    void analyze(const Event& event) {
 44      _flows->groupfillinit();
 45
 46      // Trim large-R jets and apply cuts
 47      Jets LRJ_old = apply<FastJets>(event, "LargeRJ").jetsByPt(Cuts::abseta < 4.9);
 48      Jets LRJJ;
 49      fastjet::Filter trimmer(fastjet::JetDefinition(fastjet::kt_algorithm,0.2), fastjet::SelectorPtFractionMin(0.05));
 50      for (Jet& j: LRJ_old) LRJJ.push_back(trimmer(j));
 51      Jets LRJ = iselect(LRJJ, Cuts::abseta < 2.0 && Cuts::pT > 200*GeV);
 52      if (LRJ.size() < 4) vetoEvent;
 53      LRJ = sortByPt(LRJ); // sorting for constructing sigmaM
 54
 55      // Small R jets need to pass some cuts + need to be BTagged.
 56      const Jets SRJ = apply<FastJets>(event, "SmallRJ").jetsByPt(Cuts::abseta < 2.5 && Cuts::pT > 50*GeV);
 57      Jets BT_SRJ = select(SRJ, hasBTag(Cuts::pT > 5*GeV));
 58      const int tagg = (BT_SRJ.size() == 0) ? 1 : 0; //if there are no B-TAGGED jets are present, tagg = 1.
 59
 60      // Now to find B-MATCHED Large R Jets: extra step, not useful for SR regions!
 61      const Jets BM_LRJ = selectIfAnyDeltaRLess(LRJ, BT_SRJ, 1.0);
 62
 63
 64      // Now to build observables SigmaM and deltaEta
 65      // Add mass of leading four large R jets,
 66      double sigmaM = 0.0;
 67      for (const Jet& j : head(LRJ, 4)) sigmaM += j.mass();
 68      _h_sigmaM->fill(sigmaM);
 69
 70      // Build deta between two leading large-R jets
 71      const double delta_eta = fabs(deltaEta(LRJ[0], LRJ[1]));
 72      _h_modeta->fill(delta_eta);
 73
 74      // CutFlow1
 75      if (LRJ.size() >= 4) {
 76        _flows->fillnext("CutFlow1");
 77        if (delta_eta < 1.4) {
 78          _flows->fillnext("CutFlow1");
 79          if (LRJ[0].pT() > 400*GeV) {
 80            _flows->fillnext("CutFlow1");
 81            if (sigmaM > 1000*GeV) {
 82              _flows->fillnext("CutFlow1");
 83            }
 84            if (tagg == 0) {
 85              _flows->fillnext("CutFlow1");
 86              if (sigmaM > 1000*GeV){
 87                _flows->fillnext("CutFlow1");
 88              }
 89            } //end of btagg loop
 90          } //end of pT loop
 91        } //end of delta loop
 92      } //end of Njet>4 loop
 93
 94      // CutFlow2
 95      if (LRJ.size() >= 4) {
 96        _flows->fillnext("CutFlow2");
 97        if (delta_eta < 1.4) {
 98          _flows->fillnext("CutFlow2");
 99          if (LRJ.size() >= 5) {
100            _flows->fillnext("CutFlow2");
101            if (sigmaM > 800*GeV) {
102              _flows->fillnext("CutFlow2");
103            }
104            if (tagg == 0) {
105              _flows->fillnext("CutFlow2");
106              if (sigmaM > 800*GeV) {
107                _flows->fillnext("CutFlow2");
108              }
109            } //end of btagg loop
110          } //Njet>5 loop
111        } //end of delta loop
112      } //Njet>4 loop
113
114    }
115
116
117    /// Normalise histograms etc., after the run
118    void finalize() {
119      const double expected = 36.1*crossSection()/femtobarn;
120      normalize(_h_sigmaM, expected/sumOfWeights());
121      normalize(_h_modeta, expected/sumOfWeights());
122      // _flows.scale(99.7/numEvents());
123      // @todo suppress output in reentrant mode?
124      MSG_INFO(_flows);
125    }
126
127    /// @}
128
129
130  private:
131
132    /// @name Histograms
133    Histo1DPtr _h_sigmaM, _h_modeta;
134
135    // Cutflows
136    CutflowsPtr _flows;
137
138  };
139
140
141  RIVET_DECLARE_PLUGIN(ATLAS_2018_I1667046);
142
143}