rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2017_I1495243

ttbar + jets at 13 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1495243
Status: VALIDATED
Authors:
  • Callie Bertsche
  • Judith Katzy
  • Krishna Kulkarni
  • Christian Gutschow
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • p + p -> ttbar (dileptonic, needs high statistics [ 2 million] to populate gap fractions).

Measurements of jet activity in top-quark pair events produced in proton--proton collisions are presented, using 3.2 fb$^{-1}$ of $pp$ collision data at a centre-of-mass energy of 13 TeV collected by the ATLAS experiment at the Large Hadron Collider. Events are chosen by requiring an opposite-charge $e\mu$ pair and two $b$-tagged jets in the final state. The normalised differential cross-sections of top-quark pair production are presented as functions of additional-jet multiplicity and transverse momentum, $p_\text{T}$. The fraction of signal events that do not contain additional jet activity in a given rapidity region, the gap fraction, is measured as a function of the $p_\text{T}$ threshold for additional jets, and is also presented for different invariant mass regions of the $e\mu b\bar{b}$ system. All measurements are corrected for detector effects and presented as particle-level distributions compared to predictions with different theoretical approaches for QCD radiation. While the kinematics of the jets from top-quark decays are described well, the generators show differing levels of agreement with the measurements of observables that depend on the production of additional jets.

Source code: ATLAS_2017_I1495243.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  /// @brief $t\bar{t}$ + jets at 13 TeV
 13  class ATLAS_2017_I1495243 : public Analysis {
 14  public:
 15
 16
 17    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2017_I1495243);
 18
 19
 20    void init() {
 21
 22      Cut eta_full = Cuts::abseta < 5.0 && Cuts::pT > 1.0*MeV;
 23      Cut eta_lep = Cuts::abseta < 2.5;
 24
 25      // Collect final state particles
 26      FinalState FS(eta_full);
 27
 28      // Get photons to dress leptons
 29      IdentifiedFinalState photons(FS);
 30      photons.acceptIdPair(PID::PHOTON);
 31
 32      // Projection to find the electrons
 33      IdentifiedFinalState el_id(FS);
 34      el_id.acceptIdPair(PID::ELECTRON);
 35      PromptFinalState electrons(el_id);
 36      electrons.acceptTauDecays(false);
 37      DressedLeptons dressedelectrons(photons, electrons, 0.1, Cuts::abseta< 2.5 && Cuts::pT > 25.0*GeV, true);
 38      declare(dressedelectrons, "electrons");
 39      DressedLeptons fulldressedelectrons(photons, electrons, 0.1, eta_full, true);
 40
 41      // Projection to find the muons
 42      IdentifiedFinalState mu_id(FS);
 43      mu_id.acceptIdPair(PID::MUON);
 44      PromptFinalState muons(mu_id);
 45      muons.acceptTauDecays(false);
 46      DressedLeptons dressedmuons(photons, muons, 0.1, Cuts::abseta < 2.5 && Cuts::pT > 25.0*GeV, true);
 47      declare(dressedmuons, "muons");
 48      DressedLeptons fulldressedmuons(photons, muons, 0.1, eta_full, true);
 49
 50      // Projection to find neutrinos to exclude from jets
 51      IdentifiedFinalState nu_id;
 52      nu_id.acceptNeutrinos();
 53      PromptFinalState neutrinos(nu_id);
 54      neutrinos.acceptTauDecays(false);
 55
 56      // Jet clustering
 57      VetoedFinalState vfs;
 58      vfs.addVetoOnThisFinalState(fulldressedelectrons);
 59      vfs.addVetoOnThisFinalState(fulldressedmuons);
 60      vfs.addVetoOnThisFinalState(neutrinos);
 61      FastJets jets(vfs, FastJets::ANTIKT, 0.4, JetAlg::Muons::ALL, JetAlg::Invisibles::DECAY);
 62      declare(jets, "jets");
 63
 64      // Book Histograms
 65      book(_h["bjet_pt"] , 5,1,1);
 66      book(_h["2bjet_pt"], 6,1,1);
 67      book(_h["ljet_pt"] , 7,1,1);
 68
 69      for (size_t i = 0; i < 4; ++i) {
 70        book(_h["njet"  + to_str(i)], i+1, 1, 1);
 71        book(_h["Q0"    + to_str(i)], "_Q0"    + to_str(i+ 7), refData((i>1?"d":"d0") + to_str(i+ 8) + "-x01-y01"));
 72        book(_h["MQ0"   + to_str(i)], "_MQ0"   + to_str(i+12), refData("d" + to_str(i+12) + "-x01-y01"));
 73        book(_h["Qsum"  + to_str(i)], "_Qsum"  + to_str(i+16), refData("d" + to_str(i+16) + "-x01-y01"));
 74        book(_h["MQsum" + to_str(i)], "_MQsum" + to_str(i+20), refData("d" + to_str(i+20) + "-x01-y01"));
 75        book(_s["gapFracQ0"    + to_str(i)],  8+i, 1 ,1, true);
 76        book(_s["gapFracMQ0"   + to_str(i)], 12+i, 1, 1, true);
 77        book(_s["gapFracQsum"  + to_str(i)], 16+i, 1, 1, true);
 78        book(_s["gapFracMQsum" + to_str(i)], 20+i, 1, 1, true);
 79      }
 80    }
 81
 82
 83    void analyze(const Event& event) {
 84
 85      // Get the selected objects, using the projections.
 86      Jets all_jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 25*GeV && Cuts::abseta < 2.5);
 87
 88      const vector<DressedLepton> electrons = filter_discard(apply<DressedLeptons>(event, "electrons").dressedLeptons(),
 89        [&](const DressedLepton &e) {
 90          return any(all_jets, deltaRLess(e, 0.4));
 91        });
 92
 93      const vector<DressedLepton> muons = filter_discard(apply<DressedLeptons>(event, "muons").dressedLeptons(),
 94        [&](const DressedLepton &m) {
 95          return any(all_jets, deltaRLess(m, 0.4));
 96        });
 97
 98      if (electrons.size() != 1 || muons.size() != 1)  vetoEvent;
 99      if (electrons[0].charge() == muons[0].charge())  vetoEvent;
