Rivet Analyses Reference

CMS_2017_I1598460

Triple-differential dijet pT cross section and PDF constraints at 8 TeV
Experiment: CMS (LHC)
Inspire ID: 1598460
Status: VALIDATED
Authors:
  • Klaus Rabbertz
  • Georg Sieber
References:
  • Eur.Phys.J. C77 (2017) no.11, 746
  • arXiv: 1705.02628
Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • pp collision, 8 TeV, jets with anti-kt algorithm with a jet size parameter of $R = 0.7$; only events with at least two jets up to $|y| = 5.0$ are selected and the two jets leading in $p_T$ are required to have transverse momenta greater than 50 GeV and $|y| < 3.0$.

A measurement is presented of the triple-differential dijet cross section at a centre-of-mass energy of 8 TeV using 19.7 fb$^-1$ of data collected with the CMS detector in proton-proton collisions at the LHC. The cross section is measured as a function of the average transverse momentum, half the rapidity separation, and the boost of the two leading jets in the event. The cross section is corrected for detector effects and compared to calculations in perturbative quantum chromodynamics at next-to-leading order accuracy, complemented with electroweak and nonperturbative corrections. New constraints on parton distribution functions are obtained and the inferred value of the strong coupling constant is $\alpha_s(M_{Z})$ = 0.1199 $\pm$ 0.0015 (exp)$^{+0.0031}_{-0.0020}$ (theo), where M$_{Z}$ is the mass of the Z boson.

Source code: CMS_2017_I1598460.cc
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/FastJets.hh"

namespace Rivet {


  class CMS_2017_I1598460 : public Analysis {
  public:

    /// Constructor
    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2017_I1598460);


    /// @name Analysis methods
    //@{

    /// Book histograms and initialise projections before the run
    void init() {

      const FinalState fs;
      declare(fs, "FinalState");
      FastJets fj07(fs, FastJets::ANTIKT, 0.7);
      declare(fj07, "Jets");
      /// @todo Book histograms here, e.g.:
      for (int i = 0; i < 6; i++) {
        Histo1DPtr tmp; _h_ybys.push_back(book(tmp, 2*i+1,1,1));
      }
    }


    /// Perform the per-event analysis
    void analyze(const Event& event) {
      const FastJets& fj = apply<FastJets>(event,"Jets");
      const Jets& jets = fj.jetsByPt(Cuts::pt > 50*GeV && Cuts::absrap < 5);

      // Require two jets
      if (jets.size() < 2) vetoEvent;

      // Veto events if one of two leading jets |y|>3.0, otherwise
      // the subleading jets can become the leading jets through jet selection
      if (jets[0].absrap() > 3 || jets[1].absrap() > 3) vetoEvent;

      double ystar = 0.5 * std::abs(jets[0].rap() - jets[1].rap());
      double yboost = 0.5 * std::abs(jets[0].rap() + jets[1].rap());
      double ptavg = 0.5 * (jets[0].pT() + jets[1].pT());

      // Compute index of histogram to be filled: yb0ys0 --> 0 yb0ys1 --> 1 ...
      size_t i = (size_t)yboost;
      size_t j = (size_t)ystar;
      size_t idx = j + 3*i - i*(i-1)/2;

      _h_ybys[idx]->fill(ptavg/GeV);
    }


    /// Normalise histograms etc., after the run
    void finalize() {
      scale(_h_ybys, crossSection()/sumOfWeights());
    }

    //@}


    vector<Histo1DPtr> _h_ybys;

  };


  RIVET_DECLARE_PLUGIN(CMS_2017_I1598460);

}