Rivet analyses referenceATLAS_2013_I1217863W/Z + gamma production at 7 TeVExperiment: ATLAS (LHC) Inspire ID: 1217863 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
Measurements of the differential fiducial cross sections for the production of a W or Z boson in association with a high-energy photon are measured using pp collisions at $\sqrt{s}=7 \text{TeV}$. The analysis uses a data sample with an integrated luminosity of 4.6/fb collected by the ATLAS detector during the 2011 LHC data-taking period. Events are selected using leptonic decays of the W or Z bosons with the requirement of an associated isolated photon. The default routine will consider the electron decay channel of the Z boson. Use LMODE to specify the decay channel directly. Source code: ATLAS_2013_I1217863.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#include "Rivet/Projections/VetoedFinalState.hh"
8#include "Rivet/Projections/LeadingParticlesFinalState.hh"
9
10namespace Rivet {
11
12
13 /// Electroweak Wjj production at 8 TeV
14 class ATLAS_2013_I1217863 : public Analysis {
15 public:
16
17 /// Constructor
18 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2013_I1217863);
19
20 /// @name Analysis methods
21 //@{
22
23 /// Book histograms and initialise projections before the run
24 void init() {
25
26 // Get options from the new option system
27 _mode = 2;
28 _doZ = true;
29 _doW = true;
30 if ( getOption("LMODE") == "EL" ) { _mode = 2;}
31 if ( getOption("LMODE") == "MU" ) _mode = 3;
32 if ( getOption("LMODE") == "ZEL" ) {
33 _mode = 2;
34 _doW = false;
35 }
36 if ( getOption("LMODE") == "ZMU" ) {
37 _mode = 3;
38 _doW = false;
39 }
40 if ( getOption("LMODE") == "WEL" ) {
41 _mode = 2;
42 _doZ = false;
43 }
44 if ( getOption("LMODE") == "WMU" ) {
45 _mode = 3;
46 _doZ = false;
47 }
48
49 FinalState fs;
50 declare(fs, "FS");
51
52 Cut cuts = Cuts::abseta < 2.47 && Cuts::pT > 25*GeV;
53
54 // Z finder
55 if (_doZ) {
56 ZFinder zf(fs, cuts, _mode==3? PID::MUON : PID::ELECTRON, 40.0*GeV, 1000.0*GeV, 0.1,
57 ZFinder::ChargedLeptons::PROMPT, ZFinder::ClusterPhotons::NODECAY, ZFinder::AddPhotons::NO);
58 declare(zf, "ZF");
59 }
60
61 if (_doW) {
62 // W finder for electrons and muons
63 WFinder wf(fs, cuts, _mode==3? PID::MUON : PID::ELECTRON, 0.0*GeV, 1000.0*GeV, 35.0*GeV, 0.1,
64 WFinder::ChargedLeptons::PROMPT, WFinder::ClusterPhotons::NODECAY, WFinder::AddPhotons::NO, WFinder::MassWindow::MT);
65 declare(wf, "WF");
66 }
67
68 // leading photon
69 LeadingParticlesFinalState photonfs(FinalState(Cuts::abseta < 2.37 && Cuts::pT > 15*GeV));
70 photonfs.addParticleId(PID::PHOTON);
71 declare(photonfs, "LeadingPhoton");
72
73 // jets
74 VetoedFinalState jet_fs(fs);
75 if (_doZ) { jet_fs.addVetoOnThisFinalState(getProjection<ZFinder>("ZF")); }
76 if (_doW) { jet_fs.addVetoOnThisFinalState(getProjection<WFinder>("WF")); }
77 jet_fs.addVetoOnThisFinalState(getProjection<LeadingParticlesFinalState>("LeadingPhoton"));
78 FastJets jets(jet_fs, FastJets::ANTIKT, 0.4, JetAlg::Muons::ALL, JetAlg::Invisibles::NONE);
79 declare(jets, "Jets");
80
81 // FS excluding the leading photon
82 VetoedFinalState vfs(fs);
83 vfs.addVetoOnThisFinalState(photonfs);
84 declare(vfs, "isolatedFS");
85
86
87 // Book histograms
88 if (_doZ) {
89 book(_hist_EgammaT_inclZ ,11, 1, _mode); // dSigma / dE^gamma_T for Njet >= 0
90 book(_hist_EgammaT_exclZ ,12, 1, _mode); // dSigma / dE^gamma_T for Njet = 0
91 book(_hist_Njet_EgammaT15Z ,17, 1, _mode); // dSigma / dNjet for E^gamma_T >= 15
92 book(_hist_Njet_EgammaT60Z ,18, 1, _mode); // dSigma / dNjet for E^gamma_T >= 60
93 book(_hist_mZgamma ,20, 1, _mode); // dSigma / dm^{Zgamma}
94 }
95 if (_doW){
96 book(_hist_EgammaT_inclW , 7, 1, _mode); // dSigma / dE^gamma_T for Njet >= 0
97 book(_hist_EgammaT_exclW , 8, 1, _mode); // dSigma / dE^gamma_T for Njet = 0
98 book(_hist_Njet_EgammaT15W ,15, 1, _mode); // dSigma / dNjet for E^gamma_T >= 15
99 book(_hist_Njet_EgammaT60W ,16, 1, _mode); // dSigma / dNjet for E^gamma_T >= 60
100 book(_hist_mWgammaT ,19, 1, _mode); // dSigma / dm^{Zgamma}
101 }
102
103 }
104
105
106 /// Perform the per-event analysis
107 void analyze(const Event& event) {
108
109 // retrieve leading photon
110 Particles photons = apply<LeadingParticlesFinalState>(event, "LeadingPhoton").particles();
111 if (photons.size() != 1) vetoEvent;
112 const Particle& leadingPhoton = photons[0];
113 if (leadingPhoton.Et() < 15.0*GeV) vetoEvent;
114 if (leadingPhoton.abseta() > 2.37) vetoEvent;
115
116 // check photon isolation
117 double coneEnergy(0.0);
118 Particles fs = apply<VetoedFinalState>(event, "isolatedFS").particles();
119 for (const Particle& p : fs) {
120 if ( deltaR(leadingPhoton, p) < 0.4 ) coneEnergy += p.E();
121 }
122 if (coneEnergy / leadingPhoton.E() >= 0.5 ) vetoEvent;
123
124 if (_doW) {
125 // retrieve W boson candidate
126 const WFinder& wf = apply<WFinder>(event, "WF");
127 if ( wf.bosons().size() == 1 ) {
128
129 // retrieve constituent neutrino
130 const Particle& neutrino = wf.constituentNeutrino();
131 if ( (neutrino.pT() > 35.0*GeV) ) {
132
133 // retrieve constituent lepton
134 const Particle& lepton = wf.constituentLepton();
135 if ( lepton.pT() > 25.0*GeV && lepton.abseta() < 2.47 ) {
136
137 // check photon-lepton overlap
138 if ( deltaR(leadingPhoton, lepton) > 0.7 ) {
139
140 // count jets
141 const FastJets& jetfs = apply<FastJets>(event, "Jets");
142 Jets jets = jetfs.jets(cmpMomByEt);
143 int goodJets = 0;
144 for (const Jet& j : jets) {
145 if ( !(j.Et() > 30.0*GeV) ) break;
146 if ( (j.abseta() < 4.4) && \
147 (deltaR(leadingPhoton, j) > 0.3) && \
148 (deltaR(lepton, j) > 0.3) ) ++goodJets;
149 }
150
151 double Njets = double(goodJets) + 0.5;
152 double photonEt = leadingPhoton.Et()*GeV;
153
154 const FourMomentum& lep_gamma = lepton.momentum() + leadingPhoton.momentum();
155 double term1 = sqrt(lep_gamma.mass2() + lep_gamma.pT2()) + neutrino.Et();
156 double term2 = (lep_gamma + neutrino.momentum()).pT2();
157 double mWgammaT = sqrt(term1 * term1 - term2) * GeV;
158
159 _hist_EgammaT_inclW->fill(photonEt);
160
161 _hist_Njet_EgammaT15W->fill(Njets);
162
163 if ( !goodJets ) _hist_EgammaT_exclW->fill(photonEt);
164
165 if (photonEt > 40.0*GeV) {
166 _hist_mWgammaT->fill(mWgammaT);
167 if (photonEt > 60.0*GeV) _hist_Njet_EgammaT60W->fill(Njets);
168 }
169 }
170 }
171 }
172 }
173 }
174
175 if (_doZ ){
176
177 // retrieve Z boson candidate
178 const ZFinder& zf = apply<ZFinder>(event, "ZF");
179 if ( zf.bosons().size() == 1 ) {
180 const Particle& Zboson = zf.boson();
181 if ( (Zboson.mass() > 40.0*GeV) ) {
182
183 // check charge of constituent leptons
184 const Particles& leptons = zf.constituents();
185 if (leptons.size() == 2 && leptons[0].charge() * leptons[1].charge() < 0.) {
186
187 bool lpass = true;
188 // check photon-lepton overlap
189 for (const Particle& p : leptons) {
190 if ( !(p.pT() > 25.0*GeV && p.abseta() < 2.47 && deltaR(leadingPhoton, p) > 0.7) ) lpass = false;
191 }
192 if ( lpass ) {
193
194 // count jets
195 const FastJets& jetfs = apply<FastJets>(event, "Jets");
196 Jets jets = jetfs.jets(cmpMomByEt);
197 int goodJets = 0;
198 for (const Jet& j : jets) {
199 if ( !(j.Et() > 30.0*GeV) ) break;
200 if ( (j.abseta() < 4.4) && \
201 (deltaR(leadingPhoton, j) > 0.3) && \
202 (deltaR(leptons[0], j) > 0.3) && \
203 (deltaR(leptons[1], j) > 0.3) ) ++goodJets;
204 }
205
206 double Njets = double(goodJets) + 0.5;
207 double photonEt = leadingPhoton.Et()*GeV;
208 double mZgamma = (Zboson.momentum() + leadingPhoton.momentum()).mass() * GeV;
209
210 _hist_EgammaT_inclZ->fill(photonEt);
211
212 _hist_Njet_EgammaT15Z->fill(Njets);
213
214 if ( !goodJets ) _hist_EgammaT_exclZ->fill(photonEt);
215
216 if (photonEt >= 40.0*GeV) {
217 _hist_mZgamma->fill(mZgamma);
218 if (photonEt >= 60.0*GeV) _hist_Njet_EgammaT60Z->fill(Njets);
219 }
220 }
221 }
222 }
223 }
224 }
225
226 }
227
228
229 /// Normalise histograms etc., after the run
230 void finalize() {
231
232 const double xs_fb = crossSection()/femtobarn;
233 const double sumw = sumOfWeights();
234 const double sf = xs_fb / sumw;
235
236 if (_doZ) {
237 scale(_hist_EgammaT_exclZ, sf);
238 scale(_hist_EgammaT_inclZ, sf);
239 normalize(_hist_Njet_EgammaT15Z);
240 normalize(_hist_Njet_EgammaT60Z);
241 normalize(_hist_mZgamma);
242 }
243
244 if (_doW) {
245 scale(_hist_EgammaT_exclW, sf);
246 scale(_hist_EgammaT_inclW, sf);
247 normalize(_hist_Njet_EgammaT15W);
248 normalize(_hist_Njet_EgammaT60W);
249 normalize(_hist_mWgammaT);
250 }
251
252 }
253
254 //@}
255
256 protected:
257
258 size_t _mode;
259 bool _doW;
260 bool _doZ;
261
262 private:
263
264 /// @name Histograms
265 //@{
266
267 Histo1DPtr _hist_EgammaT_inclZ;
268 Histo1DPtr _hist_EgammaT_exclZ;
269 Histo1DPtr _hist_Njet_EgammaT15Z;
270 Histo1DPtr _hist_Njet_EgammaT60Z;
271 Histo1DPtr _hist_mZgamma;
272
273 Histo1DPtr _hist_EgammaT_inclW;
274 Histo1DPtr _hist_EgammaT_exclW;
275 Histo1DPtr _hist_Njet_EgammaT15W;
276 Histo1DPtr _hist_Njet_EgammaT60W;
277 Histo1DPtr _hist_mWgammaT;
278
279 //@}
280
281 };
282
283 RIVET_DECLARE_PLUGIN(ATLAS_2013_I1217863);
284
285}
|