rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2022_I2170533

Measurements of jet multiplicity and jet transverse momentum in multijet events in proton-proton collisions at \sqrt{s} = 13 TeV
Experiment: CMS (LHC)
Inspire ID: 2170533
Status: VALIDATED
Authors:
  • cms-pag-conveners-smp@cern.ch
  • Luis Ignacio Estevez Banos
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • pp to jets at $\sqrt{s}=13$ TeV. Data collected by CMS during the year 2016.

Multijet events at large transverse momentum (pT) are measured at $\sqrt{s}$ = 13 TeV using data recorded with the CMS detector at the LHC, corresponding to an integrated luminosity of 36.3 fb-1. The multiplicity of jets with pT > 50 GeV that are produced in association with a high-pT dijet system is measured in various ranges of the pT of the jet with the highest transverse momentum and as a function of the azimuthal angle difference between the two highest pT jets in the dijet system. The differential production cross sections are measured as a function of the transverse momenta of the four highest pT jets. The measurements are compared with leading and next-to-leading order matrix element calculations supplemented with simulations of parton shower, hadronization, and multiparton interactions. In addition, the measurements are compared with next-to-leading order matrix element calculations combined with transverse-momentum dependent parton densities and transverse-momentum dependent parton shower.

Source code: CMS_2022_I2170533.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
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Projections/FinalState.hh"

namespace Rivet {

/// @brief Multijet events at 13 TeV

class CMS_2022_I2170533 : public Analysis {
   public:
    /// Constructor
    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2022_I2170533);

    /// @name Analysis methods
    ///@{

    /// Book histograms and initialise projections before the run
    void init() {
        FinalState fs;
        FastJets akt(fs, FastJets::ANTIKT, 0.4);
        declare(akt, "antikT");

        Histo1DPtr dummy;
        _h_Mult_ptmax200.add(0., 150., book(dummy, 1, 1, 1));
        _h_Mult_ptmax200.add(150., 170., book(dummy, 2, 1, 1));
        _h_Mult_ptmax200.add(170., 180., book(dummy, 3, 1, 1));
        _h_Mult_ptmax400.add(0., 150., book(dummy, 4, 1, 1));
        _h_Mult_ptmax400.add(150., 170., book(dummy, 5, 1, 1));
        _h_Mult_ptmax400.add(170., 180., book(dummy, 6, 1, 1));
        _h_Mult_ptmax800.add(0., 150., book(dummy, 7, 1, 1));
        _h_Mult_ptmax800.add(150., 170., book(dummy, 8, 1, 1));
        _h_Mult_ptmax800.add(170., 180., book(dummy, 9, 1, 1));

        // Pt of the first 4 jets
        book(_h_pT1_n1, 10, 1, 1);
        book(_h_pT2_n2, 11, 1, 1);
        book(_h_pT3_n3, 12, 1, 1);
        book(_h_pT4_n4, 13, 1, 1);
    }

    /// Perform the per-event analysis
    void analyze(const Event& event) {
        // Preselection cuts for the jets
        const Jets& jets = apply<FastJets>(event, "antikT").jetsByPt(Cuts::pT > 20 * GeV && Cuts::absrap < 3.2);  // rapidity and pt preselection |y| < 3.2 pT > 20 GeV

        int njet = 0;  // jet counting

        // cuts after preselection
        if (jets.size() < 2) vetoEvent;                                             // dijet cut
        if ((jets[0].pT() < 200.) or (jets[1].pT() < 100.)) vetoEvent;              // pt cut on 2 leading jets
        if ((fabs(jets[0].rap()) > 2.5) or (fabs(jets[1].rap()) > 2.5)) vetoEvent;  // |y| < 2.5 cut on leading jets

        double dphi = deltaPhi(jets[0].phi(), jets[1].phi()) / degree;

        for (const Jet& j : jets) {
            if (j.pT() > 50. && fabs(j.rap()) < 2.5) njet = njet + 1;  // Cuts on the extrajets
        }

        if (njet > 7) njet = 7;  // Last bin in Multiplicity is inclusive

        if (jets[0].pT() > 200 && jets[0].pT() <= 400) _h_Mult_ptmax200.fill(dphi, njet);
        if (jets[0].pT() > 400 && jets[0].pT() <= 800) _h_Mult_ptmax400.fill(dphi, njet);
        if (jets[0].pT() > 800) _h_Mult_ptmax800.fill(dphi, njet);

        // Filling dijet 1,2 pT (The events are already selected with leadin jetpT > 200 GeV and 2nd jet pT > 100 GeV)
        _h_pT1_n1->fill(jets[0].pT());
        _h_pT2_n2->fill(jets[1].pT());

        // Filling extra jets pT (jet3 and jet4)
        if (njet > 2) _h_pT3_n3->fill(jets[2].pT());
        if (njet > 3) _h_pT4_n4->fill(jets[3].pT());
    }

    /// Normalise histograms etc., after the run
    void finalize() {
        for (Histo1DPtr _h_Mult_ptmax200 : _h_Mult_ptmax200.histos()) {
            scale(_h_Mult_ptmax200, crossSection() / picobarn / sumW());
        }
        for (Histo1DPtr _h_Mult_ptmax400 : _h_Mult_ptmax400.histos()) {
            scale(_h_Mult_ptmax400, crossSection() / picobarn / sumW());
        }
        for (Histo1DPtr _h_Mult_ptmax800 : _h_Mult_ptmax800.histos()) {
            scale(_h_Mult_ptmax800, crossSection() / picobarn / sumW());
        }

        for (Histo1DPtr _h : {_h_pT1_n1, _h_pT2_n2, _h_pT3_n3, _h_pT4_n4}) {
            scale(_h, crossSection() / picobarn / sumW());
        }
    }

    ///@}

    /// @name Histograms
    ///@{

    BinnedHistogram _h_Mult_ptmax200, _h_Mult_ptmax400, _h_Mult_ptmax800;
    Histo1DPtr _h_pT1_n1, _h_pT2_n2, _h_pT3_n3, _h_pT4_n4;

    ///@}
};

RIVET_DECLARE_PLUGIN(CMS_2022_I2170533);

}  // namespace Rivet