rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2017_I1589844

$k_T$ splittings in $Z$ events at 8 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1589844
Status: VALIDATED
Authors:
  • Christian Gutschow
  • Frank Siegert
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • $pp \to Z(\to ee/\mu\mu) +$ jets at 8 TeV

A measurement of the splitting scales occurring in the $k_\text{t}$ jet-clustering algorithm is presented for final states containing a $Z$ boson. The measurement is done using 20.2 fb$^{-1}$ of proton-proton collision data collected at a centre-of-mass energy of $\sqrt{s} = 8$ TeV by the ATLAS experiment at the LHC in 2012. The measurement is based on charged-particle track information, which is measured with excellent precision in the $p_\text{T}$ region relevant for the transition between the perturbative and the non-perturbative regimes. The data distributions are corrected for detector effects, and are found to deviate from state-of-the-art predictions in various regions of the observables. If no mode specified, will try to fill both electron and muon plots. If EL or MU is specified, only the relevant plots will be filled.

Source code: ATLAS_2017_I1589844.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/IdentifiedFinalState.hh"
  5#include "Rivet/Projections/DressedLeptons.hh"
  6#include "Rivet/Projections/FastJets.hh"
  7#include "Rivet/Projections/VetoedFinalState.hh"
  8#include "Rivet/Projections/ChargedFinalState.hh"
  9
 10namespace Rivet {
 11
 12
 13  /// kT splittings in Z events at 8 TeV
 14  class ATLAS_2017_I1589844 : public Analysis {
 15  public:
 16
 17    /// @name Constructors etc.
 18    //@{
 19
 20    /// Constructors
 21    ATLAS_2017_I1589844(const string name="ATLAS_2017_I1589844",
 22                        const string ref_data="ATLAS_2017_I1589844") : Analysis(name) {
 23      setRefDataName(ref_data);
 24    }
 25
 26    //@}
 27
 28
 29    /// @name Analysis methods
 30    //@{
 31
 32    /// Book histograms and initialise projections before the run
 33    void init() {
 34
 35      // Get options from the new option system
 36      _mode = 0;
 37      if ( getOption("LMODE") == "EL" ) _mode = 1;
 38      if ( getOption("LMODE") == "MU" ) _mode = 2;
 39
 40      const FinalState fs;
 41
 42      const Cut cuts_mu = (Cuts::pT > 25*GeV) && (Cuts::abseta < 2.4);
 43      const Cut cuts_el = Cuts::pT > 25*GeV && (Cuts::abseta <= 1.37 || (Cuts::abseta >= 1.52 && Cuts::abseta < 2.47));
 44
 45      IdentifiedFinalState bare_mu(fs);
 46      bare_mu.acceptIdPair(PID::MUON);
 47      IdentifiedFinalState bare_el(fs);
 48      bare_el.acceptIdPair(PID::ELECTRON);
 49      const DressedLeptons muons(fs, bare_mu, 0.1, cuts_mu, true);
 50      const DressedLeptons elecs(fs, bare_el, 0.1, cuts_el, true);
 51      declare(muons, "muons");
 52      declare(elecs, "elecs");
 53
 54      const ChargedFinalState cfs(Cuts::abseta < 2.5 && Cuts::pT > 0.4*GeV);
 55      VetoedFinalState jet_fs(cfs);
 56      jet_fs.addVetoOnThisFinalState(muons);
 57      jet_fs.addVetoOnThisFinalState(elecs);
 58      declare(FastJets(jet_fs, FastJets::KT, 0.4), "Kt04Jets");
 59      declare(FastJets(jet_fs, FastJets::KT, 1.0), "Kt10Jets");
 60
 61      VetoedFinalState jet_fs_all(FinalState(Cuts::abseta < 2.5 && Cuts::pT > 0.4*GeV));
 62      jet_fs_all.addVetoOnThisFinalState(muons);
 63      jet_fs_all.addVetoOnThisFinalState(elecs);
 64      FastJets jetpro04_all(jet_fs_all, FastJets::KT, 0.4);
 65      jetpro04_all.useInvisibles();
 66      declare(jetpro04_all, "Kt04Jets_all");
 67      FastJets jetpro10_all(jet_fs_all, FastJets::KT, 1.0);
 68      jetpro10_all.useInvisibles();
 69      declare(jetpro10_all, "Kt10Jets_all");
 70
 71      // Histograms with data binning
 72      _ndij = 8;
 73      for (size_t i = 0; i < _ndij; ++i) {
 74        if (_mode == 0 || _mode == 1) {
 75          string label = "el_d" + to_str(i) + "_kT4";
 76          book(_h[label], i + 1, 1, 1);
 77          book(_h[label + "_all"], i + 1, 1, 5);
 78          label = "el_d" + to_str(i) + "_kT10";
 79          book(_h[label], i + 1, 1, 3);
 80          book(_h[label + "_all"], i + 1, 1, 7);
 81        }
 82        if (_mode == 0 || _mode == 2) {
 83          string label = "mu_d" + to_str(i) + "_kT4";
 84          book(_h[label], i + 1, 1, 2);
 85          book(_h[label + "_all"], i + 1, 1, 6);
 86          label = "mu_d" + to_str(i) + "_kT10";
 87          book(_h[label], i + 1, 1, 4);
 88          book(_h[label + "_all"], i + 1, 1, 8);
 89        }
 90      }
 91    }
 92
 93
 94    /// Perform the per-event analysis
 95    void analyze(const Event& e) {
 96
 97      // Check we have a Z candidate:
 98      const vector<DressedLepton>& muons = apply<DressedLeptons>(e, "muons").dressedLeptons();
 99      const vector<DressedLepton>& elecs = apply<DressedLeptons>(e, "elecs").dressedLeptons();
100
101      bool e_ok = elecs.size() == 2 && muons.empty();
102      bool m_ok = elecs.empty() && muons.size() == 2;
103      if (_mode == 0 && !e_ok && !m_ok) vetoEvent;
104      if (_mode == 1 && !e_ok )  vetoEvent;
105      if (_mode == 2 && !m_ok )  vetoEvent;
106
107      string lep_type = elecs.size()? "el_" : "mu_";
108
109      const vector<DressedLepton>& leptons = elecs.size()? elecs : muons;
110      if (leptons[0].charge()*leptons[1].charge() > 0) vetoEvent;
111
112      const double dilepton_mass = (leptons[0].momentum() + leptons[1].momentum()).mass();
113      if (!inRange(dilepton_mass, 71*GeV, 111*GeV)) vetoEvent;
114
115      // Get kT splitting scales (charged particles only)
116      const FastJets& jetpro04 = applyProjection<FastJets>(e, "Kt04Jets");
117      const shared_ptr<fastjet::ClusterSequence> seq04 = jetpro04.clusterSeq();
118      for (size_t i = 0; i < min(_ndij, (size_t)seq04->n_particles()); ++i) {
119        const double dij = sqrt(seq04->exclusive_dmerge_max(i))/GeV;
120        if (dij <= 0.0) continue;
121        const string label = lep_type + "d" + to_str(i) + "_kT4";
122        _h[label]->fill(dij);
123      }
124      const FastJets& jetpro10 = applyProjection<FastJets>(e, "Kt10Jets");
125      const shared_ptr<fastjet::ClusterSequence> seq10 = jetpro10.clusterSeq();
126      for (size_t i = 0; i < min(_ndij, (size_t)seq10->n_particles()); ++i) {
127        const double dij = sqrt(seq10->exclusive_dmerge_max(i))/GeV;
128        if (dij <= 0.0) continue;
129        const string label = lep_type + "d" + to_str(i) + "_kT10";
130        _h[label]->fill(dij);
131      }
132
133      // Get kT splitting scales (all particles)
134      const FastJets& jetpro04_all = applyProjection<FastJets>(e, "Kt04Jets_all");
135      const shared_ptr<fastjet::ClusterSequence> seq04_all = jetpro04_all.clusterSeq();
136      for (size_t i = 0; i < min(_ndij, (size_t)seq04_all->n_particles()); ++i) {
137        const double dij = sqrt(seq04_all->exclusive_dmerge_max(i))/GeV;
138        if (dij <= 0.0) continue;
139        const string label = lep_type + "d" + to_str(i) + "_kT4_all";
140        _h[label]->fill(dij);
141      }
142      const FastJets& jetpro10_all = applyProjection<FastJets>(e, "Kt10Jets_all");
143      const shared_ptr<fastjet::ClusterSequence> seq10_all = jetpro10_all.clusterSeq();
144      for (size_t i = 0; i < min(_ndij, (size_t)seq10_all->n_particles()); ++i) {
145        const double dij = sqrt(seq10_all->exclusive_dmerge_max(i))/GeV;
146        if (dij <= 0.0) continue;
147        const string label = lep_type + "d" + to_str(i) + "_kT10_all";
148        _h[label]->fill(dij);
149      }
150
151    }
152
153
154    /// Normalise histograms etc., after the run
155    void finalize() {
156      const double sf = crossSectionPerEvent();
157      for (auto& kv : _h) scale(kv.second, sf);
158    }
159
160    //@}
161
162
163  protected:
164
165    // Data members like post-cuts event weight counters go here
166    size_t _mode, _ndij;
167
168
169  private:
170
171    // Histograms
172    map<string, Histo1DPtr> _h;
173
174  };
175
176
177  RIVET_DECLARE_PLUGIN(ATLAS_2017_I1589844);
178
179}