Rivet analyses referenceATLAS_2017_I1514251Z plus jets at 13 TeVExperiment: ATLAS (LHC) Inspire ID: 1514251 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
Measurements of the production cross section of a $Z$ boson in association with jets in proton-proton collisions at $\sqrt{s}=13$ TeV are presented, using data corresponding to an integrated luminosity of 3.16 fb${}^{-1}$ collected by the ATLAS experiment at the CERN Large Hadron Collider in 2015. Inclusive and differential cross sections are measured for events containing a $Z$ boson decaying to electrons or muons and produced in association with up to seven jets with $p_\text{T}>30$ GeV and $|y|<2.5$. Predictions from different Monte Carlo generators based on leading-order and next-to-leading-order matrix elements for up to two additional partons interfaced with parton shower and fixed-order predictions at next-to-leading order and next-to-next-to-leading order are compared with the measured cross sections. Good agreement within the uncertainties is observed for most of the modelled quantities, in particular with the generators which use next-to-leading-order matrix elements and the more recent next-to-next-to-leading-order fixed-order predictions. Source code: ATLAS_2017_I1514251.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/ZFinder.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Projections/VetoedFinalState.hh"
6
7namespace Rivet {
8
9
10 /// Z + jets in pp at 13 TeV
11 class ATLAS_2017_I1514251 : public Analysis {
12 public:
13
14 /// Constructor
15 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2017_I1514251);
16
17 /// Book histograms and initialise projections before the run
18 void init() {
19
20 // Get options from the new option system
21 // default to combined.
22 _mode = 2;
23 if ( getOption("LMODE") == "EL" ) _mode = 0;
24 if ( getOption("LMODE") == "MU" ) _mode = 1;
25 if ( getOption("LMODE") == "EMU" ) _mode = 2;
26
27 const FinalState fs;
28
29 Cut cuts = (Cuts::pT > 25*GeV) && (Cuts::abseta < 2.5);
30
31 ZFinder zeefinder(fs, cuts, PID::ELECTRON, 71*GeV, 111*GeV);
32 ZFinder zmumufinder(fs, cuts, PID::MUON, 71*GeV, 111*GeV);
33 declare(zeefinder, "zeefinder");
34 declare(zmumufinder, "zmumufinder");
35
36 // Define veto FS in order to prevent Z-decay products entering the jet algorithm
37 VetoedFinalState had_fs;
38 had_fs.addVetoOnThisFinalState(zeefinder);
39 had_fs.addVetoOnThisFinalState(zmumufinder);
40 FastJets jets(had_fs, FastJets::ANTIKT, 0.4, JetAlg::Muons::ALL, JetAlg::Invisibles::DECAY);
41 declare(jets, "jets");
42
43 // individual channels
44 book(_h_Njets_excl, _mode + 1, 1, 1);
45 book(_h_Njets, _mode + 4, 1, 1);
46 book(_h_Njets_Ratio, _mode + 7, 1, 1, true);
47
48 book(_h_leading_jet_pT_eq1jet, _mode + 10, 1, 1);
49 book(_h_leading_jet_pT , _mode + 13, 1, 1);
50 book(_h_leading_jet_pT_2jet , _mode + 16, 1, 1);
51 book(_h_leading_jet_pT_3jet , _mode + 19, 1, 1);
52 book(_h_leading_jet_pT_4jet , _mode + 22, 1, 1);
53 book(_h_leading_jet_rap , _mode + 25, 1, 1);
54 book(_h_HT , _mode + 28, 1, 1);
55 book(_h_jet_dphi , _mode + 31, 1, 1);
56 book(_h_jet_mass , _mode + 34, 1, 1);
57
58 }
59
60
61
62 /// Perform the per-event analysis
63 void analyze(const Event& event) {
64
65 const ZFinder& zeefinder = apply<ZFinder>(event, "zeefinder");
66 const ZFinder& zmumufinder = apply<ZFinder>(event, "zmumufinder");
67
68 const Particles& zees = zeefinder.bosons();
69 const Particles& zmumus = zmumufinder.bosons();
70
71 //Veto Z->mumu in electron mode, and vice versa:
72 if (_mode==0 && (zees.size()!=1 || zmumus.size() ) ) vetoEvent;
73
74 if (_mode==1 && (zees.size() || zmumus.size()!=1 ) ) vetoEvent;
75
76 if (zees.size() + zmumus.size() != 1) {
77 // Running in combined mode, we did not find exactly one Z. Not good.
78 MSG_DEBUG("Did not find exactly one good Z candidate");
79 vetoEvent;
80 }
81
82 // Find the (dressed!) leptons
83 const Particles& leptons = zees.size() ? zeefinder.constituents() : zmumufinder.constituents();
84 if (leptons.size() != 2) vetoEvent;
85
86 Jets jets = apply<JetAlg>(event, "jets").jetsByPt(Cuts::pT > 30*GeV && Cuts::absrap < 2.5);
87
88 bool veto = false;
89 for(const Jet& j : jets) {
90 for(const Particle& l : leptons) { veto |= deltaR(j, l) < 0.4; }
91 }
92 if (veto) vetoEvent;
93
94 double HT=0;
95 for(const Particle& l : leptons) { HT += l.pT(); }
96
97 const size_t Njets = jets.size();
98 _h_Njets_excl->fill(Njets);
99 for(size_t i = 0; i <= Njets; ++i) { _h_Njets->fill(i); }
100
101 if (Njets < 1) vetoEvent;
102
103
104 for(size_t i = 0; i < Njets; ++i) { HT += jets[i].pT(); }
105 const double pT = jets[0].pT();
106 const double rap = jets[0].rapidity();
107
108 _h_HT->fill(HT);
109 _h_leading_jet_rap->fill(fabs(rap));
110 _h_leading_jet_pT->fill(pT);
111 if (Njets == 1) _h_leading_jet_pT_eq1jet->fill(pT);
112 if (Njets > 1) {
113 _h_leading_jet_pT_2jet->fill(pT);
114 _h_jet_dphi->fill( deltaPhi(jets[0], jets[1]));
115 _h_jet_mass->fill( (jets[0].momentum()+jets[1].momentum()).mass() );
116 }
117
118 if (Njets > 2) _h_leading_jet_pT_3jet->fill(pT);
119 if (Njets > 3) _h_leading_jet_pT_4jet->fill(pT);
120
121 }
122
123 void finalize() {
124 for (size_t i = 0; i < _h_Njets->numBins()-2; ++i) {
125 double n = _h_Njets->bin(i + 1).sumW();
126 double dN = _h_Njets->bin(i + 1).sumW2();
127 double d = _h_Njets->bin(i).sumW();
128 double dD = _h_Njets->bin(i).sumW2();
129 double r = safediv(n, d);
130 double e = sqrt( safediv(r * (1 - r), d) );
131 if ( _h_Njets->effNumEntries() != _h_Njets->numEntries() ) {
132 // use F. James's approximation for weighted events:
133 e = sqrt( safediv((1 - 2 * r) * dN + r * r * dD, d * d) );
134 }
135 _h_Njets_Ratio->point(i).setY(r, e);
136 }
137
138 // when running in combined mode, need to average to get lepton xsec
139 double normfac = crossSectionPerEvent();
140 if (_mode == 2) normfac = 0.5*normfac;
141
142 scale(_h_Njets, normfac );
143 scale(_h_Njets_excl, normfac );
144 scale(_h_HT, normfac );
145 scale(_h_leading_jet_rap, normfac );
146 scale(_h_leading_jet_pT, normfac );
147 scale(_h_leading_jet_pT_eq1jet, normfac );
148 scale(_h_leading_jet_pT_2jet, normfac );
149 scale(_h_leading_jet_pT_3jet, normfac );
150 scale(_h_leading_jet_pT_4jet, normfac );
151 scale(_h_jet_dphi, normfac );
152 scale(_h_jet_mass, normfac );
153
154 }
155
156 //@}
157
158
159 protected:
160
161 size_t _mode;
162
163
164 private:
165
166 Scatter2DPtr _h_Njets_Ratio;
167 Histo1DPtr _h_Njets;
168 Scatter2DPtr _h_Njets_excl_Ratio;
169 Histo1DPtr _h_Njets_excl;
170 Histo1DPtr _h_HT;
171 Histo1DPtr _h_leading_jet_rap;
172 Histo1DPtr _h_leading_jet_pT;
173 Histo1DPtr _h_leading_jet_pT_eq1jet;
174 Histo1DPtr _h_leading_jet_pT_2jet;
175 Histo1DPtr _h_leading_jet_pT_3jet;
176 Histo1DPtr _h_leading_jet_pT_4jet;
177 Histo1DPtr _h_jet_dphi;
178 Histo1DPtr _h_jet_mass;
179
180 };
181
182
183 RIVET_DECLARE_PLUGIN(ATLAS_2017_I1514251);
184
185}
|