rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2015_I1397635

$Wt$ at 8 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1397635
Status: VALIDATED
Authors:
  • Reinhard Schwienhorst
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • Single top $Wt$ + top--antitop production

Fiducial cross-section for $Wt+t\bar{t}$ production in events with two leptons and exactly one jet, using an integrated luminosity of 20.3 fb${}^{-1}$ of proton-proton collisions at a centre-of-mass energy of 8 TeV at the Large Hadron Collider, collected with the ATLAS detector. The cross-section for the production of a top quark and a $W$ boson is measured in a fiducial acceptance requiring two leptons with $p_\perp > 25$ GeV and $|\eta| < 2.5$, one jet with $p_\perp > 20$ GeV and $|\eta| < 2.5$, and $E^\text{miss}_\text{T} > 20$ GeV, including both $Wt$ and top-quark pair events as signal. The measured value of the fiducial cross-section is $0.85 \pm 0.01$ (stat.) $\pm 0.07$ (syst.) $\pm 0.03$ (lumi.) pb.

Source code: ATLAS_2015_I1397635.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/VetoedFinalState.hh"
  4#include "Rivet/Projections/IdentifiedFinalState.hh"
  5#include "Rivet/Projections/PromptFinalState.hh"
  6#include "Rivet/Projections/DressedLeptons.hh"
  7#include "Rivet/Projections/FastJets.hh"
  8
  9namespace Rivet {
 10
 11
 12  /// $Wt$ at 8 TeV
 13  class ATLAS_2015_I1397635 : public Analysis {
 14  public:
 15
 16    /// Constructor
 17    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2015_I1397635);
 18
 19
 20    void init() {
 21
 22      // Eta ranges
 23      Cut eta_full = Cuts::abseta < 5.0 && Cuts::pT >= 1.0*MeV;
 24      Cut eta_lep  = Cuts::abseta < 2.5;
 25
 26      // All final state particles
 27      FinalState fs(eta_full);
 28
 29      // Get photons to dress leptons
 30      IdentifiedFinalState photons(fs);
 31      photons.acceptIdPair(PID::PHOTON);
 32
 33      // Projection to find the electrons
 34      IdentifiedFinalState el_id(fs);
 35      el_id.acceptIdPair(PID::ELECTRON);
 36      PromptFinalState electrons(el_id);
 37      electrons.acceptTauDecays(true);
 38      declare(electrons, "electrons");
 39      DressedLeptons dressedelectrons(photons, electrons, 0.1, eta_lep && Cuts::pT > 25*GeV, true);
 40      declare(dressedelectrons, "dressedelectrons");
 41      DressedLeptons ewdressedelectrons(photons, electrons, 0.1, eta_full, true);
 42
 43      // Projection to find the muons
 44      IdentifiedFinalState mu_id(fs);
 45      mu_id.acceptIdPair(PID::MUON);
 46      PromptFinalState muons(mu_id);
 47      muons.acceptTauDecays(true);
 48      declare(muons, "muons");
 49      DressedLeptons dressedmuons(photons, muons, 0.1, eta_lep && Cuts::pT > 25*GeV, true);
 50      declare(dressedmuons, "dressedmuons");
 51      DressedLeptons ewdressedmuons(photons, muons, 0.1, eta_full, true);
 52
 53      // Projection to find neutrinos and produce MET
 54      IdentifiedFinalState nu_id;
 55      nu_id.acceptNeutrinos();
 56      PromptFinalState neutrinos(nu_id);
 57      neutrinos.acceptTauDecays(true);
 58      declare(neutrinos, "neutrinos");
 59
 60
 61      // Jet clustering.
 62      VetoedFinalState vfs;
 63      vfs.addVetoOnThisFinalState(ewdressedelectrons);
 64      vfs.addVetoOnThisFinalState(ewdressedmuons);
 65      vfs.addVetoOnThisFinalState(neutrinos);
 66      FastJets jets(vfs,FastJets::ANTIKT, 0.4);
 67      jets.useInvisibles();
 68      declare(jets, "jets");
 69
 70      book(_histo ,1,1,1);
 71    }
 72
 73
 74    void analyze(const Event& event) {
 75
 76      // Get the selected objects, using the projections.
 77      vector<DressedLepton> electrons = apply<DressedLeptons>(event, "dressedelectrons").dressedLeptons();
 78      vector<DressedLepton> muons = apply<DressedLeptons>(event, "dressedmuons").dressedLeptons();
 79      // also make basic event selection cuts for leptons
 80      if (electrons.empty() && muons.empty())  vetoEvent;
 81      if (electrons.size() + muons.size() != 2) vetoEvent;
 82
 83      // next selection cuts for jets
 84      const Jets jets = apply<FastJets>(event, "jets").jets(Cuts::pT>20*GeV && Cuts::abseta < 2.5, cmpMomByPt);
 85      if (jets.size() != 1) vetoEvent;
 86
 87      // and selection cuts for b-tagging
 88      Jets bjets;
 89      // Check overlap of jets/leptons.
 90      for (Jet jet : jets) {
 91        // if dR(el,jet) < 0.4 skip the event
 92        for (DressedLepton el : electrons) {
 93          if (deltaR(jet, el) < 0.4)  vetoEvent;
 94        }
 95        // if dR(mu,jet) < 0.4 skip the event
 96        for (DressedLepton mu : muons) {
 97          if (deltaR(jet, mu) < 0.4)  vetoEvent;
 98        }
 99        // Count the number of b-tags
100        // We have to check that the ghost-matched B hadrons have pT > 5 GeV
101        // By default jet.bTags() returns all B hadrons without cuts
102        bool btagged = jet.bTags(Cuts::pT >= 5*GeV).size();
103        if (btagged)  bjets += jet;
104      }
105      if (bjets.size() != 1) vetoEvent;
106
107      // Now evaluate MET selection
108      // Get the neutrinos from the event record (they have pT > 0.0 and |eta| < 4.5 at this stage
109      const Particles& neutrinos = apply<PromptFinalState>(event, "neutrinos").particlesByPt();
110      FourMomentum met;
111      for (const Particle& nu : neutrinos)  met += nu.momentum();
112      if (met.pT() <= 20*GeV)  vetoEvent;
113
114      // Make the plot
115      _histo->fill(1);
116    }
117
118
119    // Normalise histograms etc., after the run
120    void finalize() {
121      scale(_histo, crossSection() / femtobarn / sumOfWeights());
122    }
123
124
125  private:
126
127    Histo1DPtr _histo;
128
129  };
130
131
132  // Hook for the plugin system
133  RIVET_DECLARE_PLUGIN(ATLAS_2015_I1397635);
134
135}