rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2017_I1626105

Dileptonic $t\bar{t}$ at 8 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1626105
Status: VALIDATED
Authors:
  • Judith Katzy
  • Christian Gutschow
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • dileptonic top-quark pair production

This paper presents single lepton and dilepton kinematic distributions measured in dileptonic $t\bar{t}$ events produced in 20.2 fb${}^{-1}$ of $\sqrt{s} = 8$ TeV $pp$ collisions recorded by the ATLAS experiment at the LHC. Both absolute and normalised differential cross-sections are measured, using events with an opposite-charge $e\mu$ pair and one or two $b$-tagged jets. The cross-sections are measured in a fiducial region corresponding to the detector acceptance for leptons, and are compared to the predictions from a variety of Monte Carlo event generators, as well as fixed-order QCD calculations, exploring the sensitivity of the cross-sections to the gluon parton distribution function. Some of the distributions are also sensitive to the top quark pole mass, a combined fit of NLO fixed-order predictions to all the measured distributions yields a top quark mass value of $m_t^\text{pole} = 173.2 \pm 0.9 \pm 0.8 \pm 1.2$ GeV, where the three uncertainties arise from data statistics, experimental systematics, and theoretical sources.

Source code: ATLAS_2017_I1626105.cc
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/VetoedFinalState.hh"
#include "Rivet/Projections/IdentifiedFinalState.hh"
#include "Rivet/Projections/PromptFinalState.hh"
#include "Rivet/Projections/DressedLeptons.hh"
#include "Rivet/Projections/FastJets.hh"

namespace Rivet {


  /// Lepton differential ttbar analysis at 8 TeV
  class ATLAS_2017_I1626105 : public Analysis {
  public:


    DEFAULT_RIVET_ANALYSIS_CTOR(ATLAS_2017_I1626105);


    void init() {

      Cut eta_full = Cuts::abseta < 5.0 && Cuts::pT > 1.0*MeV;

      // All final state particles
      const FinalState fs;

      // Get photons to dress leptons
      IdentifiedFinalState photons(fs);
      photons.acceptIdPair(PID::PHOTON);

      // Projection to find the electrons
      PromptFinalState prompt_el(Cuts::abspid == PID::ELECTRON, true);
      DressedLeptons elecs(photons, prompt_el, 0.1, (Cuts::abseta < 2.5) && (Cuts::pT > 25*GeV));
      DressedLeptons veto_elecs(photons, prompt_el, 0.1, eta_full, false);
      declare(elecs, "elecs");

      // Projection to find the muons
      PromptFinalState prompt_mu(Cuts::abspid == PID::MUON, true);
      DressedLeptons muons(photons, prompt_mu, 0.1, (Cuts::abseta < 2.5) && (Cuts::pT > 25*GeV));
      DressedLeptons veto_muons(photons, prompt_mu, 0.1, eta_full, false);
      declare(muons, "muons");

      // Jet clustering.
      VetoedFinalState vfs;
      vfs.addVetoOnThisFinalState(veto_elecs);
      vfs.addVetoOnThisFinalState(veto_muons);
      FastJets jets(vfs, FastJets::ANTIKT, 0.4);
      jets.useInvisibles();
      declare(jets, "jets");

      // Book histograms
      bookHistos("lep_pt",       1);
      bookHistos("lep_eta",      3);
      bookHistos("dilep_pt",     5);
      bookHistos("dilep_mass",   7);
      bookHistos("dilep_rap",    9);
      bookHistos("dilep_dphi",  11);
      bookHistos("dilep_sumpt", 13);
      bookHistos("dilep_sumE",  15);
    }


    void analyze(const Event& event) {
      vector<DressedLepton> elecs = sortByPt(apply<DressedLeptons>(event, "elecs").dressedLeptons());
      vector<DressedLepton> muons = sortByPt(apply<DressedLeptons>(event, "muons").dressedLeptons());
      Jets jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 25*GeV && Cuts::abseta < 2.5);

      // Check overlap of jets/leptons.
      for (const Jet& jet : jets) {
        ifilter_discard(elecs, deltaRLess(jet, 0.4));
        ifilter_discard(muons, deltaRLess(jet, 0.4));
      }
      if (elecs.empty() || muons.empty()) vetoEvent;
      if (elecs[0].charge() == muons[0].charge()) vetoEvent;

      FourMomentum el = elecs[0].momentum();
      FourMomentum mu = muons[0].momentum();
      FourMomentum ll = elecs[0].momentum() + muons[0].momentum();

      // Fill histograms
      const double weight = event.weight();
      fillHistos("lep_pt",      el.pT()/GeV,             weight);
      fillHistos("lep_pt",      mu.pT()/GeV,             weight);
      fillHistos("lep_eta",     el.abseta(),             weight);
      fillHistos("lep_eta",     mu.abseta(),             weight);
      fillHistos("dilep_pt",    ll.pT()/GeV,             weight);
      fillHistos("dilep_mass",  ll.mass()/GeV,           weight);
      fillHistos("dilep_rap",   ll.absrap(),             weight);
      fillHistos("dilep_dphi",  deltaPhi(el, mu),        weight);
      fillHistos("dilep_sumpt", (el.pT() + mu.pT())/GeV, weight);
      fillHistos("dilep_sumE",  (el.E() + mu.E())/GeV,   weight);
    }

    void finalize() {
      // Normalize to cross-section
      const double sf = crossSection()/femtobarn/sumOfWeights();
      for (auto& hist : _h) {
        const double norm = 1.0 / hist.second->integral();
        // add overflow to last bin
        double overflow = hist.second->overflow().effNumEntries();
        hist.second->fillBin(hist.second->numBins() - 1, overflow);
        // histogram normalisation
        if (hist.first.find("norm") != string::npos)  scale(hist.second, norm);
        else  scale(hist.second, sf);
      }

    }

  private:

    /// @name Histogram helper functions
    //@{
    void bookHistos(const std::string name, unsigned int index) {
      _h[name] = bookHisto1D(index, 1, 1);
      _h["norm_" + name] = bookHisto1D(index + 1, 1, 1);
    }

    void fillHistos(const std::string name, double value, double weight) {
      _h[name]->fill(value, weight);
      _h["norm_" + name]->fill(value, weight);
    }

    map<string, Histo1DPtr> _h;
    //@}

  };

  // Declare the class as a hook for the plugin system
  DECLARE_RIVET_PLUGIN(ATLAS_2017_I1626105);
}