Rivet analyses referenceATLAS_2018_I1635273W + jets production at 8 TeVExperiment: ATLAS (LHC) Inspire ID: 1635273 Status: VALIDATED Authors:
Beam energies: (4000.0, 4000.0) GeV Run details:
This paper presents a measurement of the $W$ boson production cross section and the $W^{+}/W^{-}$ cross-section ratio, both in association with jets, in proton-proton collisions at $\sqrt{s}=8$ TeV with the ATLAS experiment at the Large Hadron Collider. The measurement is performed in final states containing one electron and missing transverse momentum using data corresponding to an integrated luminosity of 20.2 fb$^{-1}$. Differential cross sections for events with at least one or two jets are presented for a range of observables, including jet transverse momenta and rapidities, the scalar sum of transverse momenta of the visible particles and the missing transverse momentum in the event, and the transverse momentum of the $W$ boson. For a subset of the observables, the differential cross sections of positively and negatively charged $W$ bosons are measured separately. In the cross-section ratio of $W^{+}/W^{-}$ the dominant systematic uncertainties cancel out, improving the measurement precision by up to a factor of nine. The observables and ratios selected for this paper provide valuable input for the up quark, down quark, and gluon parton distribution functions of the proton. Source code: ATLAS_2018_I1635273.cc 1#include "Rivet/Analysis.hh"
2#include "Rivet/Projections/FinalState.hh"
3#include "Rivet/Projections/InvisibleFinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Projections/PromptFinalState.hh"
6#include "Rivet/Projections/DressedLeptons.hh"
7
8namespace Rivet {
9
10 /// @brief W + jets production at 8 TeV
11 class ATLAS_2018_I1635273 : public Analysis {
12 public:
13
14 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2018_I1635273);
15
16 // Book histograms and initialise projections before the run
17 void init() {
18
19 // Get options from the new option system
20 // Default uses electrons
21 // EL looks for W->ev candidate
22 // MU looks for W->mv candidate
23 _mode = 0;
24 if ( getOption("LMODE") == "EL" ) _mode = 0;
25 if ( getOption("LMODE") == "MU" ) _mode = 1;
26
27 FinalState fs;
28 Cut cuts = Cuts::pT > 25*GeV && Cuts::abseta < 2.5;
29
30 // Get photons to dress leptons
31 // (Paper says "radiated photons", but there was
32 // no promptness requirement in the analysis code)
33 FinalState photons(Cuts::abspid == PID::PHOTON);
34
35 // Get dressed leptons
36 PromptFinalState leptons(Cuts::abspid == (_mode? PID::MUON : PID::ELECTRON), false);
37 DressedLeptons dressedleptons(photons, leptons, 0.1, cuts, true);
38 declare(dressedleptons, "DressedLeptons");
39
40 // Get neutrinos for MET calculation
41 declare(InvisibleFinalState(true), "InvFS"); // true = only allow prompt invisibles
42
43 // jets
44 FastJets jets(fs, FastJets::ANTIKT, 0.4, JetAlg::Muons::NONE);
45 declare(jets, "Jets");
46
47 // book histograms
48 book(_h["N_incl_pb"], 1, 1, 1);
49 book(_h["HT_1j_fb"], 6, 1, 1);
50 book(_h["W_pt_1j_fb"], 11, 1, 1);
51 book(_h["jet_pt1_1j_fb"], 16, 1, 1);
52 book(_h["jet_y1_1j_fb"], 21, 1, 1);
53 book(_h["jet_pt2_2j_fb"], 26, 1, 1);
54 book(_h["jet_y2_2j_fb"], 28, 1, 1);
55 book(_h["DeltaRj12_2j_fb"], 30, 1, 1);
56 book(_h["jet_mass12_2j_fb"], 32, 1, 1);
57 book(_h["N_pb"], 34, 1, 1);
58 book(_h["HT_2j_fb"], 36, 1, 1);
59 book(_h["W_pt_2j_fb"], 41, 1, 1);
60 book(_h["jet_pt1_2j_fb"], 46, 1, 1);
61 book(_h["el_eta_0j_pb"], 51, 1, 1);
62 book(_h["el_eta_1j_pb"], 56, 1, 1);
63
64 book(_h["Wplus_N_incl_pb"], 3, 1, 1);
65 book(_h["Wplus_HT_1j_fb"], 8, 1, 1);
66 book(_h["Wplus_W_pt_1j_fb"], 13, 1, 1);
67 book(_h["Wplus_jet_pt1_1j_fb"], 18, 1, 1);
68 book(_h["Wplus_jet_y1_1j_fb"], 23, 1, 1);
69 book(_h["Wplus_HT_2j_fb"], 38, 1, 1);
70 book(_h["Wplus_W_pt_2j_fb"], 43, 1, 1);
71 book(_h["Wplus_jet_pt1_2j_fb"], 48, 1, 1);
72 book(_h["Wplus_el_eta_0j_pb"], 53, 1, 1);
73 book(_h["Wplus_el_eta_1j_pb"], 58, 1, 1);
74
75 book(_h["Wminus_N_incl_pb"], 3, 1, 2);
76 book(_h["Wminus_HT_1j_fb"], 8, 1, 2);
77 book(_h["Wminus_W_pt_1j_fb"], 13, 1, 2);
78 book(_h["Wminus_jet_pt1_1j_fb"],18, 1, 2);
79 book(_h["Wminus_jet_y1_1j_fb"], 23, 1, 2);
80 book(_h["Wminus_HT_2j_fb"], 38, 1, 2);
81 book(_h["Wminus_W_pt_2j_fb"], 43, 1, 2);
82 book(_h["Wminus_jet_pt1_2j_fb"],48, 1, 2);
83 book(_h["Wminus_el_eta_0j_pb"], 53, 1, 2);
84 book(_h["Wminus_el_eta_1j_pb"], 58, 1, 2);
85
86 // and ratios
87 book(_r["WplusOverWminus_N_incl_pb"], 3, 1, 3);
88 book(_r["WplusOverWminus_HT_1j_fb"], 8, 1, 3);
89 book(_r["WplusOverWminus_W_pt_1j_fb"], 13, 1, 3);
90 book(_r["WplusOverWminus_jet_pt1_1j_fb"], 18, 1, 3);
91 book(_r["WplusOverWminus_jet_y1_1j_fb"], 23, 1, 3);
92 book(_r["WplusOverWminus_HT_2j_fb"], 38, 1, 3);
93 book(_r["WplusOverWminus_W_pt_2j_fb"], 43, 1, 3);
94 book(_r["WplusOverWminus_jet_pt1_2j_fb"], 48, 1, 3);
95 book(_r["WplusOverWminus_el_eta_0j_pb"], 53, 1, 3);
96 book(_r["WplusOverWminus_el_eta_1j_pb"], 58, 1, 3);
97
98 }
99
100
101 // Perform the per-event analysis
102 void analyze(const Event& event) {
103
104 // retrieve the dressed electrons
105 const Particles& signal_leptons = apply<DressedLeptons>(event, "DressedLeptons").particlesByPt();
106 if (signal_leptons.size() != 1 ) vetoEvent;
107 const Particle& lepton = signal_leptons[0];
108
109 // calulate MET and mT
110 const Particles& invisibles = apply<InvisibleFinalState>(event, "InvFS").particles();
111 FourMomentum pMET = sum(invisibles, Kin::mom, FourMomentum()).setZ(0);
112 const double MET = pMET.pT() / GeV;
113 const double mT = sqrt( 2 * lepton.Et() / GeV * MET * (1 - cos(deltaPhi(lepton,pMET))));
114 if ( MET <= 25. ) vetoEvent;
115 if ( mT <= 40. ) vetoEvent;
116
117 // retrieve jets
118 Jets jets = apply<FastJets>(event, "Jets").jetsByPt(Cuts::pT > 30*GeV && Cuts::absrap < 4.4);
119
120 idiscardIfAnyDeltaRLess(jets, signal_leptons, 0.2);
121
122 // apply event selection on dR
123 for (const Jet& j : jets) {
124 if (deltaR(j, lepton) < 0.4) vetoEvent;
125 }
126
127 // calculate the observables
128 const double w_pt = (lepton.momentum() + pMET).pT() / GeV;
129 const size_t njets = jets.size();
130 double ST = sum(jets, Kin::pT, 0.0); // scalar pT sum of all selected jets
131 const double HT = ST + lepton.pT() / GeV + MET; //missET;
132
133 // fill W histograms
134 _h["N_pb"]->fill(njets);
135 for (size_t i = 0; i <= njets; ++i) {
136 _h["N_incl_pb"]->fill(i);
137 }
138 _h["el_eta_0j_pb"]->fill(lepton.abseta());
139
140 if (njets > 0) {
141 _h["HT_1j_fb"]->fill(HT);
142 _h["W_pt_1j_fb"]->fill(w_pt);
143 _h["jet_pt1_1j_fb"]->fill(jets[0].pT()/GeV);
144 _h["jet_y1_1j_fb"]->fill(jets[0].absrap());
145 _h["el_eta_1j_pb"]->fill(lepton.abseta());
146 }
147 if (njets > 1) {
148 _h["HT_2j_fb"]->fill(HT);
149 _h["W_pt_2j_fb"]->fill(w_pt);
150 _h["jet_pt1_2j_fb"]->fill(jets[0].pT()/GeV);
151 _h["DeltaRj12_2j_fb"]->fill(deltaR(jets[0],jets[1]));
152 _h["jet_pt2_2j_fb"]->fill(jets[1].pT()/GeV);
153 _h["jet_y2_2j_fb"]->fill(jets[1].absrap());
154 _h["jet_mass12_2j_fb"]->fill( (jets[0].mom()+jets[1].mom()).mass()/GeV);
155 }
156 // fill W+ histograms
157 if (lepton.charge() > 0) {
158 for (size_t i = 0; i <= njets; ++i) {
159 _h["Wplus_N_incl_pb"]->fill(i);
160 }
161 _h["Wplus_el_eta_0j_pb"]->fill(lepton.abseta());
162 if (njets > 0) {
163 _h["Wplus_HT_1j_fb"]->fill(HT);
164 _h["Wplus_W_pt_1j_fb"]->fill(w_pt);
165 _h["Wplus_jet_pt1_1j_fb"]->fill(jets[0].pT()/GeV);
166 _h["Wplus_jet_y1_1j_fb"]->fill(jets[0].absrap());
167 _h["Wplus_el_eta_1j_pb"]->fill(lepton.abseta());
168 if (njets > 1) {
169 _h["Wplus_HT_2j_fb"]->fill(HT);
170 _h["Wplus_W_pt_2j_fb"]->fill(w_pt);
171 _h["Wplus_jet_pt1_2j_fb"]->fill(jets[0].pT()/GeV);
172 }
173 }
174 }
175 // fill W- histograms
176 if (lepton.charge() < 0) {
177 for (size_t i = 0; i <= njets; ++i) {
178 _h["Wminus_N_incl_pb"]->fill(i);
179 }
180 _h["Wminus_el_eta_0j_pb"]->fill(lepton.abseta());
181 if (njets > 0) {
182 _h["Wminus_HT_1j_fb"]->fill(HT);
183 _h["Wminus_W_pt_1j_fb"]->fill(w_pt);
184 _h["Wminus_jet_pt1_1j_fb"]->fill(jets[0].pT()/GeV);
185 _h["Wminus_jet_y1_1j_fb"]->fill(jets[0].absrap());
186 _h["Wminus_el_eta_1j_pb"]->fill(lepton.abseta());
187 if (njets > 1) {
188 _h["Wminus_HT_2j_fb"]->fill(HT);
189 _h["Wminus_W_pt_2j_fb"]->fill(w_pt);
190 _h["Wminus_jet_pt1_2j_fb"]->fill(jets[0].pT()/GeV);
191 }
192 }
193 }
194 }
195
196
197
198 void finalize() {
199 const double scalefactor_fb = crossSection() / sumOfWeights() / femtobarn;
200 const double scalefactor_pb = crossSection() / sumOfWeights() / picobarn;
201
202 for (auto& hit : _h){
203 if (hit.first.find("_fb") != string::npos) scale(hit.second, scalefactor_fb);
204 else scale(hit.second, scalefactor_pb);
205 }
206 for (auto& rit : _r) {
207 string ratio_label = "WplusOverWminus";
208 string num_name = rit.first;
209 string denom_name = rit.first;
210 num_name.replace(rit.first.find(ratio_label),ratio_label.length(),"Wplus");
211 denom_name.replace(rit.first.find(ratio_label),ratio_label.length(),"Wminus");
212 divide(_h[num_name], _h[denom_name], rit.second);
213 }
214 }
215
216 protected:
217 // Data members like post-cuts event weight counters go here
218 size_t _mode;
219
220 private:
221 map<string, Histo1DPtr> _h;
222 map<string, Scatter2DPtr> _r;
223 };
224
225 // The hook for the plugin system
226 RIVET_DECLARE_PLUGIN(ATLAS_2018_I1635273);
227}
|