Rivet analyses referenceATLAS_2014_I1279489Measurements of electroweak production of dijets + $Z$ boson, and distributions sensitive to vector boson fusionExperiment: ATLAS (LHC) Inspire ID: 1279489 Status: VALIDATED Authors:
Beam energies: (4000.0, 4000.0) GeV Run details:
Measurements differential distributions for inclusive $Z$-boson-plus-dijet production are performed in five fiducial regions, each with different sensitivity to the electroweak contribution. Measured distributions include the differential cross section as a function of the dijet invariant mass, the differential cross section and a function of the dijet rapidity separation, the differential cross section as a function of the number of jets in the rapidity interval bounded by the two leading jets. Other measurements include the jet veto effiency as a function of the dijet invariant mass and rapdity separation, the normalized transverse momentum balance cut efficiency, and the average number of jets falling into the rapidity interval boundd by the two leading jets, as a function of dijet invariant mass and dijet rapidity separation. Source code: ATLAS_2014_I1279489.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/VetoedFinalState.hh"
5#include "Rivet/Projections/FastJets.hh"
6#include "Rivet/Projections/LeptonFinder.hh"
7
8namespace Rivet {
9
10
11 struct Plots {
12 string label;
13
14 Histo1DPtr h_dy;
15 Histo1DPtr h_mjj;
16 Histo1DPtr h_njets;
17 Histo1DPtr h_dphijj;
18 Histo1DPtr h_ptbal;
19
20 Histo1DPtr h_jetveto_mjj_veto;
21 Histo1DPtr h_jetveto_mjj_inc;
22 Histo1DPtr h_jetveto_dy_veto;
23 Histo1DPtr h_jetveto_dy_inc;
24
25 Histo1DPtr h_ptbaleff_mjj_veto;
26 Histo1DPtr h_ptbaleff_mjj_inc;
27 Histo1DPtr h_ptbaleff_dy_veto;
28 Histo1DPtr h_ptbaleff_dy_inc;
29
30 Estimate1DPtr s_jetveto_mjj;
31 Estimate1DPtr s_jetveto_dy;
32
33 Estimate1DPtr s_ptbaleff_mjj;
34 Estimate1DPtr s_ptbaleff_dy;
35
36 Profile1DPtr p_avgnjets_dy;
37 Profile1DPtr p_avgnjets_mjj;
38 };
39
40
41 struct Variables {
42
43 Variables(const Jets& jets, const Particle* lep1, const Particle* lep2) {
44 FourMomentum j1 = jets[0].momentum();
45 FourMomentum j2 = jets[1].momentum();
46 jet1pt = j1.pT();
47 jet2pt = j2.pT();
48 assert(jet1pt > jet2pt);
49
50 zpt = (lep1->mom() + lep2->mom()).pT();
51
52 deltay = fabs(j1.rapidity() - j2.rapidity());
53 mjj = (j1 + j2).mass();
54 deltaphijj = deltaPhi(j1, j2) / PI;
55
56 FourMomentum gapjet(0., 0., 0., 0.);
57 ngapjets = _getNumGapJets(jets, gapjet);
58
59 double ptbal_vec = (j1 + j2 + lep1->mom() + lep2->mom()).pT();
60 double ptbal_sc = j1.pT() + j2.pT() + lep1->pT() + lep2->pT();
61 ptbalance2 = ptbal_vec / ptbal_sc;
62
63 double ptbal3_vec = (j1 + j2 + gapjet + lep1->mom() + lep2->mom()).pT();
64 double ptbal3_sc = j1.pT() + j2.pT() + gapjet.pT() + lep1->pT() + lep2->pT();
65 ptbalance3 = ptbal3_vec / ptbal3_sc;
66
67 pass_jetveto = gapjet.pT() < 25.0*GeV;
68 pass_ptbaleff = ptbalance2 < 0.15;
69 }
70
71
72 double jet1pt;
73 double jet2pt;
74 double zpt;
75
76 double deltay;
77 double mjj;
78 double deltaphijj;
79 double ptbalance2;
80 double ptbalance3;
81 int ngapjets;
82
83 double dilepton_dr;
84
85 bool pass_jetveto;
86 bool pass_ptbaleff;
87
88
89 private:
90
91 bool _isBetween(const Jet& probe, const Jet& boundary1, const Jet& boundary2) {
92 double y_p = probe.rap();
93 double y_b1 = boundary1.rap();
94 double y_b2 = boundary2.rap();
95
96 double y_min = std::min(y_b1, y_b2);
97 double y_max = std::max(y_b1, y_b2);
98
99 if (y_p > y_min && y_p < y_max) return true;
100 else return false;
101 }
102
103 int _getNumGapJets(const Jets& jets, FourMomentum& thirdJet) {
104 if (jets.size() < 2) return 0;
105 int n_between = 0;
106 // Start loop at the 3rd hardest pT jet
107 for (size_t i = 2; i < jets.size(); ++i) {
108 // If this jet is between the boundary jets and is hard enough, increment counter
109 if (_isBetween(jets[i], jets[0], jets[1])) {
110 if (n_between == 0) thirdJet = jets[i].momentum();
111 ++n_between;
112 }
113 }
114 return n_between;
115 }
116
117 };
118
119
120
121 class ATLAS_2014_I1279489 : public Analysis {
122 public:
123
124 /// Constructor
125 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2014_I1279489);
126
127
128 /// Book histograms and initialise projections before the run
129 void init() {
130
131 FinalState fs(Cuts::abseta < 5.);
132
133 LeptonFinder dressed_electrons(0.1, Cuts::abspid == PID::ELECTRON && Cuts::abseta < 2.47 && Cuts::pT > 25*GeV);
134 declare(dressed_electrons, "DressedElectrons");
135
136 LeptonFinder dressed_muons(0.1, Cuts::abspid == PID::MUON && Cuts::abseta < 2.47 && Cuts::pT > 25*GeV);
137 declare(dressed_muons, "DressedMuons");
138
139 FastJets jets(fs, JetAlg::ANTIKT, 0.4);
140 declare(jets, "Jets");
141
142 initialisePlots(baseline_plots, "baseline");
143 initialisePlots(highpt_plots, "highpt");
144 initialisePlots(search_plots, "search");
145 initialisePlots(control_plots, "control");
146 initialisePlots(highmass_plots, "highmass");
147 }
148
149
150 void initialisePlots(Plots& plots, const string& phase_space){
151 /****************************************
152 * Plot labeling: *
153 * format = d0_-x0_-y0_ *
154 * d01 = baseline fiducial region *
155 * d02 = high-pt fiducial region *
156 * d03 = search fiducial region *
157 * d04 = control fiducial region *
158 * d05 = high-mass fiducial region *
159 * *
160 * x01 = mjj on x-axis *
161 * x02 = delta-y on x-axis *
162 * x03 = njets on x-axis *
163 * x04 = dphijj on x-axis *
164 * x05 = ptbalance on x-axis *
165 * *
166 * y01 = differential cross-section *
167 * y02 = jet veto efficiency *
168 * y03 = ptbalance efficiency *
169 * y04 = average njets *
170 ****************************************/
171 plots.label = phase_space;
172
173 if (phase_space=="baseline") {
174 book(plots.h_mjj ,1, 1, 1);
175 book(plots.h_dy ,3, 1, 1);
176
177 book(plots.h_jetveto_mjj_veto,"_jetveto_mjj_baseline_veto", refData(8,1,1));
178 book(plots.h_jetveto_mjj_inc, "_jetveto_mjj_baseline_inc", refData(8,1,1));
179 book(plots.h_jetveto_dy_veto, "_jetveto_dy_baseline_veto", refData(9,1,1));
180 book(plots.h_jetveto_dy_inc, "_jetveto_dy_baseline_inc", refData(9,1,1));
181
182 book(plots.h_ptbaleff_mjj_veto, "_ptbaleff_mjj_baseline_veto", refData(12,1,1));
183 book(plots.h_ptbaleff_mjj_inc, "_ptbaleff_mjj_baseline_inc", refData(12,1,1));
184 book(plots.h_ptbaleff_dy_veto, "_ptbaleff_dy_baseline_veto", refData(13,1,1));
185 book(plots.h_ptbaleff_dy_inc, "_ptbaleff_dy_baseline_inc", refData(13,1,1));
186
187 book(plots.s_jetveto_mjj, 8, 1, 1);
188 book(plots.s_jetveto_dy, 9, 1, 1);
189 book(plots.s_ptbaleff_mjj, 12, 1, 1);
190 book(plots.s_ptbaleff_dy, 13, 1, 1);
191
192 book(plots.p_avgnjets_mjj ,10,1,1);
193 book(plots.p_avgnjets_dy ,11,1,1);
194 }
195
196 if (phase_space=="highpt") {
197 book(plots.h_mjj , 14, 1, 1);
198 book(plots.h_dy , 16, 1, 1);
199
200 book(plots.h_jetveto_mjj_veto, "_jetveto_mjj_highpt_veto", refData(18,1,1));
201 book(plots.h_jetveto_mjj_inc, "_jetveto_mjj_highpt_inc", refData(18,1,1));
202 book(plots.h_jetveto_dy_veto, "_jetveto_dy_highpt_veto", refData(19,1,1));
203 book(plots.h_jetveto_dy_inc, "_jetveto_dy_highpt_inc", refData(19,1,1));
204
205 book(plots.h_ptbaleff_mjj_veto, "_ptbaleff_mjj_highpt_veto", refData(22,1,1));
206 book(plots.h_ptbaleff_mjj_inc, "_ptbaleff_mjj_highpt_inc", refData(22,1,1));
207 book(plots.h_ptbaleff_dy_veto, "_ptbaleff_dy_highpt_veto", refData(23,1,1));
208 book(plots.h_ptbaleff_dy_inc, "_ptbaleff_dy_highpt_inc", refData(23,1,1));
209
210 book(plots.s_jetveto_mjj, 18, 1, 1);
211 book(plots.s_jetveto_dy, 19, 1, 1);
212 book(plots.s_ptbaleff_mjj, 22, 1, 1);
213 book(plots.s_ptbaleff_dy, 23, 1, 1);
214
215 book(plots.p_avgnjets_mjj , 20,1,1);
216 book(plots.p_avgnjets_dy , 21,1,1);
217 }
218
219 if (phase_space=="search") {
220 book(plots.h_mjj , 2,1,1);
221 book(plots.h_dy , 4,1,1);
222 }
223
224 if (phase_space=="control") {
225 book(plots.h_mjj , 15,1,1);
226 book(plots.h_dy , 17,1,1);
227 }
228
229 if (phase_space=="highmass") {
230 book(plots.h_njets , 5,1,1);
231 book(plots.h_dphijj , 7,1,1);
232 book(plots.h_ptbal , 6,1,1);
233 }
234 }
235
236
237
238 /// Perform the per-event analysis
239 void analyze(const Event& event) {
240
241 // Make sure that we have a Z-candidate:
242 const Particle *lep1 = nullptr, *lep2 = nullptr;
243 //
244 const DressedLeptons& muons = apply<LeptonFinder>(event, "DressedMuons").dressedLeptons();
245 if (muons.size() == 2) {
246 const FourMomentum dimuon = muons[0].mom() + muons[1].mom();
247 if ( inRange(dimuon.mass()/GeV, 81.0, 101.0) && PID::charge3(muons[0].pid()) != PID::charge3(muons[1].pid()) ) {
248 lep1 = &muons[0];
249 lep2 = &muons[1];
250 }
251 }
252 //
253 const DressedLeptons& electrons = apply<LeptonFinder>(event, "DressedElectrons").dressedLeptons();
254 if (electrons.size() == 2) {
255 const FourMomentum dielectron = electrons[0].mom() + electrons[1].mom();
256 if ( inRange(dielectron.mass()/GeV, 81.0, 101.0) && PID::charge3(electrons[0].pid()) != PID::charge3(electrons[1].pid()) ) {
257 if (lep1 && lep2) {
258 MSG_INFO("Found Z candidates using both electrons and muons! Continuing with the muon-channel candidate");
259 } else {
260 lep1 = &electrons[0];
261 lep2 = &electrons[1];
262 }
263 }
264 }
265 // If there's no Z-candidate, we won't use this event:
266 if (!lep1 || !lep2) vetoEvent;
267
268
269 // Do lepton-jet overlap removal:
270 Jets jets = apply<FastJets>(event, "Jets").jetsByPt(Cuts::pT > 25*GeV && Cuts::absrap < 4.4);
271 idiscardIfAnyDeltaRLess(jets, muons, 0.3);
272 idiscardIfAnyDeltaRLess(jets, electrons, 0.3);
273
274 // If we don't have at least 2 good jets, we won't use this event.
275 if (jets.size() < 2) vetoEvent;
276
277
278 // Plotting, using variables and histo classes calculated by the Variables object constructor
279 Variables vars(jets, lep1, lep2);
280 bool pass_baseline = (vars.jet1pt > 55*GeV && vars.jet2pt > 45*GeV);
281 bool pass_highpt = (vars.jet1pt > 85*GeV && vars.jet2pt > 75*GeV);
282 bool pass_highmass = (pass_baseline && vars.mjj > 1000*GeV);
283 bool pass_search = (pass_baseline && vars.zpt > 20*GeV && vars.ngapjets == 0 && vars.ptbalance2 < 0.15 && vars.mjj > 250*GeV);
284 bool pass_control = (pass_baseline && vars.zpt > 20*GeV && vars.ngapjets > 0 && vars.ptbalance3 < 0.15 && vars.mjj > 250*GeV);
285 //
286 if (pass_baseline) fillPlots(vars, baseline_plots, "baseline");
287 if (pass_highpt) fillPlots(vars, highpt_plots, "highpt");
288 if (pass_highmass) fillPlots(vars, highmass_plots, "highmass");
289 if (pass_search) fillPlots(vars, search_plots, "search");
290 if (pass_control) fillPlots(vars, control_plots, "control");
291 }
292
293
294 void fillPlots(const Variables& vars, Plots& plots, string phase_space) {
295 if (phase_space == "baseline" || phase_space == "highpt" || phase_space == "search" || phase_space == "control") {
296 plots.h_dy->fill(vars.deltay);
297 plots.h_mjj->fill(vars.mjj);
298 }
299
300 if (phase_space == "baseline" || phase_space == "highpt") {
301 if (vars.pass_jetveto) {
302 plots.h_jetveto_dy_veto->fill(vars.deltay);
303 plots.h_jetveto_mjj_veto->fill(vars.mjj);
304 }
305 plots.h_jetveto_dy_inc->fill(vars.deltay);
306 plots.h_jetveto_mjj_inc->fill(vars.mjj);
307
308 if (vars.pass_ptbaleff) {
309 plots.h_ptbaleff_mjj_veto->fill(vars.mjj);
310 plots.h_ptbaleff_dy_veto->fill(vars.deltay);
311 }
312 plots.h_ptbaleff_mjj_inc->fill(vars.mjj);
313 plots.h_ptbaleff_dy_inc->fill(vars.deltay);
314
315 plots.p_avgnjets_dy->fill(vars.deltay, vars.ngapjets);
316 plots.p_avgnjets_mjj->fill(vars.mjj, vars.ngapjets);
317 }
318
319 if (phase_space == "highmass") {
320 plots.h_njets->fill(vars.ngapjets);
321 plots.h_dphijj->fill(vars.deltaphijj);
322 plots.h_ptbal->fill(vars.ptbalance2);
323 }
324 }
325
326
327 /// Normalise histograms etc., after the run
328 void finalize() {
329 finalizePlots(baseline_plots);
330 finalizePlots(highpt_plots);
331 finalizePlots(search_plots);
332 finalizePlots(control_plots);
333 finalizePlots(highmass_plots);
334 finalizeEfficiencies(baseline_plots);
335 finalizeEfficiencies(highpt_plots);
336 }
337
338 void finalizePlots(Plots& plots) {
339 if (plots.h_dy) normalize(plots.h_dy);
340 if (plots.h_mjj) normalize(plots.h_mjj);
341 if (plots.h_dphijj) normalize(plots.h_dphijj);
342 if (plots.h_njets) normalize(plots.h_njets);
343 if (plots.h_ptbal) normalize(plots.h_ptbal);
344 }
345
346 void finalizeEfficiencies(Plots& plots) {
347
348 if (plots.label != "baseline" && plots.label != "highpt") return;
349 //size_t offset = plots.label == "baseline"? 0 : 10;
350
351 if (plots.h_jetveto_mjj_veto && plots.h_jetveto_mjj_inc) {
352 divide(plots.h_jetveto_mjj_veto, plots.h_jetveto_mjj_inc, plots.s_jetveto_mjj);
353 }
354 //getEstimate1D(8+offset, 1, 1)->addAnnotation("InclusiveSumWeights", plots.h_jetveto_mjj_inc->integral());
355
356 if (plots.h_jetveto_dy_veto && plots.h_jetveto_dy_inc) {
357 divide(plots.h_jetveto_dy_veto, plots.h_jetveto_dy_inc, plots.s_jetveto_dy);
358 }
359
360 if (plots.h_ptbaleff_mjj_veto && plots.h_ptbaleff_mjj_inc) {
361 divide(plots.h_ptbaleff_mjj_veto, plots.h_ptbaleff_mjj_inc, plots.s_ptbaleff_mjj);
362 }
363
364 if (plots.h_ptbaleff_dy_veto && plots.h_ptbaleff_dy_inc) {
365 divide(plots.h_ptbaleff_dy_veto, plots.h_ptbaleff_dy_inc, plots.s_ptbaleff_dy);
366 }
367 }
368
369
370 private:
371
372 Plots baseline_plots;
373 Plots highpt_plots;
374 Plots search_plots;
375 Plots control_plots;
376 Plots highmass_plots;
377
378 };
379
380
381 RIVET_DECLARE_PLUGIN(ATLAS_2014_I1279489);
382
383}
|