rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2016_I1473674

Measurement of the differential cross sections for top quark pair production as a function of kinematic event variables at sqrt(s) = 8 TeV
Experiment: CMS (LHC)
Inspire ID: 1473674
Status: VALIDATED
Authors:
  • Markus Seidel
  • Lukas Kreczko
  • Emyr Clement
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • ttbar events at sqrt(s) = 8 TeV (inclusive or lepton+jets decay mode) Top quarks are expected in the event record to identify the decay mode.

Measurements are reported of the normalized differential cross sections for top quark pair production with respect to four kinematic event variables: the missing transverse energy; the scalar sum of the jet transverse momentum (pT); the scalar sum of the pT of all objects in the event; and the pT of leptonically decaying W bosons from top quark decays.

Source code: CMS_2016_I1473674.cc
  1#include "Rivet/Analysis.hh"
  2#include "Rivet/Projections/FinalState.hh"
  3#include "Rivet/Projections/FastJets.hh"
  4#include "Rivet/Projections/PartonicTops.hh"
  5#include "Rivet/Projections/DressedLeptons.hh"
  6#include "Rivet/Projections/IdentifiedFinalState.hh"
  7#include "Rivet/Projections/PromptFinalState.hh"
  8#include "Rivet/Projections/VetoedFinalState.hh"
  9#include "Rivet/Projections/InvMassFinalState.hh"
 10#include "Rivet/Projections/MissingMomentum.hh"
 11
 12namespace Rivet {
 13
 14
 15  class CMS_2016_I1473674 : public Analysis {
 16  public:
 17
 18    // Minimal constructor
 19    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2016_I1473674);
 20
 21
 22    // Set up projections and book histograms
 23    void init() {
 24
 25      // Complete final state
 26      FinalState fs;
 27
 28      // Parton level top quarks
 29      declare(PartonicTops(PartonicTops::DecayMode::E_MU, false), "LeptonicPartonTops");
 30      declare(PartonicTops(PartonicTops::DecayMode::HADRONIC),    "HadronicPartonTops");
 31
 32      // Projections for dressed electrons and muons
 33      IdentifiedFinalState photons(fs);
 34      photons.acceptIdPair(PID::PHOTON);
 35      //
 36      IdentifiedFinalState el_id(fs);
 37      el_id.acceptIdPair(PID::ELECTRON);
 38      PromptFinalState electrons(el_id);
 39      declare(electrons, "Electrons");
 40      DressedLeptons dressed_electrons(photons, electrons, 0.1);
 41      declare(dressed_electrons, "DressedElectrons");
 42      //
 43      IdentifiedFinalState mu_id(fs);
 44      mu_id.acceptIdPair(PID::MUON);
 45      PromptFinalState muons(mu_id);
 46      declare(muons, "Muons");
 47      DressedLeptons dressed_muons(photons, muons, 0.1);
 48      declare(dressed_muons, "DressedMuons");
 49
 50      // Projection for jets
 51      VetoedFinalState fs_jets;
 52      fs_jets.addVetoOnThisFinalState(dressed_muons);
 53      declare(FastJets(fs_jets, FastJets::ANTIKT, 0.5), "Jets");
 54
 55      // Projections for MET
 56      declare(MissingMomentum(), "MET");
 57
 58
 59      // Booking of histograms
 60      book(_hist_met ,5, 1, 1);
 61      book(_hist_ht  ,6, 1, 1);
 62      book(_hist_st  ,7, 1, 1);
 63      book(_hist_wpt ,8, 1, 1);
 64    }
 65
 66
 67    /// Per-event analysis
 68    void analyze(const Event& event) {
 69      const double weight = 1.0;
 70
 71      // Select ttbar -> lepton+jets at parton level, removing tau decays
 72      const Particles leptonicpartontops = apply<ParticleFinder>(event, "LeptonicPartonTops").particlesByPt();
 73      if (leptonicpartontops.size() != 1) vetoEvent;
 74      const Particles hadronicpartontops = apply<ParticleFinder>(event, "HadronicPartonTops").particlesByPt();
 75      if (hadronicpartontops.size() != 1) vetoEvent;
 76
 77      // Select ttbar -> lepton+jets at particle level
 78      const DressedLeptons& dressed_electrons = applyProjection<DressedLeptons>(event, "DressedElectrons");
 79      const DressedLeptons& dressed_muons = applyProjection<DressedLeptons>(event, "DressedMuons");
 80      if (dressed_electrons.dressedLeptons().size() + dressed_muons.dressedLeptons().size() != 1) vetoEvent;
 81      const FourMomentum lepton = (dressed_electrons.dressedLeptons().empty() ? dressed_muons : dressed_electrons).dressedLeptons()[0];
 82
 83      // MET
 84      const MissingMomentum& met = applyProjection<MissingMomentum>(event, "MET");
 85      _hist_met->fill(met.visibleMomentum().pT()/GeV, weight);
 86
 87      // HT and ST
 88      const FastJets& jetpro = applyProjection<FastJets>(event, "Jets");
 89      const Jets jets = jetpro.jetsByPt(20*GeV);
 90
 91      double ht = 0.0;
 92      for (const Jet& j : jets) {
 93        if (deltaR(j.momentum(), lepton) > 0.3) {
 94          ht += j.pT();
 95        }
 96      }
 97
 98      double st = ht + lepton.pT() + met.visibleMomentum().pT();
 99      _hist_ht->fill(ht/GeV, weight);
100      _hist_st->fill(st/GeV, weight);
101
102      // WPT
103      const FourMomentum w = lepton - met.visibleMomentum();
104      _hist_wpt->fill(w.pT()/GeV, weight);
105    }
106
107
108    /// Normalize histograms
109    void finalize() {
110      normalize(_hist_met);
111      normalize(_hist_ht);
112      normalize(_hist_st);
113      normalize(_hist_wpt);
114    }
115
116  private:
117    Histo1DPtr _hist_met, _hist_ht, _hist_st, _hist_wpt;
118
119  };
120
121
122  RIVET_DECLARE_PLUGIN(CMS_2016_I1473674);
123
124}