rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2014_I1326641

3-jet cross section with 7 TeV data
Experiment: ATLAS (LHC)
Inspire ID: 1326641
Status: VALIDATED
Authors:
  • Aliaksei Hrynevich
  • Pavel Starovoitov
References: Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • QCD jet production with at least three leading jets

Double-differential three-jet production cross-sections have been measured in proton-proton collisions at a centre-of-mass energy of $\sqrt{s} = 7\,$TeV using the ATLAS detector at the Large Hadron Collider. The measurements are presented as a function of the three-jet mass ($m_{jjj}$), in bins of the sum of the absolute rapidity separations between the three leading jets ($\left|Y^\ast\right|$). Invariant masses extending up to 5\,TeV are reached for $8<\left|Y^\ast\right|<10$. These measurements use a sample of data recorded using the ATLAS detector in 2011, which corresponds to an integrated luminosity of $4.51\,\text{fb}^{-1}$. Jets are identified using the anti-$k_\text{t}$ algorithm with two different jet radius parameters, $R=0.4$ and $R=0.6$.

Source code: ATLAS_2014_I1326641.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5#include "Rivet/Tools/BinnedHistogram.hh"
  6
  7namespace Rivet {
  8
  9
 10    class ATLAS_2014_I1326641 : public Analysis {
 11        public:
 12
 13            /// @name Constructors etc.
 14            //@{
 15
 16            /// Constructor
 17            ATLAS_2014_I1326641()
 18                : Analysis("ATLAS_2014_I1326641")
 19            {
 20
 21            }
 22
 23            //@}
 24
 25
 26        public:
 27
 28            /// @name Analysis methods
 29            //@{
 30
 31            /// Book histograms and initialise projections before the run
 32            void init() {
 33                //std::cout << " HELLO ANALYSIS : init " << std::'\n';
 34                const FinalState fs;
 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                double ystarBins[] = { 0.0, 2.0, 4.0, 6.0, 8.0, 10.0 };
 45
 46                size_t massDsOffset(0);
 47                for (size_t alg = 0; alg < 2; ++alg) {
 48                    for (size_t i = 0; i < 5; ++i) {
 49                        Histo1DPtr tmp;
 50                        h_trijet_Mass[alg].add(ystarBins[i], ystarBins[i+1], book(tmp, 1 + massDsOffset, 1, 1));
 51                        massDsOffset += 1;
 52                    }
 53                }
 54            }
 55
 56            /// Perform the per-event analysis
 57            void analyze(const Event& event) {
 58
 59                Jets jetAr[2];
 60                jetAr[AKT4] = apply<FastJets>(event, "AntiKT04").jetsByPt(Cuts::pT > 50*GeV);
 61                jetAr[AKT6] = apply<FastJets>(event, "AntiKT06").jetsByPt(Cuts::pT > 50*GeV);
 62
 63                const size_t nJets = 3;
 64                double ptCut[nJets] = { 150., 100., 50.};
 65
 66                // Loop over jet "radii" used in analysis
 67                for (size_t alg = 0; alg < 2; ++alg) {
 68                    // Identify 3jets
 69                    vector<FourMomentum> leadJets;
 70                    for (const Jet& jet : jetAr[alg]) {
 71                        if (jet.absrap() < 3.0 && leadJets.size() < nJets){
 72                            int filledJets = leadJets.size();
 73                            if (jet.pT() < ptCut[filledJets])  continue;
 74                            leadJets.push_back(jet.momentum());
 75                        }
 76                    }
 77
 78                    if (leadJets.size() < nJets) {
 79                        MSG_DEBUG("Could not find three suitable leading jets");
 80                        continue;
 81                    }
 82
 83                    const double y1 = leadJets[0].rapidity();
 84                    const double y2 = leadJets[1].rapidity();
 85                    const double y3 = leadJets[2].rapidity();
 86
 87                    const double yStar = fabs(y1-y2) + fabs(y2-y3) + fabs(y1-y3);
 88                    const double m = (leadJets[0] + leadJets[1] + leadJets[2]).mass();
 89                    h_trijet_Mass[alg].fill(yStar, m, 1.0);
 90                }
 91            }
 92
 93
 94            /// Normalise histograms etc., after the run
 95            void finalize() {
 96
 97                const double sf( 2.0 * crossSection() / sumOfWeights() );
 98                for (size_t alg = 0; alg < 2; ++alg) {
 99                    h_trijet_Mass[alg].scale(sf, this);
100                }
101
102            }
103
104            //@}
105
106
107        private:
108
109            // Data members like post-cuts event weight counters go here
110            enum Alg { AKT4=0, AKT6=1 };
111
112        private:
113
114            // The 3 jets mass spectrum for anti-kt 4 and anti-kt 6 jets (array index is jet type from enum above)
115            BinnedHistogram  h_trijet_Mass[2];
116    };
117
118    // The hook for the plugin system
119    RIVET_DECLARE_PLUGIN(ATLAS_2014_I1326641);
120}