rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2021_I1866118

Study of Z boson plus jets events using variables sensitive to double-parton scattering in pp collisions at 13 TeV
Experiment: CMS (LHC)
Inspire ID: 1866118
Status: VALIDATED
Authors:
  • cms-pag-conveners-smp@cern.ch
  • Rajat Gupta
  • Meena
  • Sunil Bansal
  • J.B. Singh
References:
  • CMS-SMP-20-009
  • arXiv: 2105.14511
  • JHEP 10 (2021) 176
Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • Drell-Yan events with $Z/\gamma^* \to \mu\mu$.

Double parton scattering is investigated using events with a Z boson and jets. The measurements are performed with proton-proton collision data recorded by the CMS experiment at the LHC at $\sqrt{s} = 13$ TeV, corresponding to an integrated luminosity of $35.8 \mathrm{fb}^{-1}$ collected in the year 2016. Differential cross sections of Z$+ \geq$ 1 jet and Z$+ \geq$ 2 jets are measured with transverse momentum of the jets above 20 GeV and pseudorapidity $|\eta| <$ 2.4. Several distributions with sensitivity to double parton scattering effects are measured as functions of the angle and the transverse momentum imbalance between the Z boson and the jets. The measured distributions are compared with predictions from several event generators with different hadronization models and different parameter settings for multiparton interactions. The measured distributions show a dependence on the hadronization and multiparton interaction simulation parameters, and are important input for future improvements of the simulations.

