rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2018_I1662081

Measurement of the differential cross sections of top quark pair production as a function of kinematic event variables in pp collisions at sqrt(s) = 13 TeV
Experiment: CMS (LHC)
Inspire ID: 1662081
Status: VALIDATED
Authors:
  • Doug Burns
  • Emyr Clement
  • Markus Seidel
References:
  • TOP-16-014
  • arxiv:1803.03991
  • Submitted to JHEP
Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • ttbar events at sqrt(s) = 13 TeV, lepton+jets selection at particle level

Abstract : Measurements of the differential $\mathrm{t}\overline{\mathrm{t}}$ production cross section are presented in the single-lepton decay channel, with respect to a number of global event observables. The measurements are performed with $\text{35.9 fb}^{-1}$ of proton-proton collision data collected by the CMS experiment at the LHC during 2016 at $\sqrt{s} = \text{13 TeV}$. The differential cross sections are measured at the level of long-lived generated particles in a phase space similar to that accessible by the CMS detector, and are compared to state-of-the-art leading-order and next-to-leading-order $\mathrm{t}\overline{\mathrm{t}}$ simulations. Rivet : This analysis is to be run on ${\rm t\bar{t}}$ simulation. The particle level selection requires: - exactly one electron or muon with $p_{T} > 26$ GeV and $|\eta| < 2.4$ - no additional electrons or muons with $p_{T} > 15$ GeV and $|\eta| < 2.4$ - at least thress jets with $p_{T} > 30$ GeV and $|\eta| < 2.4$, and one jet with $p_{T} > 25$ GeV and $|\eta| < 2.4$ - two of these selected jets must be tagged as originating from a b quark.

