rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2018_I1707015

Measurement of normalized differential distributions in ttbar+photon production at 13 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1707015
Status: VALIDATED
Authors:
  • Yichen Li
  • Andy Buckley
  • Leah Spiedel Johnson
  • Neil Warrack
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • top-quark pair production associated with a photon

Differential cross-sections for the production of a top-quark pair in association with a photon are measured with proton-proton collision data corresponding to an integrated luminosity of 36.1 fb^-1, collected by the ATLAS detector at the LHC in 2015 and 2016 at a centre-of-mass energy of 13 TeV. The measurements are performed in single-lepton and dilepton final states in a fiducial volume. Events with exactly one photon, one or two leptons, a channel-dependent minimum number of jets, and at least one b-jet are selected. Neural network algorithms are used to separate the signal from the backgrounds. The differential cross-sections are measured as a function of photon transverse momentum, photon absolute pseudorapidity, and angular distance between the photon and its closest lepton in both channels, as well azimuthal opening angle and absolute pseudorapidity difference between the two leptons in the dilepton channel.

Source code: ATLAS_2018_I1707015.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/ChargedFinalState.hh"
  4#include "Rivet/Projections/LeptonFinder.hh"
  5#include "Rivet/Projections/FastJets.hh"
  6#include "Rivet/Projections/PromptFinalState.hh"
  7#include "Rivet/Projections/VetoedFinalState.hh"
  8
  9namespace Rivet {
 10
 11
 12  /// @brief ttbar + gamma at 13 TeV
 13  class ATLAS_2018_I1707015 : public Analysis {
 14  public:
 15
 16
 17    // Constructor
 18    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2018_I1707015);
 19
 20
 21    // Book histograms and initialise projections before the run
 22    void init() {
 23
 24      // Set default running mode to 3 (all)
 25      _mode = 3;
 26
 27      // Get running mode
 28      if ( getOption("LMODE") == "SINGLE" )   _mode = 1;
 29      if ( getOption("LMODE") == "DILEPTON" ) _mode = 2;
 30      if ( getOption("LMODE") == "ALL" )      _mode = 3;
 31
 32      // All final state particles
 33      const FinalState fs;
 34
 35      // Charged particles for signal photon isolation
 36      ChargedFinalState cfs(fs);
 37      declare(cfs, "CFS");
 38
 39      // Signal photons
 40      PromptFinalState photons(Cuts::abspid == PID::PHOTON && Cuts::pT > 20*GeV && Cuts::abseta < 2.37, TauDecaysAs::PROMPT);
 41      declare(photons, "Photons");
 42
 43      // Leptons
 44      PromptFinalState leptons(Cuts::abspid == PID::MUON || Cuts::abspid == PID::ELECTRON, TauDecaysAs::PROMPT);
 45
 46      // Dress the leptons
 47      FinalState dressPhotons(Cuts::abspid == PID::PHOTON);
 48      Cut lepCuts = (Cuts::abseta < 2.5) && (Cuts::pT > 25*GeV);
 49      LeptonFinder dressedLeptons(leptons, dressPhotons, 0.1, lepCuts);
 50      declare(dressedLeptons, "Leptons");
 51
 52      // Jet alg input
 53      VetoedFinalState vfs(fs);
 54
 55      // Remove prompt invisibles from jet input
 56      VetoedFinalState invis_fs(fs);
 57      invis_fs.addVetoOnThisFinalState(VisibleFinalState(fs));
 58      PromptFinalState invis_pfs = PromptFinalState(invis_fs, TauDecaysAs::PROMPT);
 59      vfs.addVetoOnThisFinalState(invis_pfs);
 60
 61      // Remove prompt dressed muons (muons + associated photons) from jet input
 62      PromptFinalState muons(Cuts::abspid == PID::MUON, TauDecaysAs::PROMPT);
 63      LeptonFinder dressedmuons(muons, dressPhotons, 0.1);
 64      vfs.addVetoOnThisFinalState(dressedmuons);
 65
 66      // Jet clustering
 67      FastJets jets(vfs, JetAlg::ANTIKT, 0.4, JetMuons::ALL, JetInvisibles::ALL);
 68      declare(jets, "Jets");
 69
 70      // Book histograms
 71      if ( _mode == 1 or _mode == 3 ){  // Single lepton channel
 72        book(_h["sl_ph_pt"],   3,1,1);  // photon pT
 73        book(_h["sl_ph_eta"],  4,1,1);  // photon eta
 74        book(_h["sl_ph_l_dR"], 5,1,1);  // photon-lepton dR
 75      }
 76      if ( _mode == 2 or _mode == 3 ){  // Dilepton channel
 77        book(_h["dl_ph_pt"],   6,1,1);  // photon pT
 78        book(_h["dl_ph_eta"],  7,1,1);  // photon eta
 79        book(_h["dl_ph_l_dR"], 8,1,1);  // min photon-lepton dR
 80        book(_h["dl_l_dEta"],  9,1,1);  // lepton-lepton dEta
 81        book(_h["dl_l_dPhi"],  10,1,1); // lepton-lepton dPhi
 82      }
 83    }
 84
 85
 86    void analyze(const Event& event) {
 87
 88      // Fetch objects
 89      const DressedLeptons& leptons = apply<LeptonFinder>(event, "Leptons").dressedLeptons();
 90      Particles photons = apply<PromptFinalState>(event, "Photons").particles();
 91      ChargedFinalState charged = apply<ChargedFinalState>(event, "CFS");
 92      Jets jets = apply<JetFinder>(event, "Jets").jetsByPt(Cuts::abseta < 2.5 && Cuts::pT > 25*GeV);
 93
 94      // Immediate veto on events without one good photon
 95      if ( photons.size() != 1 ) vetoEvent;
 96      const Particle& photon = photons[0];
 97
 98      // Veto event if photon too close to a lepton
 99      for (const DressedLepton& lep : leptons) {
100        if ( deltaR(lep, photon) < 1.0 ) vetoEvent;
101      }
102
103      // Overlap removel of jets near leptons
104      for (const DressedLepton& lep : leptons) {
105        idiscard(jets, deltaRLess(lep, 0.4));
106      }
107
108      // Overlap removel of jets near isolated photon
109      double conePt = 0.0;
110      Particles photSurround = charged.particles(deltaRLess(photon, 0.3));
111      for (const Particle& p : photSurround) {
112        conePt += p.pt();
113      }
114      if ( conePt / photon.pT() < 0.1 ) idiscard(jets, deltaRLess(photon, 0.4) );
115
116      // Veto event if photon too close to good jets
117      for (const Jet& jet : jets) {
118        if ( deltaR(jet, photon) < 0.4 ) vetoEvent;
119      }
120
121      // Require at least one bjet
122      unsigned int nbjets = 0;
123      for (const Jet& jet : jets) {
124        if ( jet.bTagged(Cuts::pT > 5*GeV) ) nbjets += 1;
125      }
126      if ( nbjets == 0 ) vetoEvent;
127
128      // Fill histos in single lepton channel
129      if ( _mode == 1 || _mode == 3 ) {
130        if ( jets.size() >= 4 && leptons.size() == 1 ) {
131          _h["sl_ph_pt"]->fill(   photon.pT()/GeV );
132          _h["sl_ph_eta"]->fill(  photon.abseta() );
133          _h["sl_ph_l_dR"]->fill( deltaR(leptons[0], photon) );
134        }
135      }
136
137      // Fill histos in dilepton channel
138      if ( _mode == 2 || _mode == 3 ) {
139        if ( jets.size() >= 2 && leptons.size() == 2 ) {
140          double deltaRNew;
141          double deltaRMin = 999.0;
142          for (const DressedLepton& lep : leptons) {
143            deltaRNew = deltaR(lep, photon);
144            if ( deltaRNew < deltaRMin ) deltaRMin = deltaRNew;
145          }
146          _h["dl_ph_pt"]->fill( photon.pT()/GeV );
147          _h["dl_ph_eta"]->fill( photon.abseta() );
148          _h["dl_ph_l_dR"]->fill( deltaRMin );
149          _h["dl_l_dEta"]->fill( deltaEta(leptons[0], leptons[1]) );
150          _h["dl_l_dPhi"]->fill( deltaPhi(leptons[0], leptons[1]) );
151        }
152      }
153    }
154
155
156    void finalize() {
157
158      // Normalise histograms after the run
159      for (auto &hist : _h) {
160        normalize(hist.second, 1.0, true);
161      }
162    }
163
164
165  private:
166    int _mode;
167    map<string, Histo1DPtr> _h;
168
169  };
170
171  RIVET_DECLARE_PLUGIN(ATLAS_2018_I1707015);
172}