rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2016_I1421646

Dijet azimuthal decorrelations in $pp$ collisions at $\sqrt{s} = 8$ TeV
Experiment: CMS (LHC)
Inspire ID: 1421646
Status: VALIDATED
Authors:
  • Giannis Flouris
  • Panos Kokkas
  • Evangelos Paradas
  • Klaus Rabbertz
  • Georg Sieber
References: Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • QCD at $\sqrt{s} = 8 \text{TeV}$

A measurement of the decorrelation of azimuthal angles between the two jets with the largest transverse momenta is presented for seven regions of leading jet transverse momentum up to 2.2 TeV. The analysis is based on the proton-proton collision data collected with the CMS experiment at a centre-of-mass energy of 8 TeV corresponding to an integrated luminosity of 19.7/fb. The dijet azimuthal decorrelation is caused by the radiation of additional jets and probes the dynamics of multijet production. The results are compared to fixed-order predictions of perturbative quantum chromodynamics (QCD), and to simulations using Monte Carlo event generators that include parton showers, hadronization, and multiparton interactions. Event generators with only two outgoing high transverse momentum partons fail to describe the measurement, even when supplemented with next-to-leading-order QCD corrections and parton showers. Much better agreement is achieved when at least three outgoing partons are complemented through either next-to-leading-order predictions or parton showers. This observation emphasizes the need to improve predictions for multijet production.

Source code: CMS_2016_I1421646.cc
 1// -*- C++ -*-
 2#include "Rivet/Analysis.hh"
 3#include "Rivet/Projections/FinalState.hh"
 4#include "Rivet/Projections/FastJets.hh"
 5
 6namespace Rivet {
 7
 8
 9  /// CMS azimuthal decorrelations at 8 TeV
10  class CMS_2016_I1421646 : public Analysis {
11  public:
12
13    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2016_I1421646);
14
15
16    /// Book projections and histograms
17    void init() {
18
19      FastJets akt(FinalState(), JetAlg::ANTIKT, 0.7);
20      declare(akt, "antikT");
21
22      book(_h_deltaPhi, {200., 300., 400., 500., 700., 900., 1100., 4000.});
23      for (auto& b : _h_deltaPhi->bins()) {
24        book(b, b.index(), 1, 1);
25      }
26    }
27
28
29    /// Per-event analysis
30    void analyze(const Event & event) {
31
32      const Jets& jets = apply<JetFinder>(event, "antikT").jetsByPt(Cuts::absrap < 5.0 && Cuts::pT > 100*GeV);
33      if (jets.size() < 2) vetoEvent;
34      if (jets[0].pT() < 200*GeV) vetoEvent;
35      if (jets[0].absrap() > 2.5 || jets[1].absrap() > 2.5) vetoEvent;
36
37      const double dphi = deltaPhi(jets[0].phi(), jets[1].phi());
38      _h_deltaPhi->fill(jets[0].pT(), dphi);
39    }
40
41
42    /// Scale histograms
43    void finalize() {
44      normalize(_h_deltaPhi);
45    }
46
47
48  private:
49
50    Histo1DGroupPtr _h_deltaPhi;
51
52  };
53
54
55  RIVET_DECLARE_PLUGIN(CMS_2016_I1421646);
56
57}