Source code: CMS_2018_I1662081.cc
  1#include "Rivet/Analysis.hh"
  2#include "Rivet/Projections/FinalState.hh"
  3#include "Rivet/Projections/FastJets.hh"
  4#include "Rivet/Projections/ChargedLeptons.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/MissingMomentum.hh"
 10
 11namespace Rivet {
 12
 13
 14  class CMS_2018_I1662081 : public Analysis {
 15  public:
 16
 17    // Minimal constructor
 18    CMS_2018_I1662081()
 19      : Analysis("CMS_2018_I1662081") {}
 20
 21
 22    // Set up projections and book histograms
 23    void init() {
 24      // Complete final state
 25      FinalState fs( (Cuts::abseta < 5) and (Cuts::pT > 0.0*MeV) );
 26
 27      // Dressed leptons
 28      ChargedLeptons charged_leptons(fs);
 29      IdentifiedFinalState photons(fs);
 30      photons.acceptIdPair(PID::PHOTON);
 31
 32      PromptFinalState prompt_leptons(charged_leptons);
 33      prompt_leptons.acceptMuonDecays(true);
 34      prompt_leptons.acceptTauDecays(true);
 35      PromptFinalState prompt_photons(photons);
 36      prompt_photons.acceptMuonDecays(true);
 37      prompt_photons.acceptTauDecays(true);
 38      Cut looseLeptonCuts = Cuts::pt > 15*GeV && Cuts::abseta < 2.4;
 39
 40      DressedLeptons dressed_leptons(prompt_photons, prompt_leptons, 0.1, looseLeptonCuts, true);
 41      declare(dressed_leptons, "DressedLeptons");
 42
 43      // Projection for jets
 44      VetoedFinalState fsForJets(fs);
 45      fsForJets.addVetoOnThisFinalState(dressed_leptons);
 46      declare(FastJets(fsForJets, FastJets::ANTIKT, 0.4), "Jets");
 47
 48      // Projections for MET
 49      declare(MissingMomentum(fs), "MET");
 50
 51      // Booking of histograms
 52      book(_hist_norm_met , 4, 1, 1);
 53      book(_hist_norm_ht  , 2, 1, 1);
 54      book(_hist_norm_st  , 3, 1, 1);
 55      book(_hist_norm_wpt , 5, 1, 1);
 56      book(_hist_norm_njets , 1, 1, 1);
 57      book(_hist_norm_lpt , 6, 1, 1);
 58      book(_hist_norm_labseta , 7, 1, 1);
 59
 60      book(_hist_abs_met , 11, 1, 1);
 61      book(_hist_abs_ht  , 9, 1, 1);
 62      book(_hist_abs_st  , 10, 1, 1);
 63      book(_hist_abs_wpt , 12, 1, 1);
 64      book(_hist_abs_njets , 8, 1, 1);
 65      book(_hist_abs_lpt , 13, 1, 1);
 66      book(_hist_abs_labseta , 14, 1, 1);
 67
 68    }
 69
 70
 71    // per event analysis
 72    void analyze(const Event& event) {
 73
 74      // Lepton veto selection
 75      const DressedLeptons& dressed_leptons = applyProjection<DressedLeptons>(event, "DressedLeptons");
 76      if (dressed_leptons.dressedLeptons().size() != 1) vetoEvent;
 77
 78      // Signal lepton selection
 79      FourMomentum lepton = dressed_leptons.dressedLeptons()[0];
 80
 81      const double leptonPt = lepton.pT();
 82      const double leptonAbsEta = std::abs( lepton.eta() );
 83      if (leptonPt <= 26*GeV || leptonAbsEta >= 2.4) vetoEvent;
 84
 85      // Jet selection
 86      const FastJets& jetpro = applyProjection<FastJets>(event, "Jets");
 87      const Jets jets = jetpro.jets(Cuts::abseta < 2.4 && Cuts::pT > 20*GeV);
 88      Jets cleanedJets;
 89      unsigned int nJetsAbove30GeV = 0;
 90      unsigned int nJetsAbove20GeV = 0;
 91      unsigned int nBJetsAbove30GeV = 0;
 92      unsigned int nBJetsAbove20GeV = 0;
 93      for (const Jet& j : jets) {
 94        cleanedJets.push_back( j );
 95        ++nJetsAbove20GeV;
 96        if ( j.pT() > 30*GeV) ++nJetsAbove30GeV;
 97
 98        if ( j.bTagged() ) {
 99          ++nBJetsAbove20GeV;
100          if ( j.pT() > 30*GeV) ++nBJetsAbove30GeV;
101        }
102      }
103
104      if ( nJetsAbove30GeV < 3 || nJetsAbove20GeV < 4 ) vetoEvent;
105      if ( nBJetsAbove30GeV < 1 || nBJetsAbove20GeV < 2 ) vetoEvent;
106
107      // MET
108      const MissingMomentum& met = applyProjection<MissingMomentum>(event, "MET");
109      _hist_norm_met->fill(met.visibleMomentum().pT()/GeV);
110      _hist_abs_met->fill(met.visibleMomentum().pT()/GeV);
111
112      // HT and ST
113      double ht = 0.0;
114      for (const Jet& j : cleanedJets) ht += j.pT();
115
116      double st = ht + lepton.pT() + met.visibleMomentum().pT();
117      _hist_norm_ht->fill(ht/GeV);
118      _hist_norm_st->fill(st/GeV);
119      _hist_abs_ht->fill(ht/GeV);
120      _hist_abs_st->fill(st/GeV);
121
122      // WPT
123      FourMomentum w = lepton - met.visibleMomentum();
124      _hist_norm_wpt->fill(w.pT()/GeV);
125      _hist_abs_wpt->fill(w.pT()/GeV);
126
127      // Lepton pt and eta
128      _hist_norm_lpt->fill( leptonPt/GeV);
129      _hist_norm_labseta->fill( leptonAbsEta/GeV);
130
131      _hist_abs_lpt->fill( leptonPt/GeV);
132      _hist_abs_labseta->fill( leptonAbsEta/GeV);
133
134      // NJets
135      _hist_norm_njets->fill( cleanedJets.size());
136      _hist_abs_njets->fill( cleanedJets.size());
137
138    }
139
140
141    void finalize() {
142      normalize(_hist_norm_met);
143      normalize(_hist_norm_ht);
144      normalize(_hist_norm_st);
145      normalize(_hist_norm_wpt);
146      normalize(_hist_norm_njets);
147      normalize(_hist_norm_lpt);
148      normalize(_hist_norm_labseta);
149
150      scale(_hist_abs_met, crossSection()/picobarn / sumOfWeights());
151      scale(_hist_abs_ht, crossSection()/picobarn / sumOfWeights());
152      scale(_hist_abs_st, crossSection()/picobarn / sumOfWeights());
153      scale(_hist_abs_wpt, crossSection()/picobarn / sumOfWeights());
154      scale(_hist_abs_njets, crossSection()/picobarn / sumOfWeights());
155      scale(_hist_abs_lpt, crossSection()/picobarn / sumOfWeights());
156      scale(_hist_abs_labseta, crossSection()/picobarn / sumOfWeights());
157
158    }
159
160
161  private:
162    Histo1DPtr _hist_norm_met, _hist_norm_ht, _hist_norm_st, _hist_norm_wpt, _hist_norm_njets, _hist_norm_lpt, _hist_norm_labseta;
163    Histo1DPtr _hist_abs_met, _hist_abs_ht, _hist_abs_st, _hist_abs_wpt, _hist_abs_njets, _hist_abs_lpt, _hist_abs_labseta;
164
165  };
166
167
168  RIVET_DECLARE_PLUGIN(CMS_2018_I1662081);
169
170}