rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2014_I1303894

Differential cross-section of $W$ bosons + jets in $pp$ collisions at $\sqrt{s}=7$ TeV
Experiment: CMS (LHC)
Inspire ID: 1303894
Status: VALIDATED
Authors:
  • Darin Baumgartel
  • Emanuela Barberis
References:
  • Phys. Lett. B741 (2014) 12-37
  • https://inspirehep.net/record/1303894
  • http://arxiv.org/abs/1406.7533
Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • Run MC generators with $W$ decaying leptonically at 7 TeV CoM energy. A large number of events are required to populate the high jet multiplicity region. Suitable results can be achieved with 85M events.

A study of jet production in association with $W$ bosons has been performed, in events with the $W$ decaying to a muon. Jets are required to have $pT > 30$ GeV and $|eta| < 2.4$. Muons are required to have $pT > 25$ and $|eta| < 2.1$. Jets are only considered if they are separated from the muon by $\Delta{R} > 0.5$. Muons are dressed with photons in a cone of $0.1$ around the muon.

Source code: CMS_2014_I1303894.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5#include "Rivet/Projections/VetoedFinalState.hh"
  6#include "Rivet/Projections/PromptFinalState.hh"
  7#include "Rivet/Projections/InvMassFinalState.hh"
  8#include "Rivet/Projections/MissingMomentum.hh"
  9#include "Rivet/Projections/LeptonFinder.hh"
 10
 11namespace Rivet {
 12
 13
 14  /// @brief Differential cross-section of W bosons + jets in pp collisions at sqrt(s)=7 TeV
 15  ///
 16  /// @author Darin Baumgartel (darinb@cern.ch)
 17  ///
 18  /// Based on Rivet analysis originally created by Anil Singh (anil@cern.ch), Lovedeep Saini (lovedeep@cern.ch)
 19  class CMS_2014_I1303894 : public Analysis {
 20  public:
 21
 22    /// Constructor
 23    CMS_2014_I1303894()
 24      : Analysis("CMS_2014_I1303894")
 25    {   }
 26
 27
 28    // Book histograms and initialise projections before the run
 29    void init() {
 30      // Prompt leptons only, no test on nu flavour.
 31      // Projections
 32      const FinalState fs;
 33      declare(fs, "FS");
 34
 35      MissingMomentum missing(fs);
 36      declare(missing, "MET");
 37
 38      PromptFinalState pfs(fs);
 39      IdentifiedFinalState bareMuons(pfs);
 40      bareMuons.acceptIdPair(PID::MUON);
 41      LeptonFinder muonClusters(bareMuons, fs, -1); //, Cuts::open(), false, false);
 42      declare(muonClusters, "muonClusters");
 43
 44      VetoedFinalState jetFS(fs);
 45      jetFS.addVetoOnThisFinalState(muonClusters);
 46      jetFS.vetoNeutrinos();
 47      FastJets jetprojection(jetFS, JetAlg::ANTIKT, 0.5);
 48      declare(jetprojection, "Jets");
 49
 50      // Histograms
 51      book(_histDPhiMuJet1 ,1,1,1);
 52      book(_histDPhiMuJet2 ,2,1,1);
 53      book(_histDPhiMuJet3 ,3,1,1);
 54      book(_histDPhiMuJet4 ,4,1,1);
 55
 56      book(_histEtaJet1 ,5,1,1);
 57      book(_histEtaJet2 ,6,1,1);
 58      book(_histEtaJet3 ,7,1,1);
 59      book(_histEtaJet4 ,8,1,1);
 60
 61      book(_histHT1JetInc ,9,1,1);
 62      book(_histHT2JetInc ,10,1,1);
 63      book(_histHT3JetInc ,11,1,1);
 64      book(_histHT4JetInc ,12,1,1);
 65
 66      book(_histJet30MultExc  ,13,1,1);
 67      book(_histJet30MultInc  ,14,1,1);
 68
 69      book(_histPtJet1 ,15,1,1);
 70      book(_histPtJet2 ,16,1,1);
 71      book(_histPtJet3 ,17,1,1);
 72      book(_histPtJet4 ,18,1,1);
 73
 74      // Counters
 75      book(_n_1jet, "n_1jet");
 76      book(_n_2jet, "n_2jet");
 77      book(_n_3jet, "n_3jet");
 78      book(_n_4jet, "n_4jet");
 79      book(_n_inclusivebinsummation, "n_inclusivebinsummation");
 80    }
 81
 82
 83    void analyze(const Event& event) {
 84      // Get the dressed muon
 85      const LeptonFinder& muonClusters = apply<LeptonFinder>(event, "muonClusters");
 86      int nmu = muonClusters.dressedLeptons().size();
 87      if (nmu < 1) vetoEvent;
 88      DressedLepton dressedmuon = muonClusters.dressedLeptons()[0];
 89      if (dressedmuon.momentum().abseta() > 2.1) vetoEvent;
 90      if (dressedmuon.momentum().pT() < 25.0*GeV) vetoEvent;
 91
 92      // Check that the muon and neutrino are not decay products of tau
 93      if (dressedmuon.bareLepton().hasAncestorWith(Cuts::pid ==  PID::TAU)) vetoEvent;
 94      if (dressedmuon.bareLepton().hasAncestorWith(Cuts::pid == -PID::TAU)) vetoEvent;
 95
 96      // Get the missing momentum
 97      const MissingMomentum& met = apply<MissingMomentum>(event, "MET");
 98      const double ptmet = met.visibleMomentum().pT();
 99      const double phimet = (-met.visibleMomentum()).phi();
100
101      // Calculate MET and MT(mu,MET), and remove events with MT < 50 GeV
102      const double ptmuon = dressedmuon.pT();
103      const double phimuon = dressedmuon.phi();
104      const double mt_mumet = sqrt(2*ptmuon*ptmet*(1.0 - cos(phimet-phimuon)));
105
106      // Remove events in MT < 50 region
107      if (mt_mumet < 50*GeV) vetoEvent;
108
109      // Loop over jets and fill pt/eta/phi quantities in vectors
110      const Jets& jets_filtered = apply<FastJets>(event, "Jets").jetsByPt();
111      vector<float> finaljet_pT_list, finaljet_eta_list, finaljet_phi_list;
112      double htjets = 0.0;
113      for (size_t ii = 0; ii < jets_filtered.size(); ++ii) {
114        // Jet pT/eta/phi
115        double jet_pt = jets_filtered[ii].pT();
116        double jet_eta = jets_filtered[ii].eta();
117        double jet_phi = jets_filtered[ii].phi();
118
119        // Kinemetic cuts for jet acceptance
120        if (fabs(jet_eta) > 2.4) continue;
121        if (jet_pt < 30.0*GeV) continue;
122        if (deltaR(dressedmuon, jets_filtered[ii]) < 0.5) continue;
123
124        // Add jet to jet list and increases the HT variable
125        finaljet_pT_list.push_back(jet_pt);
126        finaljet_eta_list.push_back(jet_eta);
127        finaljet_phi_list.push_back(jet_phi);
128        htjets += fabs(jet_pt);
129      }
130
131
132      // Filling of histograms:
133      // Fill as many jets as there are into the exclusive jet multiplicity
134      if (!finaljet_pT_list.empty())
135        _histJet30MultExc->fill(finaljet_pT_list.size());
136
137      for (size_t ij = 0; ij < finaljet_pT_list.size(); ++ij) {
138        _histJet30MultInc->fill(ij+1);
139        _n_inclusivebinsummation->fill();
140      }
141
142      if (finaljet_pT_list.size() >= 1) {
143        _histPtJet1->fill(finaljet_pT_list[0]);
144        _histEtaJet1->fill(fabs(finaljet_eta_list[0]));
145        _histDPhiMuJet1->fill(deltaPhi(finaljet_phi_list[0], phimuon));
146        _histHT1JetInc->fill(htjets);
147        _n_1jet->fill();
148      }
149
150      if (finaljet_pT_list.size() >= 2) {
151        _histPtJet2->fill(finaljet_pT_list[1]);
152        _histEtaJet2->fill(fabs(finaljet_eta_list[1]));
153        _histDPhiMuJet2->fill(deltaPhi(finaljet_phi_list[1], phimuon));
154        _histHT2JetInc->fill(htjets);
155        _n_2jet->fill();
156      }
157
158      if (finaljet_pT_list.size() >= 3) {
159        _histPtJet3->fill(finaljet_pT_list[2]);
160        _histEtaJet3->fill(fabs(finaljet_eta_list[2]));
161        _histDPhiMuJet3->fill(deltaPhi(finaljet_phi_list[2], phimuon));
162        _histHT3JetInc->fill(htjets);
163        _n_3jet->fill();
164      }
165
166      if (finaljet_pT_list.size() >=4 ) {
167        _histPtJet4->fill(finaljet_pT_list[3]);
168        _histEtaJet4->fill(fabs(finaljet_eta_list[3]));
169        _histDPhiMuJet4->fill(deltaPhi(finaljet_phi_list[3], phimuon));
170        _histHT4JetInc-> fill(htjets);
171        _n_4jet->fill();
172      }
173
174    }
175
176
177    // Finalize the histograms.
178    void finalize() {
179
180      const double inclusive_cross_section = crossSection()/picobarn;
181      const double norm_1jet_histo = inclusive_cross_section*dbl(*_n_1jet)/sumOfWeights();
182      const double norm_2jet_histo = inclusive_cross_section*dbl(*_n_2jet)/sumOfWeights();
183      const double norm_3jet_histo = inclusive_cross_section*dbl(*_n_3jet)/sumOfWeights();
184      const double norm_4jet_histo = inclusive_cross_section*dbl(*_n_4jet)/sumOfWeights();
185      const double norm_incmultiplicity = inclusive_cross_section*dbl(*_n_inclusivebinsummation)/sumOfWeights();
186
187      normalize(_histJet30MultExc, norm_1jet_histo);
188      normalize(_histJet30MultInc, norm_incmultiplicity);
189
190      normalize(_histPtJet1, norm_1jet_histo);
191      normalize(_histHT1JetInc, norm_1jet_histo);
192      normalize(_histEtaJet1, norm_1jet_histo);
193      normalize(_histDPhiMuJet1, norm_1jet_histo);
194
195      normalize(_histPtJet2, norm_2jet_histo);
196      normalize(_histHT2JetInc, norm_2jet_histo);
197      normalize(_histEtaJet2, norm_2jet_histo);
198      normalize(_histDPhiMuJet2, norm_2jet_histo);
199
200      normalize(_histPtJet3, norm_3jet_histo);
201      normalize(_histHT3JetInc, norm_3jet_histo);
202      normalize(_histEtaJet3, norm_3jet_histo);
203      normalize(_histDPhiMuJet3, norm_3jet_histo);
204
205      normalize(_histPtJet4, norm_4jet_histo);
206      normalize(_histHT4JetInc, norm_4jet_histo);
207      normalize(_histEtaJet4, norm_4jet_histo);
208      normalize(_histDPhiMuJet4, norm_4jet_histo);
209    }
210
211
212  private:
213
214    Histo1DPtr  _histJet30MultExc, _histJet30MultInc;
215    Histo1DPtr  _histPtJet1, _histPtJet2, _histPtJet3, _histPtJet4;
216    Histo1DPtr  _histEtaJet1, _histEtaJet2, _histEtaJet3, _histEtaJet4;
217    Histo1DPtr  _histDPhiMuJet1, _histDPhiMuJet2, _histDPhiMuJet3, _histDPhiMuJet4;
218    Histo1DPtr  _histHT1JetInc, _histHT2JetInc, _histHT3JetInc, _histHT4JetInc;
219
220    CounterPtr _n_1jet, _n_2jet, _n_3jet, _n_4jet, _n_inclusivebinsummation;
221
222  };
223
224
225  RIVET_DECLARE_PLUGIN(CMS_2014_I1303894);
226
227}