Source code: CMS_2021_I1866118.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/ZFinder.hh"
  5#include "Rivet/Projections/FinalState.hh"
  6#include "Rivet/Projections/ChargedFinalState.hh"
  7#include "Rivet/Projections/VetoedFinalState.hh"
  8#include "Rivet/Projections/FastJets.hh"
  9
 10namespace Rivet {
 11
 12  /// @brief Study of Z boson plus jets events using variables sensitive to double-parton scattering in pp collisions at 13 TeV
 13  class CMS_2021_I1866118 : public Analysis {
 14  public:
 15    /// Constructor
 16    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2021_I1866118);
 17
 18    /// Initialization
 19    void init() {
 20      const FinalState fs;
 21
 22      Cut cut = Cuts::abseta < 2.4 && Cuts::pT > 27 * GeV;
 23      ZFinder zmumufinder(fs, cut, PID::MUON, 70 * GeV, 110 * GeV);
 24      declare(zmumufinder, "zmumufinder");
 25
 26      // Define veto FS in order to prevent Z-decay products entering the jet algorithm
 27
 28      VetoedFinalState had_fs;
 29      had_fs.addVetoOnThisFinalState(zmumufinder);
 30      FastJets jets(had_fs, FastJets::ANTIKT, 0.4);
 31      jets.useInvisibles();
 32      declare(jets, "jets");
 33
 34      book(h_dphi_Z1J_cn, 1, 1, 1);
 35      book(h_reldpt_Z1J_cn, 2, 1, 1);
 36      book(h_dphi_Zdijet_Z2J_cn, 3, 1, 1);
 37      book(h_reldpt_Zdijet_Z2J_cn, 4, 1, 1);
 38      book(h_reldpt_j1j2_Z2J_cn, 5, 1, 1);
 39
 40      book(h_dphi_Z1J_sc, 6, 1, 1);
 41      book(h_reldpt_Z1J_sc, 7, 1, 1);
 42      book(h_dphi_Zdijet_Z2J_sc, 8, 1, 1);
 43      book(h_reldpt_Zdijet_Z2J_sc, 9, 1, 1);
 44      book(h_reldpt_j1j2_Z2J_sc, 10, 1, 1);
 45    }
 46
 47    /// Perform the per-event analysis
 48    void analyze(const Event& event) {
 49      const ZFinder& zmumufinder = apply<ZFinder>(event, "zmumufinder");
 50      const Particles& zmumus = zmumufinder.bosons();
 51      if (zmumus.size() != 1) {
 52        vetoEvent;
 53      }
 54
 55      // Find the (dressed!) leptons
 56      const Particles& leptons = zmumufinder.constituents();
 57      if (leptons.size() != 2)
 58        vetoEvent;
 59
 60      Jets jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 20 * GeV && Cuts::abseta < 2.4);
 61      idiscardIfAnyDeltaRLess(jets, leptons, 0.4);
 62
 63      const size_t Njets = jets.size();
 64
 65      if (Njets < 1)
 66        vetoEvent;
 67
 68      if (Njets >= 1) {
 69        h_dphi_Z1J_cn->fill(deltaPhi(zmumus[0], jets[0]));
 70        h_reldpt_Z1J_cn->fill((zmumus[0] + jets[0].momentum()).pt() / (zmumus[0].pt() + jets[0].pt()));
 71
 72        h_dphi_Z1J_sc->fill(deltaPhi(zmumus[0], jets[0]));
 73        h_reldpt_Z1J_sc->fill((zmumus[0] + jets[0].momentum()).pt() / (zmumus[0].pt() + jets[0].pt()));
 74      }
 75
 76      if (Njets >= 2) {
 77        FourMomentum dij = jets[0].momentum() + jets[1].momentum();
 78
 79        h_dphi_Zdijet_Z2J_cn->fill(deltaPhi(zmumus[0], dij));
 80        h_reldpt_Zdijet_Z2J_cn->fill((zmumus[0] + dij).pt() / (zmumus[0].pt() + dij.pt()));
 81        h_reldpt_j1j2_Z2J_cn->fill(dij.pt() / (jets[0].pt() + jets[1].pt()));
 82
 83        h_dphi_Zdijet_Z2J_sc->fill(deltaPhi(zmumus[0], dij));
 84        h_reldpt_Zdijet_Z2J_sc->fill((zmumus[0] + dij).pt() / (zmumus[0].pt() + dij.pt()));
 85        h_reldpt_j1j2_Z2J_sc->fill(dij.pt() / (jets[0].pt() + jets[1].pt()));
 86      }
 87    }
 88
 89    void normalizeToSum(Histo1DPtr hist) {
 90      double sum = 0.;
 91      for (size_t i = 0; i < hist->numBins(); ++i) {
 92        sum += hist->bin(i).height();
 93      }
 94      scale(hist, 1. / sum);
 95    }
 96
 97    /// Normalise histograms etc., after the run
 98    void finalize() {
 99      double norm = (sumOfWeights() != 0) ? crossSection()/picobarn/sumOfWeights() : 1.0;
100
101      scale(h_dphi_Z1J_cn, norm);
102      scale(h_reldpt_Z1J_cn, norm);
103      scale(h_dphi_Zdijet_Z2J_cn, norm);
104      scale(h_reldpt_Zdijet_Z2J_cn, norm);
105      scale(h_reldpt_j1j2_Z2J_cn, norm);
106
107      normalizeToSum(h_dphi_Z1J_sc);
108      normalizeToSum(h_reldpt_Z1J_sc);
109      normalizeToSum(h_dphi_Zdijet_Z2J_sc);
110      normalizeToSum(h_reldpt_Zdijet_Z2J_sc);
111      normalizeToSum(h_reldpt_j1j2_Z2J_sc);
112    }
113
114  private:
115    /// @name Histogram objects
116    //@{
117    Histo1DPtr h_dphi_Z1J_cn;
118    Histo1DPtr h_reldpt_Z1J_cn;
119    Histo1DPtr h_dphi_Zdijet_Z2J_cn;
120    Histo1DPtr h_reldpt_Zdijet_Z2J_cn;
121    Histo1DPtr h_reldpt_j1j2_Z2J_cn;
122
123    Histo1DPtr h_dphi_Z1J_sc;
124    Histo1DPtr h_reldpt_Z1J_sc;
125    Histo1DPtr h_dphi_Zdijet_Z2J_sc;
126    Histo1DPtr h_reldpt_Zdijet_Z2J_sc;
127    Histo1DPtr h_reldpt_j1j2_Z2J_sc;
128
129    //@}
130  };
131
132  // Hook for the plugin system
133  RIVET_DECLARE_PLUGIN(CMS_2021_I1866118);
134
135}  // namespace Rivet