Rivet analyses referenceCMS_2024_I2780732Measurement of Rdphi observable for the determination of the strong coupling at 13 TeVExperiment: CMS (LHC) Inspire ID: 2780732 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
A measurement is presented of the ratio observable $R_{\Delta\phi}(p_\mathrm{T})$ that provides a measure of the azimuthal correlations among jets with large transverse momentum $p_\mathrm{T}$. The $R_{\Delta\phi}(p_\mathrm{T})$ variable is measured in multijet events over the $p_\mathrm{T} = 360\mathord{-}3170 \mathrm{GeV}$ range based on data collected by the CMS experiment in proton-proton collisions at a centre-of-mass energy of 13 TeV, corresponding to an integrated luminosity of $134 \mathrm{fb}^{-1}$. The results are compared with predictions from Monte Carlo parton-shower event generator simulations, as well as with fixed-order perturbative quantum chromodynamics (pQCD) predictions at next-to-leading-order (NLO) accuracy obtained with different parton distribution functions (PDFs) and corrected for nonperturbative and electroweak effects. Data and theory agree within uncertainties. From the comparison of the measured distribution with the pQCD prediction obtained with the NNPDF3.1 NLO PDFs, the strong coupling constant at the Z boson mass scale is $\alpha_{\text{S}}(m_\text{Z}) = 0.1177 \pm 0.0013\, \text{(exp)} _{-0.0073}^{+0.0116} \, \text{(theo)} = 0.1177_{-0.0074}^{+0.0117}$, where the total uncertainty is dominated by the scale dependence of the fixed-order predictions. A test of the running of $\alpha_S(Q)$ in the TeV region shows no deviation from the expected NLO pQCD behaviour. Source code: CMS_2024_I2780732.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 /// @brief Rdphi observable for the determination of the strong coupling at 13 TeV
9 class CMS_2024_I2780732 : public Analysis {
10 public:
11
12 /// Constructor
13 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2024_I2780732);
14
15 /// Book histograms and initialise projections before the run
16 void init() {
17
18 // Initialise and register projections
19 FinalState fs;
20 FastJets akt(fs, JetAlg::ANTIKT, 0.7);
21 declare(akt, "antikT");
22
23 // Book histograms
24 // Denominator and numerator of Rdphi observable
25 book(_h_Inclusive_Denominator, "_Inclusive_Denominator", refData(2, 1, 1));
26 book(_h_Dphi_Numerator, "_Inclusive_Numerator", refData(2, 1, 1));
27 book(_h_Rdphi, 2, 1, 1);
28
29 }
30
31
32 /// Perform the per-event analysis
33 void analyze(const Event& event) {
34
35 // Select jets with pt > 50 GeV and absolute rapidity < 2.5
36 const Jets& jets = apply<FastJets>(event, "antikT").jetsByPt(Cuts::absrap < 2.5 && Cuts::pT > 50.*GeV);
37
38 // Loop over all jets and fill the denominator (inclusive jets)
39 for (unsigned int i=0; i<jets.size(); ++i) {
40 _h_Inclusive_Denominator->fill(jets[i].pT()/GeV);
41
42 // Loop looking for neighboring jets for each jet
43 for (unsigned int j=0; j<jets.size(); ++j) {
44 if (i == j) continue;
45 const double dphi = deltaPhi(jets[i], jets[j]);
46
47 // Select neighboring jets and fill the numerator
48 if ( (dphi>2.*pi/3.) && (dphi<7.*pi/8.) && (jets[j].pT() > 100.*GeV) ) {
49 _h_Dphi_Numerator->fill(jets[i].pT()/GeV);
50 }
51 }
52 }
53
54 }
55
56
57 /// Calculation of Rdphi observable
58 void finalize() {
59
60 divide(_h_Dphi_Numerator, _h_Inclusive_Denominator, _h_Rdphi);
61
62 }
63
64
65 /// @name Histograms
66 Histo1DPtr _h_Dphi_Numerator, _h_Inclusive_Denominator;
67 Estimate1DPtr _h_Rdphi;
68
69
70 };
71
72
73 RIVET_DECLARE_PLUGIN(CMS_2024_I2780732);
74
75}
|