rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2017_I1610623

Measurements of differential cross sections for the associated production of a W boson and jets in proton-proton collisions at sqrt(s) = 13 TeV
Experiment: CMS (LHC)
Inspire ID: 1610623
Status: VALIDATED
Authors:
  • Emanuela Barberis
  • Apichart Hortiangtham
  • Kadir Ocalan
References:
  • DOI:10.1103/PhysRevD.96.072005
  • arXiv: 1707.05979
  • CMS-SMP-16-005
Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • 13 TeV $pp \to W+jets$.

Differential cross sections for a W boson produced in association with jets are measured in a data sample of proton-proton collisions at a center-of-mass energy of 13 TeV recorded with the CMS detector and corresponding to an integrated luminosity of 2.2 inverse femtobarns. The W bosons are identified through their decay into a muon and a neutrino. The cross sections are reported as functions of jet multiplicity, jet transverse momenta, jet rapidity $|y|$, the scalar sum of jet transverse momenta ($H_T$), and angular correlations between the muon and each jet for different jet multiplicities. The cross sections are measured in the fiducial region defined by a muon with $p_T > 25$ GeV and pseudorapidity $|\eta| < 2.4$, and by a transverse mass between the muon and the missing transverse energy $M_T > 50$ GeV. Jets are reconstructed using the anti-kT algorithm with a distance parameter R = 0.4, and only jets with $p_T > 30$ GeV, $|y| < 2.4$, and a separation of $\Delta R > 0.4$ from the muon are considered. In addition, the differential cross section is measured as a function of the angular distance (DeltaR) between the muon and the closest jet for events with one or more jets, requiring jets with $p_T > 100$ GeV and the leading jet with $p_T > 300$ GeV.

