rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2015_I1310737

Jet multiplicity and differential cross-sections of $Z$+jets events in $pp$ at $\sqrt{s} = 7$ TeV
Experiment: CMS (LHC)
Inspire ID: 1310737
Status: VALIDATED
Authors:
  • Fabio Cossutti
  • Chiara La Licata
References:
  • Phys.Rev. D91 (2015) 052008
  • http://dx.doi.org/10.1103/PhysRevD.91.052008
  • arxiv:1408.3104
  • http://inspirehep.net/record/1310737
Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • Run MC generators with Z decaying leptonically into both electrons and muons at 7 TeV CoM energy. Order of 5 million unweighted events can give a reasonable global comparison, but precision in the high jet multiplicity region/high jet pt may require substantially larger samples or statistical enhancement of high jet multiplicities.

Measurements of differential cross sections are presented for the production of a Z boson and at least one hadronic jet in proton-proton collisions at $\sqrt{s}=7$ TeV, recorded by the CMS detector, using a data sample corresponding to an integrated luminosity of 4.9 $\text{fb}^{-1}$. The jet multiplicity distribution is measured for up to six jets. The differential cross sections are measured as a function of jet transverse momentum and pseudorapidity for the four highest transverse momentum jets. The distribution of the scalar sum of jet transverse momenta is also measured as a function of the jet multiplicity. The measurements are compared with theoretical predictions at leading and next-to-leading order in perturbative QCD. Cuts: First two leading electrons or muons with $p_T > 20$ GeV and $|\eta| < 2.4$ Dilepton invariant mass in the [71,111] GeV range Jets $p_T > 30$ GeV and $|\eta| < 2.4$ $\Delta{R}($lepton,jet$) > 0.5$