100
101      Jets bjets, extrajets;
102      for (Jet j : all_jets) {
103        size_t b_tagged = j.bTags(Cuts::pT > 5*GeV).size();
104        if (bjets.size() < 2 && b_tagged)  bjets += j;
105        else  extrajets += j;
106      }
107      if (bjets.size() < 2)  vetoEvent;
108
109      double bjetpt = bjets[0].pt();
110      if (bjetpt > 250*GeV)  bjetpt = 275*GeV;
111      _h["bjet_pt"]->fill(bjetpt);
112
113      double b2jetpt = bjets[1].pt();
114      if (b2jetpt > 150*GeV)  b2jetpt = 175*GeV;
115      _h["2bjet_pt"]->fill(b2jetpt);
116
117      if (extrajets.size()) {
118        double ljetpt = extrajets[0].pt();
119        if (ljetpt > 250*GeV)  ljetpt = 275*GeV;
120        _h["ljet_pt"]->fill(ljetpt);
121      }
122
123      double Memubb = (electrons[0].momentum() + muons[0].momentum() + bjets[0].momentum() + bjets[1].momentum()).mass();
124      vector<double> leadpt = { 0., 0., 0., 0. }, ptsum = { 0., 0., 0., 0. };
125      vector<size_t> njetcount = { 0, 0, 0, 0 };
126      for (size_t i = 0; i < extrajets.size(); ++i) {
127        double absrap = extrajets[i].absrap(), pt = extrajets[i].pT();
128        if (pt > 25*GeV)  ++njetcount[0];
129        if (pt > 40*GeV)  ++njetcount[1];
130        if (pt > 60*GeV)  ++njetcount[2];
131        if (pt > 80*GeV)  ++njetcount[3];
132
133        if (absrap < 0.8 && pt > leadpt[0])  leadpt[0] = pt;
134        else if (absrap > 0.8 && absrap < 1.5 && pt > leadpt[1])  leadpt[1] = pt;
135        else if (absrap > 1.5 && absrap < 2.1 && pt > leadpt[2])  leadpt[2] = pt;
136        if (absrap < 2.1 && pt > leadpt[3])  leadpt[3] = pt;
137
138        if (absrap < 0.8)  ptsum[0] += pt;
139        else if (absrap > 0.8 && absrap < 1.5)  ptsum[1] += pt;
140        else if (absrap > 1.5 && absrap < 2.1)  ptsum[2] += pt;
141        if (absrap < 2.1)  ptsum[3] += pt;
142      }
143
144
145      for (size_t i = 0; i < 4; ++i) {
146        size_t cutoff = i? 3 : 4;
147        if (njetcount[i] > cutoff)  njetcount[i] = cutoff;
148        _h["njet" + to_str(i)]->fill(njetcount[i]);
149
150        if (leadpt[i] > 305*GeV)  leadpt[i] = 305*GeV;
151        _h["Q0" + to_str(i)]->fill(leadpt[i]);
152
153        if (ptsum[i] > 505*GeV)  ptsum[i] = 505*GeV;
154        _h["Qsum" + to_str(i)]->fill(ptsum[i]);
155      }
156
157
158      for (size_t i = 0; i < 4; ++i) {
159        if (i == 0 && !(Memubb < 300*GeV))  continue;
160        if (i == 1 && !(Memubb > 300*GeV && Memubb < 425*GeV))  continue;
161        if (i == 2 && !(Memubb > 425*GeV && Memubb < 600*GeV))  continue;
162        if (i == 3 && !(Memubb > 600*GeV))  continue;
163        _h["MQ0"   + to_str(i)]->fill(leadpt[3]);
164        _h["MQsum" + to_str(i)]->fill(ptsum[3]);
165      }
166    }
167
168
169    void constructGapFraction(Scatter2DPtr out, Histo1DPtr in) {
170      bool hasWeights = in->effNumEntries() != in->numEntries();
171      double denW  = in->sumW();
172      double denW2 = in->sumW2();
173      size_t nEnd = out->numPoints();
174
175      for (size_t i = 0; i < nEnd; ++i) {
176          double numW = in->sumW(), numW2 = in->sumW2();
177          for (size_t j = i; j < nEnd; ++j) {
178            numW  -= in->bin(j).sumW();
179            numW2 -= in->bin(j).sumW2();
180          }
181          double yval = safediv(numW, denW);
182          double yerr = sqrt(safediv(yval * (1 - yval), denW));
183          if (hasWeights) { // use F. James's approximation for weighted events
184            yerr = sqrt( safediv((1 - 2 * yval) * numW2 + yval * yval * denW2, denW * denW) );
185          }
186          out->point(i).setY(yval, yerr);
187      }
188    }
189
190
191    void finalize() {
192
193      // Build gap fraction plots
194      for (size_t i = 0; i < 4; ++i) {
195        constructGapFraction(_s["gapFracQ0"    + to_str(i)], _h["Q0"    + to_str(i)]);
196        constructGapFraction(_s["gapFracMQ0"   + to_str(i)], _h["MQ0"   + to_str(i)]);
197        constructGapFraction(_s["gapFracQsum"  + to_str(i)], _h["Qsum"  + to_str(i)]);
198        constructGapFraction(_s["gapFracMQsum" + to_str(i)], _h["MQsum" + to_str(i)]);
199      }
200
201      // Normalize to cross-section
202      for (map<string, Histo1DPtr>::iterator hit = _h.begin(); hit != _h.end(); ++hit) {
203        if (hit->first.find("jet") != string::npos)  normalize(hit->second);
204      }
205    }
206
207
208  private:
209
210    /// @name Histogram helper functions
211    map<string, Histo1DPtr> _h;
212    map<string, Scatter2DPtr> _s;
213  };
214
215
216  // Declare the class as a hook for the plugin system
217  RIVET_DECLARE_PLUGIN(ATLAS_2017_I1495243);
218
219
220}