rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2014_I1304688

Measurement of jet multiplicity and transverse momentum spectra in top events using full 7 TeV ATLAS dataset
Experiment: ATLAS (LHC)
Inspire ID: 1304688
Status: VALIDATED
Authors:
  • W. H. Bell
  • A. Grohsjean
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • ttbar events with at least one lepton in the ttbar final state at 7 TeV, i.e. both semileptonic and dileptonic decays should be enabled. The tau decay channels also count as leptonic.

Measurement of the differential t$\bar{\mathrm t}$ production cross-section in 7 TeV proton-proton collisions in the single-lepton channel from ATLAS. The data comprise the full 2011 data sample corresponding to an integrated luminosity of $4.6\,\mathrm{fb}^{-1}$. The differential cross-sections are measured as a function of the jet multiplicity for up to eight jets using jet transverse momentum thresholds of 25, 40, 60, and 80 GeV, and as a function of jet transverse momentum up to the fifth leading jet. The results after background subtraction are corrected for all detector effects, within a kinematic range closely matched to the experimental acceptance.

Source code: ATLAS_2014_I1304688.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/VetoedFinalState.hh"
  4#include "Rivet/Projections/InvisibleFinalState.hh"
  5#include "Rivet/Projections/PromptFinalState.hh"
  6#include "Rivet/Projections/LeptonFinder.hh"
  7#include "Rivet/Projections/FastJets.hh"
  8
  9namespace Rivet {
 10
 11
 12  /// @brief ATLAS 7 TeV jets in ttbar events analysis
 13  ///
 14  /// @author W. H. Bell <W.Bell@cern.ch>
 15  /// @author A. Grohsjean <alexander.grohsjean@desy.de>
 16  class ATLAS_2014_I1304688 : public Analysis {
 17  public:
 18
 19    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2014_I1304688);
 20
 21    void init() {
 22      // Eta ranges
 23      /// @todo 1 MeV? Really?
 24      Cut eta_full = Cuts::abseta < 5.0 && Cuts::pT > 1.0*MeV;
 25      Cut eta_lep = Cuts::abseta < 2.5;
 26
 27      // Get photons to dress leptons
 28      FinalState photons(eta_full && Cuts::abspid == PID::PHOTON);
 29
 30      // Projection to find the electrons
 31      PromptFinalState electrons(eta_full && Cuts::abspid == PID::ELECTRON, TauDecaysAs::PROMPT);
 32      LeptonFinder dressedelectrons(electrons, photons, 0.1, eta_lep && Cuts::pT > 25*GeV);
 33      declare(dressedelectrons, "dressedelectrons");
 34      LeptonFinder vetodressedelectrons(electrons, photons, 0.1, eta_lep && Cuts::pT >= 15*GeV);
 35      declare(vetodressedelectrons, "vetodressedelectrons");
 36      LeptonFinder ewdressedelectrons(electrons, photons, 0.1, eta_full);
 37      declare(ewdressedelectrons, "ewdressedelectrons");
 38
 39      // Projection to find the muons
 40      PromptFinalState muons(eta_full && Cuts::abspid == PID::MUON, TauDecaysAs::PROMPT);
 41      LeptonFinder dressedmuons(muons, photons, 0.1, eta_lep && Cuts::pT >= 25*GeV);
 42      declare(dressedmuons, "dressedmuons");
 43      LeptonFinder vetodressedmuons(muons, photons, 0.1, eta_lep && Cuts::pT >= 15*GeV);
 44      declare(vetodressedmuons, "vetodressedmuons");
 45      LeptonFinder ewdressedmuons(muons, photons, 0.1, eta_full);
 46      declare(ewdressedmuons, "ewdressedmuons");
 47
 48      // Projection to find neutrinos and produce MET
 49      InvisibleFinalState neutrinos(OnlyPrompt::YES, TauDecaysAs::PROMPT);
 50      declare(neutrinos, "neutrinos");
 51
 52      // Jet clustering.
 53      VetoedFinalState vfs;
 54      vfs.addVetoOnThisFinalState(ewdressedelectrons);
 55      vfs.addVetoOnThisFinalState(ewdressedmuons);
 56      vfs.addVetoOnThisFinalState(neutrinos);
 57      FastJets jets(vfs, JetAlg::ANTIKT, 0.4);
 58      jets.useInvisibles();
 59      declare(jets, "jets");
 60
 61      // Book histograms
 62      for (size_t i = 0; i < pTcuts.size(); ++i) {
 63        const string name = "mult_"+std::to_string(i);
 64        book(_s[name], i+1, 1, 1);
 65      }
 66      for (size_t i = 0; i < njets; ++i) {
 67        const string name = "jet_"+std::to_string(i);
 68        book(_h[name], i+5, 1, 1);
 69      }
 70    }
 71
 72
 73    void analyze(const Event& event) {
 74
 75      if (_edges.empty()) {
 76        _edges.resize(_s.size());
 77        for (size_t i = 0; i < _edges.size(); ++i) {
 78          const string name = "mult_"+std::to_string(i);
 79          _edges[i] = _s[name]->xEdges();
 80        }
 81      }
 82
 83      // Get the selected objects, using the projections.
 84      const DressedLeptons dressedelectrons = sortByPt(apply<LeptonFinder>(event, "dressedelectrons").dressedLeptons());
 85      const DressedLeptons vetodressedelectrons = apply<LeptonFinder>(event, "vetodressedelectrons").dressedLeptons();
 86
 87      const DressedLeptons dressedmuons = sortByPt(apply<LeptonFinder>(event, "dressedmuons").dressedLeptons());
 88      const DressedLeptons vetodressedmuons = apply<LeptonFinder>(event, "vetodressedmuons").dressedLeptons();
 89
 90      if (dressedelectrons.empty() && dressedmuons.empty())  vetoEvent;
 91      if (dressedelectrons.size()) {
 92        if (vetodressedelectrons.size() > 1 || vetodressedmuons.size())  vetoEvent;
 93      }
 94      if (dressedmuons.size()) {
 95        if (vetodressedmuons.size() > 1 || vetodressedelectrons.size())  vetoEvent;
 96      }
 97
 98      Jets jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 25*GeV && Cuts::abseta < 2.5);
 99      if (jets.size() < 3)  vetoEvent;
100
101
102      // Calculate the missing ET, using the prompt neutrinos only (really?)
103      /// @todo Why not use MissingMomentum?
104      FourMomentum pmet;
105      for (const Particle& p : apply<InvisibleFinalState>(event, "neutrinos").particlesByPt())  pmet += p.momentum();
106      const double met_et = pmet.pT();
107      const double met_phi = pmet.phi();
108      if (met_et <= 30*GeV)  vetoEvent;
109      if (dressedelectrons.size()) {
110        if (transMass(dressedelectrons[0].pT(), dressedelectrons[0].phi(), met_et, met_phi) <= 35*GeV)  vetoEvent;
111      }
112      if (dressedmuons.size()) {
113        if (transMass(dressedmuons[0].pT(), dressedmuons[0].phi(), met_et, met_phi) <= 35*GeV)  vetoEvent;
114      }
115
116      // Check overlap of jets/leptons.
117      size_t jet_ntag = 0;
118      bool overlap = false;
119      for (size_t j1 = 0; j1 < jets.size(); ++j1) {
120        const Jet& jet = jets[j1];
121        // If dR(el,jet) < 0.4 skip the event
122        for (const DressedLepton& el : dressedelectrons) {
123          if (deltaR(jet, el) < 0.4)  overlap = true;
124        }
125        // If dR(mu,jet) < 0.4 skip the event
126        for (const DressedLepton& mu : dressedmuons) {
127          if (deltaR(jet, mu) < 0.4)  overlap = true;
128        }
129        // If dR(jet,jet) < 0.5 skip the event
130        for (size_t j2 = j1+1; j2 < jets.size(); ++j2) {
131          const Jet& jet2 = jets[j2];
132          if (deltaR(jet, jet2) < 0.5)  overlap = true;
133        }
134        // Count the number of b-tags
135        if (jet.bTags().size())  ++jet_ntag;
136      }
137
138      // Remove events with object overlap
139      if (overlap || !jet_ntag)  vetoEvent;
140
141      // Count the jet multiplicity for 25, 40, 60 and 80GeV
142      vector<size_t> jet_n; jet_n.resize(pTcuts.size());
143      for (size_t i = 0; i < pTcuts.size(); ++i) {
144        jet_n[i] = countJets(jets, i);
145        const string name = "mult_" + std::to_string(i);
146        const string& edge = _edges[i][jet_n[i]-3];
147        _s[name]->fill(edge);
148      }
149
150      // Check if the additional pT threshold requirements are passed
151      const bool pass_jetPt = jets.size() > 1 && jets[0].pT() > 50*GeV && jets[1].pT() > 35*GeV;
152      if (!pass_jetPt)  vetoEvent;
153
154
155      // Fill histograms
156      for (size_t i = 0; i < njets; ++i) {
157        const string name = "jet_" + std::to_string(i);
158        if (i > 1 && jet_n[0] <= i)  continue;
159        _h[name]->fill(jets[i].pT()/GeV);
160      }
161    }
162
163
164    void finalize() {
165      // Normalize to cross-section x 0.5 to average lepton flavours
166      const double norm = 0.5*crossSection()/picobarn/sumOfWeights();
167      scale(_h, norm);
168      scale(_s, norm);
169    }
170
171
172
173  private:
174
175
176    /// @name Physics object helper functions
177    /// @{
178
179    double transMass(double ptLep, double phiLep, double met, double phiMet) const {
180      return sqrt(2.0*ptLep*met*(1 - cos(phiLep-phiMet)));
181    }
182
183    size_t countJets(const Jets& jets, size_t thresh) const {
184      size_t jet_n = 0;
185      for (const Jet& jet : jets) {
186        if (jet.pT() > pTcuts[thresh])  ++jet_n;
187      }
188      const vector<size_t> ncutoff{8,7,6,5};
189      return min(jet_n, ncutoff[thresh]);
190    }
191
192    /// @}
193
194
195  private:
196
197    /// @name Objects that are used by the event selection decisions
198    /// @{
199
200    map<string, Histo1DPtr> _h;
201    map<string, BinnedHistoPtr<string>> _s;
202    vector<vector<string>> _edges;
203    const size_t njets = 5;
204    const vector<int> pTcuts{25,40,60,80};
205    /// @}
206
207  };
208
209
210
211  RIVET_DECLARE_PLUGIN(ATLAS_2014_I1304688);
212
213}