Source code: CMS_2015_I1310737.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FastJets.hh"
  4#include "Rivet/Projections/PromptFinalState.hh"
  5#include "Rivet/Projections/ChargedFinalState.hh"
  6#include "Rivet/Projections/InvMassFinalState.hh"
  7#include "Rivet/Projections/VisibleFinalState.hh"
  8#include "Rivet/Projections/VetoedFinalState.hh"
  9#include "Rivet/Projections/DileptonFinder.hh"
 10
 11namespace Rivet {
 12
 13
 14  class CMS_2015_I1310737 : public Analysis {
 15  public:
 16
 17    /// Constructor
 18    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2015_I1310737);
 19
 20
 21    /// Book histograms and initialise projections before the run
 22    void init() {
 23
 24      DileptonFinder zeeFinder(91.2*GeV, 0.1, Cuts::abseta < 2.4 && Cuts::pT > 20*GeV &&
 25                               Cuts::abspid == PID::ELECTRON, Cuts::massIn(71*GeV, 111*GeV));
 26      declare(zeeFinder, "ZeeFinder");
 27
 28      DileptonFinder zmumuFinder(91.2*GeV, 0.1, Cuts::abseta < 2.4 && Cuts::pT > 20*GeV &&
 29                                 Cuts::abspid == PID::MUON, Cuts::massIn(71*GeV, 111*GeV));
 30      declare(zmumuFinder, "ZmumuFinder");
 31
 32      VisibleFinalState visfs;
 33      VetoedFinalState jetConstits(visfs);
 34      jetConstits.addVetoOnThisFinalState(zeeFinder);
 35      jetConstits.addVetoOnThisFinalState(zmumuFinder);
 36
 37      FastJets akt05Jets(jetConstits, JetAlg::ANTIKT, 0.5);
 38      declare(akt05Jets, "AntiKt05Jets");
 39
 40
 41      book(_h_excmult_jets_tot ,1, 1, 1);
 42      book(_h_incmult_jets_tot ,2, 1, 1);
 43      book(_h_leading_jet_pt_tot ,3, 1, 1);
 44      book(_h_second_jet_pt_tot ,4, 1, 1);
 45      book(_h_third_jet_pt_tot ,5, 1, 1);
 46      book(_h_fourth_jet_pt_tot ,6, 1, 1);
 47      book(_h_leading_jet_eta_tot ,7, 1, 1);
 48      book(_h_second_jet_eta_tot ,8, 1, 1);
 49      book(_h_third_jet_eta_tot ,9, 1, 1);
 50      book(_h_fourth_jet_eta_tot ,10, 1, 1);
 51      book(_h_ht1_tot ,11, 1, 1);
 52      book(_h_ht2_tot ,12, 1, 1);
 53      book(_h_ht3_tot ,13, 1, 1);
 54      book(_h_ht4_tot ,14, 1, 1);
 55    }
 56
 57
 58    /// Perform the per-event analysis
 59    void analyze(const Event& event) {;
 60
 61      const DileptonFinder& zeeFS = apply<DileptonFinder>(event, "ZeeFinder");
 62      const DileptonFinder& zmumuFS = apply<DileptonFinder>(event, "ZmumuFinder");
 63
 64      const Particles& zees = zeeFS.bosons();
 65      const Particles& zmumus = zmumuFS.bosons();
 66
 67      // We did not find exactly one Z. No good.
 68      if (zees.size() + zmumus.size() != 1) {
 69        MSG_DEBUG("Did not find exactly one good Z candidate");
 70        vetoEvent;
 71      }
 72
 73      // Find the (dressed!) leptons
 74      const Particles& dressedLeptons = zees.size() ? zeeFS.constituents() : zmumuFS.constituents();
 75
 76      // Cluster jets
 77      // NB. Veto has already been applied on leptons and photons used for dressing
 78      const FastJets& fj = apply<FastJets>(event, "AntiKt05Jets");
 79      const Jets& jets = fj.jetsByPt(Cuts::abseta < 2.4 && Cuts::pT > 30*GeV);
 80
 81      // Perform lepton-jet overlap and HT calculation
 82      double ht = 0;
 83      Jets goodjets;
 84      for (const Jet& j : jets) {
 85        // Decide if this jet is "good", i.e. isolated from the leptons
 86        /// @todo Nice use-case for any() and a C++11 lambda
 87        bool overlap = false;
 88        for (const Particle& l : dressedLeptons) {
 89          if (deltaR(j, l) < 0.5) {
 90            overlap = true;
 91            break;
 92          }
 93        }
 94
 95        // Fill HT and good-jets collection
 96        if (overlap) continue;
 97        goodjets.push_back(j);
 98        ht += j.pT();
 99      }
100
101      // We don't care about events with no isolated jets
102      if (goodjets.empty()) {
103        MSG_DEBUG("No jets in event");
104        vetoEvent;
105      }
106
107
108      /////////////////
109
110
111      // Weight to be used for histo filling
112      const double w = 0.5 * 1.0;
113
114      // Fill jet number integral histograms
115      _h_excmult_jets_tot->fill(goodjets.size(), w);
116      /// @todo Could be better computed by toIntegral transform on exclusive histo
117      for (size_t iJet = 1; iJet <= goodjets.size(); iJet++ )
118        _h_incmult_jets_tot->fill(iJet, w);
119
120      // Fill leading jet histograms
121      const Jet& j1 = goodjets[0];
122      _h_leading_jet_pt_tot->fill(j1.pT()/GeV, w);
123      _h_leading_jet_eta_tot->fill(j1.abseta(), w);
124      _h_ht1_tot->fill(ht/GeV, w);
125
126      // Fill 2nd jet histograms
127      if (goodjets.size() < 2) return;
128      const Jet& j2 = goodjets[1];
129      _h_second_jet_pt_tot->fill(j2.pT()/GeV, w);
130      _h_second_jet_eta_tot->fill(j2.abseta(), w);
131      _h_ht2_tot->fill(ht/GeV, w);
132
133      // Fill 3rd jet histograms
134      if (goodjets.size() < 3) return;
135      const Jet& j3 = goodjets[2];
136      _h_third_jet_pt_tot->fill(j3.pT()/GeV, w);
137      _h_third_jet_eta_tot->fill(j3.abseta(), w);
138      _h_ht3_tot->fill(ht/GeV, w);
139
140      // Fill 4th jet histograms
141      if (goodjets.size() < 4) return;
142      const Jet& j4 = goodjets[3];
143      _h_fourth_jet_pt_tot->fill(j4.pT()/GeV, w);
144      _h_fourth_jet_eta_tot->fill(j4.abseta(), w);
145      _h_ht4_tot->fill(ht/GeV, w);
146    }
147
148
149    /// Normalise histograms etc., after the run
150    void finalize() {
151
152      const double norm = (sumOfWeights() != 0) ? crossSection()/picobarn/sumOfWeights() : 1.0;
153
154      scale(_h_excmult_jets_tot, norm );
155      scale(_h_incmult_jets_tot, norm );
156      scale(_h_leading_jet_pt_tot, norm );
157      scale(_h_second_jet_pt_tot, norm );
158      scale(_h_third_jet_pt_tot, norm );
159      scale(_h_fourth_jet_pt_tot, norm );
160      scale(_h_leading_jet_eta_tot, norm );
161      scale(_h_second_jet_eta_tot, norm );
162      scale(_h_third_jet_eta_tot, norm );
163      scale(_h_fourth_jet_eta_tot, norm );
164      scale(_h_ht1_tot, norm );
165      scale(_h_ht2_tot, norm );
166      scale(_h_ht3_tot, norm );
167      scale(_h_ht4_tot, norm );
168    }
169
170
171  private:
172
173    /// @name Histograms
174
175    Histo1DPtr _h_excmult_jets_tot,  _h_incmult_jets_tot;
176    Histo1DPtr _h_leading_jet_pt_tot, _h_second_jet_pt_tot, _h_third_jet_pt_tot, _h_fourth_jet_pt_tot;
177    Histo1DPtr _h_leading_jet_eta_tot, _h_second_jet_eta_tot, _h_third_jet_eta_tot, _h_fourth_jet_eta_tot;
178    Histo1DPtr _h_ht1_tot, _h_ht2_tot, _h_ht3_tot, _h_ht4_tot;
179
180  };
181
182
183  RIVET_DECLARE_PLUGIN(CMS_2015_I1310737);
184
185
186}