Source code: CMS_2017_I1610623.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/WFinder.hh"
  7
  8namespace Rivet {
  9
 10
 11  class CMS_2017_I1610623 : public Analysis {
 12  public:
 13
 14    /// Constructor
 15    RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2017_I1610623);
 16
 17
 18    /// @name Analysis methods
 19    //@{
 20
 21    /// Book histograms and initialise projections before the run
 22    void init() {
 23
 24      // Initialise and register projections
 25      FinalState fs;
 26      WFinder wfinder_mu(fs, Cuts::abseta < 2.4 && Cuts::pT > 0*GeV, PID::MUON, 0*GeV, 1000000*GeV, 0*GeV, 0.1, WFinder::ChargedLeptons::PROMPT, WFinder::ClusterPhotons::NODECAY, WFinder::AddPhotons::YES, WFinder::MassWindow::MT);
 27      //WFinder wfinder_mu(fs, Cuts::abseta < 2.4 && Cuts::pT > 0*GeV, PID::MUON, 0*GeV, 1000000*GeV, 0*GeV, 0.1, WFinder::ChargedLeptons::PROMPT, WFinder::ClusterPhotons::NODECAY, WFinder::AddPhotons::NO, WFinder::MassWindow::MT);
 28      declare(wfinder_mu, "WFinder_mu");
 29
 30      // Define veto FS
 31      VetoedFinalState vfs;
 32      vfs.addVetoOnThisFinalState(wfinder_mu);
 33      vfs.addVetoPairId(PID::MUON);
 34      vfs.vetoNeutrinos();
 35
 36      FastJets fastjets(vfs, FastJets::ANTIKT, 0.4);
 37      declare(fastjets, "Jets");
 38
 39      //-------------
 40      book(_hist_Mult_exc      ,"d01-x01-y01");
 41      book(_hist_inc_WJetMult  ,"d02-x01-y01");
 42
 43      //-------------
 44      book(_hist_JetPt1j ,"d03-x01-y01");
 45      book(_hist_JetPt2j ,"d04-x01-y01");
 46      book(_hist_JetPt3j ,"d05-x01-y01");
 47      book(_hist_JetPt4j ,"d06-x01-y01");
 48
 49      //-------------
 50      book(_hist_JetRap1j ,"d07-x01-y01");
 51      book(_hist_JetRap2j ,"d08-x01-y01");
 52      book(_hist_JetRap3j ,"d09-x01-y01");
 53      book(_hist_JetRap4j ,"d10-x01-y01");
 54
 55      //-------------
 56      book(_hist_Ht_1j ,"d11-x01-y01");
 57      book(_hist_Ht_2j ,"d12-x01-y01");
 58      book(_hist_Ht_3j ,"d13-x01-y01");
 59      book(_hist_Ht_4j ,"d14-x01-y01");
 60
 61      //-------------
 62      book(_hist_dphij1mu_1j , "d15-x01-y01");
 63      book(_hist_dphij2mu_2j , "d16-x01-y01");
 64      book(_hist_dphij3mu_3j , "d17-x01-y01");
 65      book(_hist_dphij4mu_4j , "d18-x01-y01");
 66
 67      //-------------
 68      book(_hist_dRmuj_1j , "d19-x01-y01");
 69
 70    }
 71
 72    // define function used for filiing inc Njets histo
 73    void _fill(Histo1DPtr& _histJetMult, std::vector<FourMomentum>& finaljet_list){
 74      _histJetMult->fill(0);
 75      for (size_t i=0 ; i<finaljet_list.size() ; ++i) {
 76        if (i==6) break;
 77        _histJetMult->fill(i+1);  // inclusive multiplicity
 78      }
 79    }
 80
 81    /// Perform the per-event analysis
 82    void analyze(const Event& event) {
 83
 84      /// @todo Do the event by event analysis here
 85      const WFinder& wfinder_mu = applyProjection<WFinder>(event, "WFinder_mu");
 86
 87      if (wfinder_mu.bosons().size() != 1) {
 88        vetoEvent;
 89      }
 90
 91      if (wfinder_mu.bosons().size() == 1) {
 92
 93        const FourMomentum lepton0 = wfinder_mu.constituentLepton().momentum();
 94        const FourMomentum neutrino = wfinder_mu.constituentNeutrino().momentum();
 95        double WmT = wfinder_mu.mT();
 96
 97        if (WmT < 50.0*GeV) vetoEvent;
 98
 99        double pt0 = lepton0.pT();
100        double eta0 = lepton0.eta();
101
102        if ( (fabs(eta0) > 2.4) || (pt0 < 25.0*GeV) ) vetoEvent;
103
104        // Obtain the jets.
105        vector<FourMomentum> finaljet_list;
106        vector<FourMomentum> jet100_list;
107        double HT = 0.0;
108
109        // loop over jets in an event, pushback in finaljet_list collection
110        for (const Jet& j : applyProjection<FastJets>(event, "Jets").jetsByPt(30.0*GeV)) {
111          const double jrap = j.momentum().rap();
112          const double jpt = j.momentum().pT();
113          if ( (fabs(jrap) < 2.4) && (deltaR(lepton0, j.momentum()) > 0.4) ) {
114            if(jpt > 30.0*GeV) {
115              finaljet_list.push_back(j.momentum());
116              HT += j.momentum().pT();
117            }
118            if(jpt > 100.0*GeV) {
119              jet100_list.push_back(j.momentum());
120            }
121          }
122        } // end looping over jets
123
124        //---------------------- FILL HISTOGRAMS ------------------
125
126        // Multiplicity exc plot.
127        _hist_Mult_exc->fill(finaljet_list.size());
128
129        // Multiplicity inc plot.
130        _fill(_hist_inc_WJetMult, finaljet_list);
131
132        // dRmuj plot.
133        double mindR(99999);
134        if(jet100_list.size()>=1) {
135          for (unsigned ji = 0; ji < jet100_list.size(); ji++){
136            double dr_(9999);
137            dr_ = fabs(deltaR(lepton0, jet100_list[ji]));
138            if (dr_ < mindR){
139              mindR = dr_;
140            }
141          }
142          if(jet100_list[0].pT() > 300.0*GeV){
143            _hist_dRmuj_1j->fill(mindR);
144          }
145        }
146
147        if(finaljet_list.size()>=1) {
148          _hist_JetPt1j->fill(finaljet_list[0].pT());
149          _hist_JetRap1j->fill(fabs(finaljet_list[0].rap()));
150          _hist_Ht_1j->fill(HT);
151          _hist_dphij1mu_1j->fill(deltaPhi(finaljet_list[0].phi(), lepton0.phi()));
152        }
153
154        if(finaljet_list.size()>=2) {
155          _hist_JetPt2j->fill(finaljet_list[1].pT());
156          _hist_JetRap2j->fill(fabs(finaljet_list[1].rap()));
157          _hist_Ht_2j->fill(HT);
158          _hist_dphij2mu_2j->fill(deltaPhi(finaljet_list[1].phi(), lepton0.phi()));
159        }
160
161        if(finaljet_list.size()>=3) {
162          _hist_JetPt3j->fill(finaljet_list[2].pT());
163          _hist_JetRap3j->fill(fabs(finaljet_list[2].rap()));
164          _hist_Ht_3j->fill(HT);
165          _hist_dphij3mu_3j->fill(deltaPhi(finaljet_list[2].phi(), lepton0.phi()));
166        }
167
168        if(finaljet_list.size()>=4) {
169          _hist_JetPt4j->fill(finaljet_list[3].pT());
170          _hist_JetRap4j->fill(fabs(finaljet_list[3].rap()));
171          _hist_Ht_4j->fill(HT);
172          _hist_dphij4mu_4j->fill(deltaPhi(finaljet_list[3].phi(), lepton0.phi()));
173        }
174      } // close the Wboson loop
175
176    } //void loop
177
178
179    /// Normalise histograms etc., after the run
180    void finalize() {
181      const double crossec = !std::isnan(crossSectionPerEvent()) ? crossSection() : 61526.7*picobarn;
182      if (std::isnan(crossSectionPerEvent())){
183        MSG_INFO("No valid cross-section given, using NNLO xsec calculated by FEWZ " << crossec/picobarn << " pb");
184      }
185
186      scale(_hist_Mult_exc, crossec/picobarn/sumOfWeights());
187      scale(_hist_inc_WJetMult, crossec/picobarn/sumOfWeights());
188
189      scale(_hist_JetPt1j, crossec/picobarn/sumOfWeights());
190      scale(_hist_JetPt2j, crossec/picobarn/sumOfWeights());
191      scale(_hist_JetPt3j, crossec/picobarn/sumOfWeights());
192      scale(_hist_JetPt4j, crossec/picobarn/sumOfWeights());
193
194      scale(_hist_JetRap1j, crossec/picobarn/sumOfWeights());
195      scale(_hist_JetRap2j, crossec/picobarn/sumOfWeights());
196      scale(_hist_JetRap3j, crossec/picobarn/sumOfWeights());
197      scale(_hist_JetRap4j, crossec/picobarn/sumOfWeights());
198
199      scale(_hist_Ht_1j, crossec/picobarn/sumOfWeights());
200      scale(_hist_Ht_2j, crossec/picobarn/sumOfWeights());
201      scale(_hist_Ht_3j, crossec/picobarn/sumOfWeights());
202      scale(_hist_Ht_4j, crossec/picobarn/sumOfWeights());
203
204      scale(_hist_dphij1mu_1j, crossec/picobarn/sumOfWeights());
205      scale(_hist_dphij2mu_2j, crossec/picobarn/sumOfWeights());
206      scale(_hist_dphij3mu_3j, crossec/picobarn/sumOfWeights());
207      scale(_hist_dphij4mu_4j, crossec/picobarn/sumOfWeights());
208
209      scale(_hist_dRmuj_1j, crossec/picobarn/sumOfWeights());
210
211    }
212
213    //@}
214
215
216  private:
217
218    /// @name Histograms
219    //@{
220    Histo1DPtr _hist_Mult_exc;
221    Histo1DPtr _hist_inc_WJetMult;
222
223    Histo1DPtr _hist_JetPt1j;
224    Histo1DPtr _hist_JetPt2j;
225    Histo1DPtr _hist_JetPt3j;
226    Histo1DPtr _hist_JetPt4j;
227
228    Histo1DPtr _hist_JetRap1j;
229    Histo1DPtr _hist_JetRap2j;
230    Histo1DPtr _hist_JetRap3j;
231    Histo1DPtr _hist_JetRap4j;
232
233    Histo1DPtr _hist_Ht_1j;
234    Histo1DPtr _hist_Ht_2j;
235    Histo1DPtr _hist_Ht_3j;
236    Histo1DPtr _hist_Ht_4j;
237
238    Histo1DPtr _hist_dphij1mu_1j;
239    Histo1DPtr _hist_dphij2mu_2j;
240    Histo1DPtr _hist_dphij3mu_3j;
241    Histo1DPtr _hist_dphij4mu_4j;
242
243    Histo1DPtr _hist_dRmuj_1j;
244    //@}
245
246  };
247
248
249  // The hook for the plugin system
250  RIVET_DECLARE_PLUGIN(CMS_2017_I1610623);
251
252
253}