rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2018_I1705857

ttbb at 13 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1705857
Status: VALIDATED
Authors:
  • Stewart Swift
  • Tom Neep
  • Christian Gutschow
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • ttbb production at 13 TeV in dilepton and lepton+jets channel

This paper presents measurements of $t\bar{t}$ production in association with additional $b$-jets in $pp$ collisions at the LHC at a centre-of-mass energy of 13 TeV. The data were recorded with the ATLAS detector and correspond to an integrated luminosity of 36.1 fb$^{-1}$. Fiducial cross-section measurements are performed in the dilepton and lepton-plus-jets $t\bar{t}$ decay channels. Results are presented at particle level in the form of inclusive cross-sections of $t\bar{t}$ final states with three and four $b$-jets as well as differential cross-sections as a function of global event properties and properties of $b$-jet pairs. The measured inclusive fiducial cross-sections generally exceed the $t\bar{t}b\bar{b}$ predictions from various next-to-leading-order matrix element calculations matched to a parton shower but are compatible within the total uncertainties. The experimental uncertainties are smaller than the uncertainties in the predictions. Comparisons of state-of-the-art theoretical predictions with the differential measurements are shown and good agreement with data is found for most of them.

Source code: ATLAS_2018_I1705857.cc
  1#include "Rivet/Analysis.hh"
  2#include "Rivet/Projections/FastJets.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/DressedLeptons.hh"
  5#include "Rivet/Projections/PromptFinalState.hh"
  6#include "Rivet/Projections/VetoedFinalState.hh"
  7
  8namespace Rivet {
  9
 10
 11/// @brief ttbb at 13 TeV
 12class ATLAS_2018_I1705857 : public Analysis {
 13 public:
 14   /// Constructor
 15   RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2018_I1705857);
 16
 17    void init() {
 18      // Eta ranges
 19      Cut eta_full = (Cuts::abseta < 5.0);
 20      // Lepton cuts
 21      Cut lep_cuts25 = (Cuts::abseta < 2.5) && (Cuts::pT >= 25*GeV);
 22      // All final state particles
 23      FinalState fs(eta_full);
 24
 25      // Get photons to dress leptons
 26      PromptFinalState photons(eta_full && Cuts::abspid == PID::PHOTON, true);
 27
 28      // Projection to find the electrons
 29      PromptFinalState electrons(eta_full && Cuts::abspid == PID::ELECTRON, true);
 30
 31      // Projection to find the muons
 32      PromptFinalState muons(eta_full && Cuts::abspid == PID::MUON, true);
 33
 34      DressedLeptons dressedelectrons25(photons, electrons, 0.1, lep_cuts25, true);
 35      DressedLeptons dressedmuons25(photons, muons, 0.1, lep_cuts25, true);
 36
 37      declare(dressedelectrons25, "elecs");
 38      declare(dressedmuons25, "muons");
 39
 40      // From here on we are just setting up the jet clustering
 41      IdentifiedFinalState nu_id;
 42      nu_id.acceptNeutrinos();
 43      PromptFinalState neutrinos(nu_id);
 44      neutrinos.acceptTauDecays(true);
 45
 46      PromptFinalState jet_photons(eta_full && Cuts::abspid == PID::PHOTON, false);
 47      DressedLeptons all_dressed_electrons(jet_photons, electrons, 0.1, eta_full, true);
 48      DressedLeptons all_dressed_muons(jet_photons, muons, 0.1, eta_full, true);
 49
 50      VetoedFinalState vfs(fs);
 51      vfs.addVetoOnThisFinalState(all_dressed_electrons);
 52      vfs.addVetoOnThisFinalState(all_dressed_muons);
 53      vfs.addVetoOnThisFinalState(neutrinos);
 54
 55      FastJets jets(vfs, FastJets::ANTIKT, 0.4, JetAlg::Muons::DECAY, JetAlg::Invisibles::DECAY);
 56      declare(jets, "jets");
 57
 58      // fiducial cross-section histogram
 59      book(_histograms["fid_xsec"], 1, 1, 1);
 60      book(_histograms["fid_xsec_no_ttX"], 2, 1, 1);
 61
 62      book(_histograms["nbjets_emu"], 3, 1, 1);
 63      book(_histograms["nbjets_emu_no_ttX"], 4, 1, 1);
 64
 65      // HT
 66      book_hist("ht_emu", 3);
 67      book_hist("ht_had_emu", 4);
 68
 69      book_hist("ht_ljets", 5);
 70      book_hist("ht_had_ljets", 6);
 71
 72      // b-jet pTs
 73      book_hist("lead_bjet_pt_emu",    7);
 74      book_hist("sublead_bjet_pt_emu", 8);
 75      book_hist("third_bjet_pt_emu",   9);
 76
 77      book_hist("lead_bjet_pt_ljets",    10);
 78      book_hist("sublead_bjet_pt_ljets", 11);
 79      book_hist("third_bjet_pt_ljets",   12);
 80      book_hist("fourth_bjet_pt_ljets",  13);
 81
 82      // leading bb pair
 83      book_hist("m_bb_leading_emu",  14);
 84      book_hist("pt_bb_leading_emu", 15);
 85      book_hist("dR_bb_leading_emu", 16);
 86
 87      book_hist("m_bb_leading_ljets",  17);
 88      book_hist("pt_bb_leading_ljets", 18);
 89      book_hist("dR_bb_leading_ljets", 19);
 90
 91      // closest bb pair
 92      book_hist("m_bb_closest_emu",  20);
 93      book_hist("pt_bb_closest_emu", 21);
 94      book_hist("dR_bb_closest_emu", 22);
 95
 96      book_hist("m_bb_closest_ljets",  23);
 97      book_hist("pt_bb_closest_ljets", 24);
 98      book_hist("dR_bb_closest_ljets", 25);
 99    }
100
101
102    void analyze(const Event& event) {
103
104      vector<DressedLepton> leptons;
105      for (auto &lep : apply<DressedLeptons>(event, "muons").dressedLeptons()) { leptons.push_back(lep); }
106      for (auto &lep : apply<DressedLeptons>(event, "elecs").dressedLeptons()) { leptons.push_back(lep); }
107
108      const Jets jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 25*GeV && Cuts::abseta < 2.5);
109      for (const auto& jet : jets) {
110        ifilter_discard(leptons, [&](const DressedLepton& lep) { return deltaR(jet, lep) < 0.4; });
111      }
112
113      Jets bjets;
114      for (const Jet& jet : jets) {
115        if (jet.bTagged(Cuts::pT >= 5*GeV))  bjets += jet;
116      }
117
118      size_t njets = jets.size();
119      size_t nbjets = bjets.size();
120
121      // Evaluate basic event selection
122      bool pass_ljets = (leptons.size() == 1 && leptons[0].pT() > 27*GeV);
123
124      bool pass_emu =
125        // 2 leptons > 27 GeV
126        (leptons.size() == 2) &&
127        (leptons[0].pT() > 27*GeV && leptons[1].pT() > 27*GeV) &&
128        // emu events
129        ((leptons[0].abspid() == 11 && leptons[1].abspid() == 13) ||
130         (leptons[0].abspid() == 13 && leptons[1].abspid() == 11)) &&
131        // opposite charge
132        (leptons[0].charge() != leptons[1].charge());
133
134      // If we don't have exactly 1 or 2 leptons then veto the event
135      if (!pass_emu && !pass_ljets)  vetoEvent;
136
137      if (pass_emu) {
138        if (nbjets >= 2)  fill("nbjets_emu", nbjets - 1);
139        if (nbjets >= 3)  fill("fid_xsec", 1);
140        if (nbjets >= 4)  fill("fid_xsec", 2);
141      }
142
143      if (pass_ljets) {
144        if (nbjets >= 3 && njets >= 5)  fill("fid_xsec", 3);
145        if (nbjets >= 4 && njets >= 6)  fill("fid_xsec", 4);
146      }
147
148      if (pass_emu && (nbjets < 3 || njets < 3))    vetoEvent;
149      if (pass_ljets && (nbjets < 4 || njets < 6))  vetoEvent;
150
151      double hthad = sum(jets, pT, 0.0);
152      double ht = sum(leptons, pT, hthad);
153
154      FourMomentum jsum = bjets[0].momentum() + bjets[1].momentum();
155      double dr_leading = deltaR(bjets[0], bjets[1]);
156
157      size_t ind1 = 0, ind2 = 0; double mindr = 999.;
158      for (size_t i = 0; i < bjets.size(); ++i) {
159        for (size_t j = 0; j < bjets.size(); ++j) {
160          if (i == j)  continue;
161          double dr = deltaR(bjets[i], bjets[j]);
162          if (dr < mindr) {
163            ind1 = i;
164            ind2 = j;
165            mindr = dr;
166          }
167        }
168      }
169
170      FourMomentum bb_closest = bjets[ind1].momentum() + bjets[ind2].momentum();
171      double dr_closest = deltaR(bjets[ind1], bjets[ind2]);
172
173      if (pass_ljets) {
174        // b-jet pTs
175        fill("lead_bjet_pt_ljets", bjets[0].pT()/GeV);
176        fill("sublead_bjet_pt_ljets", bjets[1].pT()/GeV);
177        fill("third_bjet_pt_ljets", bjets[2].pT()/GeV);
178
179        if (nbjets >= 4)  fill("fourth_bjet_pt_ljets", bjets[3].pT()/GeV);
180
181        // HT
182        fill("ht_ljets", ht/GeV);
183        fill("ht_had_ljets", hthad/GeV);
184
185        // leading bb pair
186        fill("m_bb_leading_ljets", jsum.mass()/GeV);
187        fill("pt_bb_leading_ljets", jsum.pT()/GeV);
188        fill("dR_bb_leading_ljets", dr_leading);
189
190        // closest bb pair
191        fill("m_bb_closest_ljets", bb_closest.mass()/GeV);
192        fill("pt_bb_closest_ljets", bb_closest.pT()/GeV);
193        fill("dR_bb_closest_ljets", dr_closest);
194      }
195      if (pass_emu) {
196        // b-jet pTs
197        fill("lead_bjet_pt_emu", bjets[0].pT()/GeV);
198        fill("sublead_bjet_pt_emu", bjets[1].pT()/GeV);
199        fill("third_bjet_pt_emu", bjets[2].pT()/GeV);
200
201        // HT
202        fill("ht_emu", ht/GeV);
203        fill("ht_had_emu", hthad/GeV);
204
205        // leading bb pair
206        fill("m_bb_leading_emu", jsum.mass()/GeV);
207        fill("pt_bb_leading_emu", jsum.pT()/GeV);
208        fill("dR_bb_leading_emu", dr_leading);
209
210        // closest bb pair
211        fill("m_bb_closest_emu", bb_closest.mass()/GeV);
212        fill("pt_bb_closest_emu", bb_closest.pT()/GeV);
213        fill("dR_bb_closest_emu", dr_closest);
214      }
215    }
216
217    void finalize() {
218      // Normalise all histograms
219      const double sf = crossSection() / femtobarn / sumOfWeights();
220      for (auto const& h : _histograms) {
221        scale(h.second, sf);
222        if (h.first.find("fid_xsec") != string::npos)  continue;
223        normalize(h.second, 1.0);
224      }
225    }
226
227    void fill(const string& name, const double value) {
228      _histograms[name]->fill(value);
229      _histograms[name + "_no_ttX"]->fill(value);
230    }
231
232    void book_hist(const std::string& name, unsigned int d) {
233      // ttX substracted histograms are even numbers
234      book(_histograms[name], (d * 2) - 1, 1, 1);
235      book(_histograms[name + "_no_ttX"], d * 2, 1, 1);
236    }
237
238    private:
239      map<std::string, Histo1DPtr> _histograms;
240
241  };
242
243  // The hook for the plugin system
244  RIVET_DECLARE_PLUGIN(ATLAS_2018_I1705857);
245}