rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2018_I1677498

WWbb at 13 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1677498
Status: VALIDATED
Authors:
  • Christian Herwig
  • Christian Gutschow
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • dileptonic WWbb production at 13 TeV

This Letter presents a normalized differential cross-section measurement in a fiducial phase-space region where interference effects between top-quark pair production and associated production of a single top quark with a $W$ boson and a $b$-quark are significant. Events with exactly two leptons ($ee$, $\mu\mu$, or $e\mu$) and two $b$-tagged jets that satisfy a multiparticle invariant mass requirement are selected from 36.1 fb$^{-1}$ of proton-proton collision data taken at $\sqrt{s}=13$ TeV with the ATLAS detector at the LHC in 2015 and 2016. The results are compared with predictions from simulations using various strategies for the interference. The standard prescriptions for interference modeling are significantly different from each other but are within $2\sigma$ of the data. State-of-the-art predictions that naturally incorporate interference effects provide the best description of the data in the measured region of phase space most sensitive to these effects. These results provide an important constraint on interference models and will guide future model development and tuning.

Source code: ATLAS_2018_I1677498.cc
  1#include "Rivet/Analysis.hh"
  2#include "Rivet/Projections/FinalState.hh"
  3#include "Rivet/Projections/PromptFinalState.hh"
  4#include "Rivet/Projections/LeptonFinder.hh"
  5#include "Rivet/Projections/FastJets.hh"
  6
  7namespace Rivet {
  8
  9
 10  /// WWbb at 13 TeV
 11  class ATLAS_2018_I1677498 : public Analysis {
 12  public:
 13
 14    /// Constructor
 15    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2018_I1677498);
 16
 17
 18    /// Book cuts and projections
 19    void init() {
 20
 21      // All final state particles
 22      FinalState fs(Cuts::abseta < 5.0);
 23
 24      PromptFinalState photons(Cuts::abspid == PID::PHOTON, TauDecaysAs::PROMPT);
 25      PromptFinalState bare_el(Cuts::abspid == PID::ELECTRON, TauDecaysAs::PROMPT);
 26      LeptonFinder elecs(bare_el, photons, 0.1, Cuts::pT > 7*GeV && Cuts::abseta < 2.47);
 27      declare(elecs, "elecs");
 28
 29      PromptFinalState bare_mu(Cuts::abspid == PID::MUON, TauDecaysAs::PROMPT);
 30      LeptonFinder muons(bare_mu, photons, 0.1, Cuts::pT > 6*GeV && Cuts::abseta < 2.5);
 31      declare(muons, "muons");
 32
 33      FastJets jets(fs, JetAlg::ANTIKT, 0.4, JetMuons::NONE, JetInvisibles::NONE);
 34      declare(jets, "jets");
 35
 36      book(_h, 3, 1, 1);
 37    }
 38
 39
 40    void analyze(const Event& event) {
 41
 42      // Identify bjets (signal), light jets (OR) and soft b-jets (veto)
 43      size_t soft_bjets = 0;
 44      Jets bjets, lightjets;
 45      for (const Jet& jet : apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 5*GeV)) {
 46        bool isBjet = jet.bTagged(Cuts::pT > 5*GeV);
 47        if (isBjet) soft_bjets += 1;
 48        if (jet.abseta() < 2.5) {
 49          if ( isBjet && jet.pT() > 25*GeV) bjets += jet;
 50          if (!isBjet && jet.pT() > 20*GeV) lightjets += jet;
 51        }
 52      }
 53      if (soft_bjets != 2) vetoEvent;
 54      if (bjets.size() != 2) vetoEvent;
 55
 56      // Get dressed leptons
 57      DressedLeptons leptons;
 58      for (auto& lep : apply<LeptonFinder>(event, "muons").dressedLeptons()) { leptons.push_back(lep); }
 59      for (auto& lep : apply<LeptonFinder>(event, "elecs").dressedLeptons()) { leptons.push_back(lep); }
 60
 61      // 1. Find which light jets survive OR
 62      for (const auto& lep : leptons) {
 63        idiscard(lightjets, [&](const Jet& jet) {
 64          return deltaR(jet, lep) < 0.2 && (lep.abspid() == PID::ELECTRON || lep.pT()/jet.pT() > 0.7);
 65        });
 66      }
 67
 68      // 2. Find which leptons survive the OR and apply signal selection
 69      for (const auto& jet : (lightjets + bjets)) {
 70        idiscard(leptons, [&](const DressedLepton& lep) {
 71          return lep.pT() < 28*GeV || deltaR(jet, lep) < min(0.4, 0.04+10*GeV/lep.pT());
 72        });
 73      }
 74
 75      if (leptons.size() != 2) vetoEvent;
 76      std::sort(leptons.begin(), leptons.end(), cmpMomByPt);
 77
 78      // Z veto
 79      const size_t nEl = count(leptons, [](const DressedLepton& lep) { return  lep.abspid() == PID::ELECTRON; });
 80      const double mll = (leptons[0].mom() + leptons[1].mom()).mass();
 81      if (nEl != 1 && !(fabs(mll - 91*GeV) > 15*GeV && mll > 10*GeV)) vetoEvent;
 82
 83      const double m00 = (leptons[0].mom() + bjets[0].mom()).mass();
 84      const double m10 = (leptons[1].mom() + bjets[0].mom()).mass();
 85      const double m01 = (leptons[0].mom() + bjets[1].mom()).mass();
 86      const double m11 = (leptons[1].mom() + bjets[1].mom()).mass();
 87      double minimax = min( max(m00,m11), max(m10,m01) );
 88      minimax = min(minimax, 419.); // put overflow in the last bin
 89      _h->fill(minimax/GeV);
 90    }
 91
 92
 93    /// Finalise
 94    void finalize() {
 95      normalize(_h);
 96    }
 97
 98
 99  private:
100
101    /// Histogram
102    Histo1DPtr _h;
103
104  };
105
106
107  RIVET_DECLARE_PLUGIN(ATLAS_2018_I1677498);
108
109}