Rivet analyses referenceATLAS_2017_I1627873EW Zjj using early Run-2 dataExperiment: ATLAS (LHC) Inspire ID: 1627873 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
The cross-section for the production of two jets in association with a leptonically decaying Z boson ($Zjj$) is measured in proton-proton collisions at a centre-of-mass energy of 13 TeV, using data recorded with the ATLAS detector at the Large Hadron Collider, corresponding to an integrated luminosity of 3.2 fbi$^{-1}$. The electroweak $Zjj$ cross-section is extracted in a fiducial region chosen to enhance the electroweak contribution relative to the dominant Drell–Yan $Zjj$ process, which is constrained using a data-driven approach. The measured fiducial electroweak cross-section is $\sigma_{EWZjj}$ = 119$\pm$16(stat.)$\pm$20(syst.)$\pm$2(lumi.) fb for dijet invariant mass greater than 250 GeV, and 34.2$\pm$5.8(stat.)$\pm$5.5(syst.)$\pm$0.7(lumi.) fb for dijet invariant mass greater than 1 TeV. Standard Model predictions are in agreement with the measurements. The inclusive $Zjj$ cross-section is also measured in six different fiducial regions with varying contributions from electroweak and Drell-Yan $Zjj$ production. Source code: ATLAS_2017_I1627873.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FastJets.hh"
4#include "Rivet/Projections/FinalState.hh"
5#include "Rivet/Projections/PromptFinalState.hh"
6#include "Rivet/Projections/VetoedFinalState.hh"
7#include "Rivet/Projections/DressedLeptons.hh"
8
9namespace Rivet {
10
11
12 /// @brief EW Zjj using early Run-2 data
13 class ATLAS_2017_I1627873 : public Analysis {
14 public:
15
16 /// Constructor
17 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2017_I1627873);
18
19
20 /// Book histograms and initialise projections before the run
21 void init() {
22
23 _mode = 0;
24 if ( getOption("TYPE") == "EW_ONLY" ) _mode = 1;
25
26 FinalState fs(Cuts::abseta < 5.0);
27
28 FinalState photon_fs(Cuts::abspid == PID::PHOTON);
29
30 PromptFinalState electron_fs(Cuts::abspid == PID::ELECTRON);
31 PromptFinalState muon_fs(Cuts::abspid == PID::MUON);
32
33 DressedLeptons dressed_electrons(photon_fs, electron_fs, 0.1, Cuts::abseta < 2.47 && Cuts::pT > 25*GeV);
34 declare(dressed_electrons, "DressedElectrons");
35
36 DressedLeptons dressed_muons(photon_fs, muon_fs, 0.1, Cuts::abseta < 2.47 && Cuts::pT > 25*GeV);
37 declare(dressed_muons, "DressedMuons");
38
39 VetoedFinalState remfs(fs);
40 remfs.addVetoOnThisFinalState(dressed_electrons);
41 remfs.addVetoOnThisFinalState(dressed_muons);
42
43 FastJets jets(remfs, FastJets::ANTIKT, 0.4, JetAlg::Muons::ALL, JetAlg::Invisibles::ALL);
44 declare(jets, "Jets");
45
46 if (_mode) book(_h["zjj-ew"], 3, 1, 1);
47 else book(_h["zjj"], 2, 1, 1);
48 }
49
50
51 /// Perform the per-event analysis
52 void analyze(const Event& event) {
53
54 const Jets& jets = apply<FastJets>(event, "Jets").jetsByPt(Cuts::pT > 25*GeV && Cuts::abseta < 4.4);
55 vector<DressedLepton> electrons = apply<DressedLeptons>(event, "DressedElectrons").dressedLeptons();
56 vector<DressedLepton> muons = apply<DressedLeptons>(event, "DressedMuons").dressedLeptons();
57
58 // Overlap Removal
59 idiscardIfAnyDeltaRLess(electrons, jets, 0.4);
60 idiscardIfAnyDeltaRLess(muons, jets, 0.4);
61
62 Particle lep1, lep2;
63 if (electrons.size() == 2 && muons.empty()) {
64 lep1 = electrons[0]; lep2 = electrons[1];
65 if ( lep1.charge3() == lep2.charge3() ) vetoEvent;
66 }
67 else if (electrons.empty() && muons.size() == 2) {
68 lep1 = muons[0]; lep2 = muons[1];
69 if (lep1.charge3() == lep2.charge3()) vetoEvent;
70 }
71 else vetoEvent;
72
73 if (jets.size() < 2) vetoEvent;
74
75 const FourMomentum dilepton = lep1.mom()+lep2.mom();
76 if ( !inRange(dilepton.mass(), 81.0*GeV, 101.0*GeV) ) vetoEvent;
77
78 const double jet1pt = jets[0].pT();
79 const double jet2pt = jets[1].pT();
80 const double mjj = (jets[0].mom() + jets[1].mom()).mass();
81 const double zpt = (lep1.mom() + lep2.mom()).pT();
82
83 size_t ngapjets = 0;
84 Jet thirdjet;
85 for (size_t i = 2; i < jets.size(); ++i) {
86 const Jet j = jets[i];
87 if (_isBetween(j, jets[0], jets[1])) {
88 if (!ngapjets) thirdjet = j;
89 ++ngapjets;
90 }
91 }
92
93 const double ptbal_vec = (jets[0].mom() + jets[1].mom() + lep1.mom() + lep2.mom()).pT();
94 const double ptbal_sc = jets[0].pT() + jets[1].pT() + lep1.pT() + lep2.pT();
95 const double ptbalance2 = ptbal_vec / ptbal_sc;
96
97 const double ptbal3_vec = (jets[0].mom() + jets[1].mom() + thirdjet.mom() + lep1.mom() + lep2.mom()).pT();
98 const double ptbal3_sc = jets[0].pT() + jets[1].pT() + thirdjet.pT() + lep1.pT() + lep2.pT();
99 const double ptbalance3 = ptbal3_vec / ptbal3_sc;
100
101
102
103 //categories: baseline, high-PT, EW-enriched, QCD-enriched, high-mass, EW-enriched and high-mass
104 if(!(jet1pt > 55*GeV && jet2pt > 45*GeV)) vetoEvent;
105
106 if (_mode) {
107 if (zpt > 20.0*GeV && !ngapjets && ptbalance2 < 0.15 && mjj > 250.0*GeV) _h["zjj-ew"]->fillBin(0);
108 if (zpt > 20.0*GeV && !ngapjets && ptbalance2 < 0.15 && mjj > 1000.0*GeV) _h["zjj-ew"]->fillBin(1);
109 }
110 else {
111 _h["zjj"]->fillBin(0);
112 if (jet1pt > 85.0*GeV && jet2pt > 75.0*GeV) _h["zjj"]->fillBin(1);
113 if (zpt > 20.0*GeV && ngapjets == 0 && ptbalance2 < 0.15 && mjj > 250.0*GeV) _h["zjj"]->fillBin(2);
114 if (zpt > 20.0*GeV && ngapjets && ptbalance3 < 0.15 && mjj > 250.0*GeV) _h["zjj"]->fillBin(3);
115 if (mjj > 1000.0*GeV) _h["zjj"]->fillBin(4);
116 if (zpt > 20.0*GeV && !ngapjets && ptbalance2 < 0.15 && mjj > 1000.0*GeV) _h["zjj"]->fillBin(2);
117 }
118
119 }
120
121
122
123 /// Normalise histograms etc., after the run
124 void finalize() {
125
126 double factor = crossSection()/femtobarn/sumOfWeights();
127 scale(_h, factor);
128 }
129
130 bool _isBetween(const Jet probe, const Jet boundary1, const Jet boundary2) {
131 double y_p = probe.rapidity();
132 double y_b1 = boundary1.rapidity();
133 double y_b2 = boundary2.rapidity();
134
135 double y_min = std::min(y_b1, y_b2);
136 double y_max = std::max(y_b1, y_b2);
137
138 if (y_p > y_min && y_p < y_max) return true;
139 else return false;
140 }
141
142 ///@}
143
144 private:
145
146 size_t _mode;
147
148 /// @name Histograms
149 ///@{
150 map<string, Histo1DPtr> _h;
151 ///@}
152
153 };
154
155
156 RIVET_DECLARE_PLUGIN(ATLAS_2017_I1627873);
157
158}
|