rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2012_CONF_2012_103

High jet-multiplicity + MET squark and gluino search
Experiment: ATLAS (LHC)
Status: OBSOLETE
Authors:
  • Peter Richardson
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • BSM signal events at 8000 GeV.

Search for SUSY using events with 6 or more jets in association with missing transverse momentum produced in proton-proton collisions at a centre-of-mass energy of 8 TeV. The data sample has a total integrated luminosity of 5.8 fb$^{-1}$. Distributions in the $W$ and top control regions are not produced, while in addition to the plots from the paper the count of events in the different signal regions is included. The analysis is identical to the previous 7 TeV paper.

Source code: ATLAS_2012_CONF_2012_103.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/VetoedFinalState.hh"
  7#include "Rivet/Projections/IdentifiedFinalState.hh"
  8#include "Rivet/Projections/FastJets.hh"
  9#include "Rivet/Tools/RivetMT2.hh"
 10
 11namespace Rivet {
 12
 13
 14  class ATLAS_2012_CONF_2012_103 : public Analysis {
 15  public:
 16
 17    /// Constructor
 18    ATLAS_2012_CONF_2012_103()
 19      : Analysis("ATLAS_2012_CONF_2012_103")
 20    {    }
 21
 22
 23    /// @name Analysis methods
 24    //@{
 25
 26    /// Book histograms and initialise projections before the run
 27    void init() {
 28
 29      // projection to find the electrons
 30      IdentifiedFinalState elecs(Cuts::abseta < 2.47 && Cuts::pT > 20*GeV);
 31      elecs.acceptIdPair(PID::ELECTRON);
 32      declare(elecs, "elecs");
 33
 34      // projection to find the muons
 35      IdentifiedFinalState muons(Cuts::abseta < 2.4 && Cuts::pT > 10*GeV);
 36      muons.acceptIdPair(PID::MUON);
 37      declare(muons, "muons");
 38
 39      // for pTmiss
 40      declare(VisibleFinalState(Cuts::abseta < 4.9), "vfs");
 41
 42      VetoedFinalState vfs;
 43      vfs.addVetoPairId(PID::MUON);
 44
 45      /// Jet finder
 46      declare(FastJets(vfs, FastJets::ANTIKT, 0.4), "AntiKtJets04");
 47
 48      /// Book histograms
 49      book(_etmiss_HT_7j55 ,"etmiss_HT_7j55", 8, 0., 16.);
 50      book(_etmiss_HT_8j55 ,"etmiss_HT_8j55", 8, 0., 16.);
 51      book(_etmiss_HT_9j55 ,"etmiss_HT_9j55", 8, 0., 16.);
 52      book(_etmiss_HT_6j80 ,"etmiss_HT_6j80", 8, 0., 16.);
 53      book(_etmiss_HT_7j80 ,"etmiss_HT_7j80", 8, 0., 16.);
 54      book(_etmiss_HT_8j80 ,"etmiss_HT_8j80", 8, 0., 16.);
 55
 56      book(_hist_njet55 ,"hist_njet55", 4, 5.5, 9.5);
 57      book(_hist_njet80 ,"hist_njet80", 4, 4.5, 8.5);
 58
 59      book(_count_7j55 ,"count_7j55", 1, 0., 1.);
 60      book(_count_8j55 ,"count_8j55", 1, 0., 1.);
 61      book(_count_9j55 ,"count_9j55", 1, 0., 1.);
 62      book(_count_6j80 ,"count_6j80", 1, 0., 1.);
 63      book(_count_7j80 ,"count_7j80", 1, 0., 1.);
 64      book(_count_8j80 ,"count_8j80", 1, 0., 1.);
 65
 66    }
 67
 68
 69    /// Perform the per-event analysis
 70    void analyze(const Event& event) {
 71      const double weight = 1.0;
 72
 73      // get the jet candidates
 74      Jets cand_jets;
 75      for (const Jet& jet :
 76               apply<FastJets>(event, "AntiKtJets04").jetsByPt(20.0*GeV) ) {
 77        if ( fabs( jet.eta() ) < 2.8 ) {
 78          cand_jets.push_back(jet);
 79        }
 80      }
 81
 82      // candidate muons
 83      Particles cand_mu =
 84        apply<IdentifiedFinalState>(event, "muons").particlesByPt();
 85
 86      // candidate electrons
 87      Particles cand_e  =
 88        apply<IdentifiedFinalState>(event, "elecs").particlesByPt();
 89
 90      // resolve jet/lepton ambiguity
 91      Jets recon_jets;
 92      for ( const Jet& jet : cand_jets ) {
 93        // candidates after |eta| < 2.8
 94        if ( fabs( jet.eta() ) >= 2.8 ) continue;
 95        bool away_from_e = true;
 96        for ( const Particle & e : cand_e ) {
 97          if ( deltaR(e.momentum(),jet.momentum()) <= 0.2 ) {
 98            away_from_e = false;
 99            break;
100          }
101        }
102        if ( away_from_e ) recon_jets.push_back( jet );
103      }
104
105      // only keep electrons more than R=0.4 from jets
106      Particles recon_e;
107      for ( const Particle & e : cand_e ) {
108        bool away = true;
109        for ( const Jet& jet : recon_jets ) {
110          if ( deltaR(e.momentum(),jet.momentum()) < 0.4 ) {
111            away = false;
112            break;
113          }
114        }
115        if ( away )
116          recon_e.push_back( e );
117      }
118
119      // only keep muons more than R=0.4 from jets
120      Particles recon_mu;
121      for ( const Particle & mu : cand_mu ) {
122        bool away = true;
123        for ( const Jet& jet : recon_jets ) {
124          if ( deltaR(mu.momentum(),jet.momentum()) < 0.4 ) {
125            away = false;
126            break;
127          }
128        }
129        if ( away )
130          recon_mu.push_back( mu );
131      }
132
133      // pTmiss
134      Particles vfs_particles =
135        apply<VisibleFinalState>(event, "vfs").particles();
136      FourMomentum pTmiss;
137      for ( const Particle & p : vfs_particles ) {
138        pTmiss -= p.momentum();
139      }
140      double eTmiss = pTmiss.pT();
141
142      // now only use recon_jets, recon_mu, recon_e
143
144      // reject events with electrons and muons
145      if ( ! ( recon_mu.empty() && recon_e.empty() ) ) {
146        MSG_DEBUG("Charged leptons left after selection");
147        vetoEvent;
148      }
149
150      // calculate H_T
151      double HT=0;
152      for ( const Jet& jet : recon_jets ) {
153        if ( jet.pT() > 40 * GeV )
154          HT += jet.pT() ;
155      }
156
157      // number of jets
158      unsigned int njet55=0, njet80=0;
159      for (unsigned int ix=0;ix<recon_jets.size();++ix) {
160        if(recon_jets[ix].pT()>80.*GeV) ++njet80;
161        if(recon_jets[ix].pT()>55.*GeV) ++njet55;
162      }
163
164      double ratio = eTmiss/sqrt(HT);
165
166      if(ratio>4.) {
167        if(njet55>9) njet55 = 9;
168        if(njet80>8) njet80 = 8;
169        _hist_njet55->fill(njet55,weight);
170        _hist_njet80->fill(njet80,weight);
171        // 7j55
172        if(njet55>=7)
173          _count_7j55->fill( 0.5, weight);
174        // 8j55
175        if(njet55>=8)
176          _count_8j55->fill( 0.5, weight) ;
177        // 8j55
178        if(njet55==9)
179          _count_9j55->fill( 0.5, weight) ;
180        // 6j80
181        if(njet80>=6)
182          _count_6j80->fill( 0.5, weight) ;
183        // 7j80
184        if(njet80>=7)
185          _count_7j80->fill( 0.5, weight) ;
186        // 8j80
187        if(njet80==8)
188          _count_8j80->fill( 0.5, weight) ;
189      }
190
191      if(njet55>=7)
192        _etmiss_HT_7j55->fill( ratio, weight);
193      // 8j55
194      if(njet55>=8)
195        _etmiss_HT_8j55->fill( ratio, weight) ;
196      // 8j55
197      if(njet55>=9)
198        _etmiss_HT_9j55->fill( ratio, weight) ;
199      // 6j80
200      if(njet80>=6)
201        _etmiss_HT_6j80->fill( ratio, weight) ;
202      // 7j80
203      if(njet80>=7)
204        _etmiss_HT_7j80->fill( ratio, weight) ;
205      // 8j80
206      if(njet80>=8)
207        _etmiss_HT_8j80->fill( ratio, weight) ;
208
209    }
210
211    //@}
212
213    void finalize() {
214      double norm = crossSection()/femtobarn*5.8/sumOfWeights();
215
216      scale(_etmiss_HT_7j55,2.*norm);
217      scale(_etmiss_HT_8j55,2.*norm);
218      scale(_etmiss_HT_9j55,2.*norm);
219      scale(_etmiss_HT_6j80,2.*norm);
220      scale(_etmiss_HT_7j80,2.*norm);
221      scale(_etmiss_HT_8j80,2.*norm);
222
223      scale(_hist_njet55,norm);
224      scale(_hist_njet80,norm);
225
226      scale(_count_7j55,norm);
227      scale(_count_8j55,norm);
228      scale(_count_9j55,norm);
229      scale(_count_6j80,norm);
230      scale(_count_7j80,norm);
231      scale(_count_8j80,norm);
232    }
233
234  private:
235
236    /// @name Histograms
237    //@{
238    Histo1DPtr _etmiss_HT_7j55;
239    Histo1DPtr _etmiss_HT_8j55;
240    Histo1DPtr _etmiss_HT_9j55;
241    Histo1DPtr _etmiss_HT_6j80;
242    Histo1DPtr _etmiss_HT_7j80;
243    Histo1DPtr _etmiss_HT_8j80;
244
245    Histo1DPtr _hist_njet55;
246    Histo1DPtr _hist_njet80;
247
248    Histo1DPtr _count_7j55;
249    Histo1DPtr _count_8j55;
250    Histo1DPtr _count_9j55;
251    Histo1DPtr _count_6j80;
252    Histo1DPtr _count_7j80;
253    Histo1DPtr _count_8j80;
254    //@}
255
256  };
257
258  // The hook for the plugin system
259  RIVET_DECLARE_PLUGIN(ATLAS_2012_CONF_2012_103);
260
261}