|
Rivet analyses reference
ATLAS_2022_I2614196
Zy+jets at 13 TeV
Experiment: ATLAS (LHC)
Inspire ID: 2614196
Status: VALIDATED
Authors:
- Lorenzo Rossini
- Christian Gutschow
References:
Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
- pp -> Zy + jets at 13 TeV
Differential cross-section measurements of $Z\gamma$ production in association with hadronic jets are presented, using the full 139 fb${}^{-1}$ dataset of $\sqrt{s}$=13 TeV proton-proton collisions collected by the ATLAS detector during Run 2 of the LHC. Distributions are measured using events in which the $Z$ boson decays leptonically and the photon is usually radiated from an initial-state quark. Measurements are made in both one and two observables, including those sensitive to the hard scattering in the event and others which probe additional soft and collinear radiation. Different Standard Model predictions, from both parton-shower Monte Carlo simulation and fixed-order QCD calculations, are compared with the measurements. In general, good agreement is observed between data and predictions from MATRIX and MiNNLOPS, as well as next-to-leading-order predictions from MadGraph5_aMC@NLO and Sherpa.
Source code:
ATLAS_2022_I2614196.cc
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254 | // -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Projections/VetoedFinalState.hh"
#include "Rivet/Projections/DressedLeptons.hh"
#include "Rivet/Projections/PromptFinalState.hh"
#include "Rivet/Projections/InvisibleFinalState.hh"
namespace Rivet {
/// @brief Zy+jets at 13 TeV
class ATLAS_2022_I2614196 : public Analysis {
public:
/// Constructor
RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2022_I2614196);
/// @name Analysis methods
//@{
/// Book histograms and initialise projections before the run
void init() {
// Prompt photons
const PromptFinalState photon_fs(Cuts::abspid == PID::PHOTON && Cuts::pT > 30*GeV && Cuts::abseta < 2.37);
declare(photon_fs, "Photons");
// Prompt leptons
const PromptFinalState bareelectron_fs(Cuts::abspid == PID::ELECTRON);
const PromptFinalState baremuon_fs(Cuts::abspid == PID::MUON);
// Dressed leptons
const FinalState allphoton_fs(Cuts::abspid == PID::PHOTON);
const Cut leptoncut = Cuts::pT > 25*GeV && Cuts::abseta < 2.47;
const DressedLeptons dressedelectron_fs(allphoton_fs, bareelectron_fs, 0.1, leptoncut, true); // use *all* photons for lepton dressing
const DressedLeptons dressedmuon_fs(allphoton_fs, baremuon_fs, 0.1, leptoncut, true); // use *all* photons for lepton dressing
declare(dressedelectron_fs, "Electrons");
declare(dressedmuon_fs, "Muons");
// FS excluding the leading photon
VetoedFinalState vfs;
vfs.addVetoOnThisFinalState(photon_fs);
vfs.addVetoOnThisFinalState(dressedmuon_fs);
vfs.addVetoOnThisFinalState(InvisibleFinalState());
declare(vfs, "isolatedFS");
VetoedFinalState hadrons(FinalState(Cuts::abseta < 4.4));
hadrons.addVetoOnThisFinalState(dressedelectron_fs);
hadrons.addVetoOnThisFinalState(dressedmuon_fs);
declare(hadrons, "hadrons");
//Jets
FastJets jets(hadrons, FastJets::ANTIKT, 0.4, JetAlg::Muons::ALL, JetAlg::Invisibles::DECAY);
declare(jets, "jets");
// Histograms
book(_h["pTll"], 1, 1, 1);
book(_h["DiffpTll_pTy"], 2, 1, 1);
book(_h["SumpTll_pTy"], 3, 1, 1);
book(_h["DeltaRll"], 4, 1, 1);
book(_h["NJetsMix"], 5, 1, 1);
book(_h["pTJet1"], 6, 1, 1);
book(_h["pTJet2"], 7, 1, 1);
book(_h["RatiopTJet12"], 8, 1, 1);
book(_h["mjj"], 9, 1, 1);
book(_h["mllyj"], 10, 1, 1);
book(_h["HT"], 11, 1, 1);
book(_h["pTysqrtHT"], 12, 1, 1);
book(_h["DeltaPhiJetY"], 13, 1, 1);
book(_h["pTllyj"], 14, 1, 1);
// The binning on HepData is not directly usable for the following two.
// (Moreover, the bin edges have insufficient numerical precision.)
book(_h["phi_CS_Slice1"], 51, 1, 1);
book(_h["phi_CS_Slice2"], 51, 1, 2);
book(_h["phi_CS_Slice3"], 51, 1, 3);
book(_h["phi_CS_Slice4"], 51, 1, 4);
book(_h["phi_CS_Slice5"], 51, 1, 5);
book(_h["ctheta_CS_Slice1"], 52, 1, 1);
book(_h["ctheta_CS_Slice2"], 52, 1, 2);
book(_h["ctheta_CS_Slice3"], 52, 1, 3);
book(_h["ctheta_CS_Slice4"], 52, 1, 4);
book(_h["ctheta_CS_Slice5"], 52, 1, 5);
book(_h["RatiopTlly_mlly_Slice1"], 17, 1, 1);
book(_h["RatiopTlly_mlly_Slice2"], 18, 1, 1);
book(_h["RatiopTlly_mlly_Slice3"], 19, 1, 1);
book(_h["DiffpTll_pTy_Slice1"], 20, 1, 1);
book(_h["DiffpTll_pTy_Slice2"], 21, 1, 1);
book(_h["DiffpTll_pTy_Slice3"], 22, 1, 1);
book(_h["pTllyj_Slice1"], 23, 1, 1);
book(_h["pTllyj_Slice2"], 24, 1, 1);
book(_h["pTllyj_Slice3"], 25, 1, 1);
}
/// Perform the per-event analysis
void analyze(const Event& event) {
// const double weight = event.weight();
// Get objects
Particles electrons = apply<DressedLeptons>(event, "Electrons").particlesByPt();
Particles muons = apply<DressedLeptons>(event, "Muons").particlesByPt();
Particles photons = apply<PromptFinalState>(event, "Photons").particlesByPt();
if (photons.empty()) vetoEvent;
if (electrons.size()<2 && muons.size()<2) vetoEvent;
Particles lep;
// Sort the dressed leptons by pt
if (electrons.size() >= 2) {
lep.push_back(electrons[0]);
lep.push_back(electrons[1]);
}
else {
lep.push_back(muons[0]);
lep.push_back(muons[1]);
}
if (lep[0].momentum().pT() < 30*GeV) vetoEvent;
const double mll = (lep[0].momentum() + lep[1].momentum()).mass();
if (mll < 40*GeV) vetoEvent;
Particles selectedPh;
Particles fs = apply<VetoedFinalState>(event, "isolatedFS").particles();
for (const Particle& ph : photons) {
// check photon isolation
double coneEnergy(0.0);
for (const Particle& p : fs) {
if ( deltaR(ph, p) < 0.2 ) coneEnergy += p.Et();
}
if (coneEnergy / ph.momentum().pT() > 0.07 ) continue;
if (any(electrons, deltaRLess(ph, 0.4))) continue;
if (any(muons, deltaRLess(ph, 0.4))) continue;
selectedPh.push_back(ph);
}
if(selectedPh.size()<1) vetoEvent;
const double mlly = (lep[0].momentum() + lep[1].momentum() + selectedPh[0].momentum()).mass();
if (mll + mlly <= 182*GeV) vetoEvent;
// Get jets
// abseta
const Cut jetscut = (Cuts::pT > 30*GeV && Cuts::absrap < 2.5) || (Cuts::pT > 50*GeV && Cuts::absrapIn(2.5, 4.5));
const Cut jetscut25 = (Cuts::pT > 25*GeV && Cuts::abseta < 2.5) || (Cuts::pT > 50*GeV && Cuts::absetaIn(2.5, 4.4));
Jets jetsMix = apply<FastJets>(event, "jets").jetsByPt(jetscut);
Jets jetsMix25 = apply<FastJets>(event, "jets").jetsByPt(jetscut25);
iselectIfAnyDeltaRLess(jetsMix, photons, 0.4);
iselectIfAnyDeltaRLess(jetsMix, electrons, 0.2);
iselectIfAnyDeltaRLess(jetsMix25, photons, 0.4);
iselectIfAnyDeltaRLess(jetsMix25, electrons, 0.2);
// Jet multiplicity
size_t njets = jetsMix.size();
_h["NJetsMix"]->fill(njets);
double pTJet1 = njets>0 ? jetsMix[0].pT(): -1.;
double pTJet2 = njets>1 ? jetsMix[1].pT(): -1.;
double mjj = njets>1 ? (jetsMix[0].mom() + jetsMix[1].mom()).mass() : -1.;
double mllyj = njets>0 ? (lep[0].mom() + lep[1].mom() + selectedPh[0].mom() + jetsMix[0].mom()).mass() : -1.;
double HT = lep[0].pT() + lep[1].pT() + selectedPh[0].pT();
HT = sum(jetsMix25, Kin::pT, HT);
double pTy = selectedPh[0].pT();
double DeltaPhijy = njets>0 ? fabs(deltaPhi(jetsMix[0], selectedPh[0])) : -1.;
double DeltaRll = fabs(deltaR(lep[0], lep[1]));
double pTll = (lep[0].momentum() + lep[1].momentum()).pT();
double pTlly = (lep[0].mom() + lep[1].mom() + selectedPh[0].mom()).pT();
double Ratio_pTlly_mlly = pTlly/mlly;
double Diff_pTll_pTy = pTll-pTy;
double Sum_pTll_pTy = pTll+pTy;
double pTllyj = njets>0 ? (lep[0].mom() + lep[1].mom() + selectedPh[0].mom() + jetsMix[0].mom()).pT() : -1.;
// Fill plots
_h["pTJet1"]->fill(pTJet1/GeV);
_h["pTJet2"]->fill(pTJet2/GeV);
_h["RatiopTJet12"]->fill(pTJet2/pTJet1);
_h["mllyj"]->fill(mllyj/GeV);
_h["mjj"]->fill(mjj/GeV);
_h["HT"]->fill(HT/GeV);
_h["pTysqrtHT"]->fill((pTy/GeV)/sqrt(HT/GeV));
_h["DeltaPhiJetY"]->fill(DeltaPhijy);
_h["DeltaRll"]->fill(DeltaRll);
_h["pTll"]->fill(pTll/GeV);
if (mlly > 125*GeV && mlly < 200*GeV) _h["RatiopTlly_mlly_Slice1"]->fill(Ratio_pTlly_mlly);
else if (mlly > 200*GeV && mlly < 300*GeV) _h["RatiopTlly_mlly_Slice2"]->fill(Ratio_pTlly_mlly);
else if (mlly > 300*GeV) _h["RatiopTlly_mlly_Slice3"]->fill(Ratio_pTlly_mlly);
_h["SumpTll_pTy"]->fill(Sum_pTll_pTy/GeV);
_h["DiffpTll_pTy"]->fill(Diff_pTll_pTy/GeV);
if (Sum_pTll_pTy < 200*GeV) _h["DiffpTll_pTy_Slice1"]->fill(Diff_pTll_pTy/GeV);
else if (Sum_pTll_pTy > 200*GeV && Sum_pTll_pTy < 300*GeV) _h["DiffpTll_pTy_Slice2"]->fill(Diff_pTll_pTy/GeV);
else if (Sum_pTll_pTy > 300*GeV) _h["DiffpTll_pTy_Slice3"]->fill(Diff_pTll_pTy/GeV);
_h["pTllyj"]->fill(pTllyj/GeV);
if (pTlly < 50*GeV) _h["pTllyj_Slice1"]->fill(pTllyj/GeV);
else if (pTlly > 50*GeV && pTlly < 75*GeV) _h["pTllyj_Slice2"]->fill(pTllyj/GeV);
else if (pTlly > 75*GeV) _h["pTllyj_Slice3"]->fill(pTllyj/GeV);
// polarization variables
FourMomentum lep1, lep2;
if (lep[0].pid() > 0 && lep[1].pid() < 0) {
lep1 = lep[0].mom();
lep2 = lep[1].mom();
}
else if (lep[0].pid() < 0 && lep[1].pid() > 0) {
lep1 = lep[1].mom();
lep2 = lep[0].mom();
}
else {
MSG_DEBUG("Same sign lepton!");
vetoEvent;
}
FourMomentum Dilepton = (lep1+lep2);
FourMomentum lep1_CS;
LorentzTransform boost;
boost.setBetaVec(-Dilepton.betaVec());
lep1_CS = boost.transform(lep1);
double ctheta_CS, phi_CS;
ctheta_CS = cos(lep1_CS.p3().theta());
phi_CS = lep1_CS.p3().azimuthalAngle();
vector<double> pTll_binning = {0, 35, 60, 90, 135, 2500};
int pTll_binnum = std::lower_bound(pTll_binning.begin(), pTll_binning.end(), pTll) - pTll_binning.begin();
if (mll > 80*GeV && mll<100*GeV) {
string s = "_Slice" + to_string(pTll_binnum);
_h["ctheta_CS" + s]->fill(ctheta_CS);
_h["phi_CS" + s]->fill(phi_CS);
}
}
// end of analysis
/// Normalise histograms etc., after the run
void finalize() {
const double sf = crossSection()/femtobarn/sumOfWeights();
scale(_h, sf);
}
private:
/// Histograms
map<string, Histo1DPtr> _h;
};
// The hook for the plugin system
RIVET_DECLARE_PLUGIN(ATLAS_2022_I2614196);
}
|
|