rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2012_I1082936

Inclusive jet and dijet cross sections at 7 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1082936
Status: VALIDATED
Authors:
  • Holger Schulz
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • QCD jet production with a minimum leading jet pT of 30 GeV and minimum second jet pT of 20 GeV at 7 TeV.

Inclusive jet and dijet cross sections have been measured in proton-proton collisions at a centre-of-mass energy of 7 TeV using the ATLAS detector at the Large Hadron Collider. The cross sections were measured using jets clustered with the anti-kT algorithm with parameters R=0.4 and R=0.6. These measurements are based on the 2010 data sample, consisting of a total integrated luminosity of 37 inverse picobarns. Inclusive jet double-differential cross sections are presented as a function of jet transverse momentum, in bins of jet rapidity. Dijet double-differential cross sections are studied as a function of the dijet invariant mass, in bins of half the rapidity separation of the two leading jets. The measurements are performed in the jet rapidity range $|y|<4.4$, covering jet transverse momenta from 20 GeV to 1.5 TeV and dijet invariant masses from 70 GeV to 5 TeV. This is the successor analysis of ATLAS_2010_S8817804

Source code: ATLAS_2012_I1082936.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FastJets.hh"
  4#include "Rivet/Tools/BinnedHistogram.hh"
  5#include "Rivet/Projections/FinalState.hh"
  6
  7namespace Rivet {
  8
  9
 10  class ATLAS_2012_I1082936 : public Analysis {
 11  public:
 12
 13    /// @name Constructors etc.
 14    //@{
 15
 16    /// Constructor
 17    ATLAS_2012_I1082936()
 18      : Analysis("ATLAS_2012_I1082936")
 19    {
 20    }
 21
 22    //@}
 23
 24
 25  public:
 26
 27    /// @name Analysis methods
 28    //@{
 29
 30    /// Book histograms and initialise projections before the run
 31    void init() {
 32
 33      const FinalState fs;
 34      declare(fs,"FinalState");
 35
 36      FastJets fj04(fs,  FastJets::ANTIKT, 0.4);
 37      fj04.useInvisibles();
 38      declare(fj04, "AntiKT04");
 39
 40      FastJets fj06(fs,  FastJets::ANTIKT, 0.6);
 41      fj06.useInvisibles();
 42      declare(fj06, "AntiKT06");
 43
 44
 45      // Histogram booking copied from the previous analysis
 46      double ybins[] = { 0.0, 0.3, 0.8, 1.2, 2.1, 2.8, 3.6, 4.4 };
 47      double ystarbins[] = { 0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.4};
 48
 49      size_t ptDsOffset(0), massDsOffset(2);
 50      for (size_t alg = 0; alg < 2; ++alg) {
 51        for (size_t i = 0; i < 7; ++i) {
 52          Histo1DPtr tmp; 
 53          _pThistos[alg].add(ybins[i], ybins[i+1], book(tmp, 1 + ptDsOffset, 1, i+1));
 54        }
 55        ptDsOffset += 1;
 56
 57        for (size_t i = 0; i < 9; ++i) {
 58          Histo1DPtr tmp;
 59          _mass[alg].add(ystarbins[i], ystarbins[i+1], book(tmp, 1 + massDsOffset, 1, i+1));
 60        }
 61        massDsOffset += 1;
 62      }
 63    }
 64
 65    /// Perform the per-event analysis
 66    void analyze(const Event& event) {
 67      Jets jetAr[2];
 68      jetAr[AKT6] = apply<FastJets>(event, "AntiKT06").jetsByPt(20*GeV);
 69      jetAr[AKT4] = apply<FastJets>(event, "AntiKT04").jetsByPt(20*GeV);
 70
 71      // Loop over jet "radii" used in analysis
 72      for (size_t alg = 0; alg < 2; ++alg) {
 73        // Identify dijets
 74        vector<FourMomentum> leadjets;
 75        for (const Jet& jet : jetAr[alg]) {
 76          const double pT = jet.pT();
 77          const double absy = jet.absrap();
 78          _pThistos[alg].fill(absy, pT/GeV);
 79
 80          if (absy < 4.4 && leadjets.size() < 2) {
 81            if (leadjets.empty() && pT < 30*GeV) continue;
 82            leadjets.push_back(jet.momentum());
 83          }
 84        }
 85        // Make sure we have a leading jet with pT >30 GeV and a second to leading jet with pT>20 GeV
 86        if (leadjets.size() < 2) {
 87          MSG_DEBUG("Could not find two suitable leading jets");
 88          continue;
 89        }
 90
 91        const double y1 = leadjets[0].rapidity();
 92        const double y2 = leadjets[1].rapidity();
 93        const double ystar = fabs(y1-y2)/2.;
 94        const double m = (leadjets[0] + leadjets[1]).mass();
 95        // Fill mass histogram
 96        _mass[alg].fill(ystar, m/TeV);
 97      }
 98    }
 99
100
101    /// Normalise histograms etc., after the run
102    void finalize() {
103      for (size_t alg = 0; alg < 2; ++alg) {
104        // factor 0.5 needed because it is differential in dy and not d|y|
105        _pThistos[alg].scale(0.5*crossSectionPerEvent()/picobarn, this);
106        _mass[alg].scale(crossSectionPerEvent()/picobarn, this);
107}
108    }
109
110    //@}
111
112
113  private:
114
115    // Data members like post-cuts event weight counters go here
116
117    enum Alg { AKT4=0, AKT6=1 };
118
119  private:
120
121    /// The inclusive pT spectrum for akt6 and akt4 jets (array index is jet type from enum above)
122    BinnedHistogram _pThistos[2];
123
124    /// The di-jet mass spectrum binned in rapidity for akt6 and akt4 jets (array index is jet type from enum above)
125    BinnedHistogram _mass[2];
126  };
127
128  // The hook for the plugin system
129  RIVET_DECLARE_PLUGIN(ATLAS_2012_I1082936);
130
131}