rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

LHCB_2016_I1454404

Measurement of forward W and Z boson production in association with jets at LHCb
Experiment: LHCB (LHC)
Inspire ID: 1454404
Status: VALIDATED
Authors:
  • Abbie Jane Chadwick
  • Stephen Farry
References:
  • 10.1007/JHEP05(2016)131
  • arXiv: hep-ex/1605.00951
  • LHCb-PAPER-2016-011, CERN-EP-2016-092
Beams: p+ p+
Beam energies: (4000.0, 4000.0) GeV
Run details:
  • p+ p+ -> W/Z + jet + X

Measurements are made of forward W and Z production in association with jets in the forward region. Muons and jets are required to have transverse momentum in excess of 20 Gev, and to have a pseudorapidity between 2.0 and 4.5 for muons, and between 2.2 and 4.2 for jets. A single muon is required in the case of $W$ production, and two opposite sign muons with a combined invariant mass of between 60 and 120 GeV are required for $W$ production. The leptons and jets are required to be separated by a radius of 0.5 in ($\eta$, $\phi$) space.

Source code: LHCB_2016_I1454404.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/FastJets.hh"
  5#include "Rivet/Projections/WFinder.hh"
  6#include "Rivet/Projections/ZFinder.hh"
  7
  8namespace Rivet {
  9
 10
 11  /// @brief Measurement of forward W and Z boson production with jets in pp collisions at 8 TeV
 12  class LHCB_2016_I1454404 : public Analysis {
 13  public:
 14
 15    /// Constructor
 16    RIVET_DEFAULT_ANALYSIS_CTOR(LHCB_2016_I1454404);
 17
 18
 19    /// @name Analysis methods
 20    //@{
 21
 22    /// Book histograms and initialise projections before the run
 23    void init() {
 24
 25      _mode = 0;
 26      string mode = getOption("MODE");
 27      if (mode == "ALL" ) _mode = 0;
 28      else if (mode == "WpJET") _mode = 1;
 29      else if (mode == "WmJET") _mode = 2;
 30      else if (mode == "ZJET") _mode = 3;
 31      else if (mode == "WJET") _mode = 4;
 32
 33
 34      const Cut muSel = Cuts::eta >= 2.0 && Cuts::eta <= 4.5 && Cuts::pT > 20*GeV;
 35
 36      // Z boson
 37      ZFinder zfinder(FinalState(), muSel, PID::MUON, 60*GeV, 120*GeV);
 38      declare(zfinder, "ZFinder");
 39
 40      // W boson
 41      WFinder wfinder(FinalState(), muSel, PID::MUON, 0*GeV ,500*GeV ,0.);
 42      declare(wfinder, "WFinder");
 43
 44      // Jet Z
 45      FastJets jetproZ(zfinder.remainingFinalState(), FastJets::ANTIKT, 0.5);
 46      declare(jetproZ, "JetsZ");
 47
 48      // Jet W
 49      FastJets jetproW(wfinder.remainingFinalState(), FastJets::ANTIKT, 0.5);
 50      declare(jetproW, "JetsW");
 51
 52      // Book histograms
 53      /////////
 54      if (_mode == 0 || _mode == 1 || _mode == 4) {
 55        book(_h_wpj, 1, 1, 1);
 56        book(_h_eta_wpj, 4, 1, 1);
 57        book(_h_etaj_wpj, 5, 1, 1);
 58        book(_h_ptj_wpj, 6, 1, 1);
 59      }
 60      /////////
 61      if (_mode == 0 || _mode == 2 || _mode == 4) {
 62        book(_h_wmj, 1, 1, 2);
 63        book(_h_eta_wmj, 4, 1, 2);
 64        book(_h_etaj_wmj, 5, 1, 2);
 65        book(_h_ptj_wmj, 6, 1, 2);
 66      }
 67
 68      /////////
 69      if (_mode == 0 || _mode == 3) {
 70        book(_h_zj, 1, 1, 3);
 71        book(_h_yz_zj, 7, 1, 1);
 72        book(_h_etaj_zj, 8, 1, 1);
 73        book(_h_ptj_zj, 9, 1, 1);
 74        book(_h_dphi_zj, 10, 1, 1);
 75      }
 76
 77      if (_mode == 0 ){
 78        book(_h_rwz, 2,1,1);
 79        book(_h_rwpz, 2,1,2);
 80        book(_h_rwmz, 2,1,3);
 81      }
 82
 83      if (_mode == 0 || _mode == 4){
 84        book(_h_rwpm, 2,1,4);
 85        book(_h_aw, 3,1,1);
 86        // this is a temporary histogram to construct rwz later
 87        book(_h_wj, "_temp_wj", refData(1,1,1));
 88      }
 89    }
 90
 91    /// Perform the per-event analysis
 92    void analyze(const Event& event) {
 93      const Cut jetSel = Cuts::eta >= 2.2 && Cuts::eta <= 4.2 && Cuts::pT > 20*GeV;
 94      if (_mode == 0 || _mode == 3) {
 95
 96        //////////////////////////////////////////////////////////
 97        ///////////////ZFinder Muon //////////////////////////////
 98        //////////////////////////////////////////////////////////
 99
100        const ZFinder& zfinder = apply<ZFinder>(event, "ZFinder");
101        if (zfinder.bosons().size() ==1){
102          const Particles muon = zfinder.constituentLeptons(); //zfinder.constituents()?
103          const Particles Z = zfinder.bosons();
104          const FourMomentum Zmom = Z[0].momentum();
105          const Jets jetsZ = apply<FastJets>(event, "JetsZ").jetsByPt(jetSel);
106          const Jets cleanedJetsZ = filter_discard(jetsZ, [&](const Jet& j) {return any(muon, deltaRLess(j, 0.5)); });
107
108          if (cleanedJetsZ.size() > 0 && cleanedJetsZ.at(0).pT() > 20*GeV) {
109            const double yZ = Zmom.rap(); //histogram 7
110            const double etaj = cleanedJetsZ[0].eta(); //histogram 8
111            const double ptj  = cleanedJetsZ[0].pT()/GeV; //histogram 9
112            double dphi_tmp = abs(Zmom.phi() - cleanedJetsZ[0].phi());
113            const double dphi = dphi_tmp < Rivet::pi ? dphi_tmp : Rivet::twopi - dphi_tmp;
114            _h_zj->fill(sqrtS()/GeV);
115            _h_dphi_zj->fill(dphi);
116            _h_yz_zj->fill(yZ); // boson rapidity vs diff cross section
117            _h_etaj_zj->fill(etaj); // jet pseudorapidity vs diff cross section
118            _h_ptj_zj->fill(ptj); //jet transverse momentum vs diff cross section
119          }
120        }
121      }
122
123      if (_mode == 0 || _mode == 1 || _mode == 2  || _mode == 4) {
124        //////////////////////////////////////////////////////////
125        ///////////////WFinder Muon //////////////////////////////
126        //////////////////////////////////////////////////////////
127
128        const WFinder& wfinder = apply<WFinder>(event, "WFinder");
129        if (wfinder.bosons().size() == 1) {
130          const Particles Muons = wfinder.constituentLeptons();
131          const FourMomentum muonmom = Muons[0].momentum();
132          const Jets jetsW = apply<FastJets>(event, "JetsW").jetsByPt(jetSel);
133
134          const Jets cleanedJetsW = filter_discard(jetsW, [&](const Jet& j) {return any(Muons, deltaRLess(j, 0.5)); });
135
136          if (cleanedJetsW.size() > 0 && cleanedJetsW.at(0).pT() > 20*GeV) {
137            const double etaj = cleanedJetsW[0].eta(); //histogram 5
138            const double etamu = muonmom.eta(); //histogram 4
139            if( (_mode == 0 || _mode == 1 || _mode == 4) && Muons[0].charge() > 0) {
140              //fill with W related analysis
141              if (_mode != 1 ) _h_wj->fill(sqrtS()/GeV); // don't need this for single charge case
142              _h_wpj->fill(sqrtS()/GeV);
143              _h_eta_wpj->fill(etamu); // W+ Jet  muon pseudorapidity vs diff cross section
144              _h_etaj_wpj->fill(etaj); // W+ Jet jet pseudorapidity vs diff cross section
145              _h_ptj_wpj->fill(cleanedJetsW[0].pT()/GeV); // W+ Jet jet transverse momentum vs diff cross section
146            }
147            else if( (_mode == 0 || _mode == 2 || _mode == 4 ) && Muons[0].charge() < 0) {
148              //fill with W related analysis
149              if (_mode != 2) _h_wj->fill(sqrtS()/GeV);
150              _h_wmj->fill(sqrtS()/GeV); // don't need this for single charge case
151              _h_eta_wmj->fill(etamu); // W- Jet  muon pseudorapidity vs diff cross section
152              _h_etaj_wmj->fill(etaj); // W- Jet jet pseudorapidity vs diff cross section
153              _h_ptj_wmj->fill(cleanedJetsW[0].pT()/GeV); // W+ Jet jet transverse momentum vs diff cross section
154            }
155          }
156        }
157
158      }
159    }
160
161    /// Normalise histograms etc., after the run
162    void finalize() {
163
164      double scalefactor = crossSection()/picobarn/sumOfWeights();
165
166      if(_mode == 0 || _mode == 1 || _mode == 4) {
167        scale({_h_wpj, _h_eta_wpj, _h_etaj_wpj, _h_ptj_wpj}, scalefactor);
168      }
169      if(_mode == 0 || _mode == 2 || _mode == 4) {
170        scale({_h_wmj, _h_eta_wmj, _h_etaj_wmj, _h_ptj_wmj}, scalefactor);
171      }
172      if(_mode == 0 || _mode == 3) {
173        scale({_h_zj, _h_yz_zj, _h_etaj_zj, _h_ptj_zj, _h_dphi_zj}, scalefactor);
174      }
175      if (_mode == 0 ) {
176        scale(_h_wj, scalefactor); // need to scale this for consistency
177        divide(_h_wpj, _h_zj, _h_rwpz);
178        divide(_h_wmj, _h_zj, _h_rwmz);
179        divide(_h_wj, _h_zj, _h_rwz);
180
181      }
182      if (_mode == 0 || _mode == 4) {
183        divide(_h_wpj, _h_wmj, _h_rwpm);
184        asymm(_h_wpj, _h_wmj, _h_aw);
185      }
186
187    }
188
189  protected:
190
191    Log& log = getLog(); // in Analysis or Projection
192    size_t _mode;
193
194    /// histograms
195    Histo1DPtr _h_wpj, _h_wmj, _h_wj, _h_zj;
196    Scatter2DPtr _h_rwz, _h_rwpz, _h_rwmz, _h_rwpm, _h_aw;
197    Histo1DPtr _h_eta_wpj, _h_eta_wmj, _h_etaj_wpj, _h_etaj_wmj, _h_ptj_wpj, _h_ptj_wmj;
198    Histo1DPtr _h_yz_zj, _h_etaj_zj, _h_ptj_zj, _h_dphi_zj;
199  };
200
201
202  RIVET_DECLARE_PLUGIN(LHCB_2016_I1454404);
203
204}