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