Rivet analyses referenceLHCB_2016_I1454404Measurement of forward W and Z boson production in association with jets at LHCbExperiment: LHCB (LHC) Inspire ID: 1454404 Status: VALIDATED Authors:
Beam energies: (4000.0, 4000.0) GeV Run details:
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}
|