Rivet analyses referenceATLAS_2019_I1720442Inclusive 4-lepton lineshape at 13 TeVExperiment: ATLAS (LHC) Inspire ID: 1720442 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
A measurement of the four-lepton invariant mass spectrum is made with the ATLAS detector, using an integrated luminosity of 36.1 fb of proton-proton collisions at TeV delivered by the Large Hadron Collider. The differential cross-section is measured for events containing two same-flavour opposite-sign lepton pairs. It exhibits a rich structure, with different mass regions dominated in the Standard Model by single boson production, Higgs boson production, and boson pair production, and non-negligible interference effects at high invariant masses. The measurement is compared with state-of-the-art Standard Model calculations, which are found to be consistent with the data. These calculations are used to interpret the data in terms of and subprocesses, and to place constraints on a possible contribution from physics beyond the Standard Model. Source code: ATLAS_2019_I1720442.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/PromptFinalState.hh"
5#include "Rivet/Projections/LeptonFinder.hh"
6
7namespace Rivet {
8
9
10 /// ATLAS 4-lepton lineshape at 13 TeV
11 class ATLAS_2019_I1720442 : public Analysis {
12 public:
13
14 /// Constructor
15 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2019_I1720442);
16
17 void init() {
18
19 PromptFinalState photons(Cuts::abspid == PID::PHOTON);
20 PromptFinalState elecs(Cuts::abspid == PID::ELECTRON);
21 PromptFinalState muons(Cuts::abspid == PID::MUON);
22
23 // Selection
24 Cut el_fid_sel = (Cuts::abseta < 2.47) && (Cuts::pT > 7*GeV);
25 Cut mu_fid_sel = (Cuts::abseta < 2.7) && (Cuts::pT > 5*GeV);
26
27 LeptonFinder dressed_elecs(elecs, photons, 0.005, el_fid_sel);
28 declare(dressed_elecs, "elecs");
29
30 LeptonFinder dressed_muons(muons, photons, 0.005, mu_fid_sel);
31 declare(dressed_muons, "muons");
32
33 // Book histos
34 book(_h["m4l_inclusive"], 1,1,1);
35
36 book(_h["m4l_ptslice1"], 2,1,1);
37 book(_h["m4l_ptslice2"], 3,1,1);
38 book(_h["m4l_ptslice3"], 4,1,1);
39 book(_h["m4l_ptslice4"], 5,1,1);
40
41 book(_h["m4l_rapidityslice1"], 6,1,1);
42 book(_h["m4l_rapidityslice2"], 7,1,1);
43 book(_h["m4l_rapidityslice3"], 8,1,1);
44 book(_h["m4l_rapidityslice4"], 9,1,1);
45
46 book(_h["m4l_4mu"], 12,1,1);
47 book(_h["m4l_4e"], 13,1,1);
48 book(_h["m4l_2e2mu"], 14,1,1);
49 }
50
51
52 /// @brief Generic dilepton candidate
53 /// @todo Move into the Rivet core?
54 struct Dilepton : public ParticlePair {
55 Dilepton() { }
56 Dilepton(ParticlePair _particlepair) : ParticlePair(_particlepair) {
57 assert(first.abspid() == second.abspid());
58 }
59 FourMomentum mom() const { return first.momentum() + second.momentum(); }
60 operator FourMomentum() const { return mom(); }
61 static bool cmppT(const Dilepton& lx, const Dilepton& rx) { return lx.mom().pT() < rx.mom().pT(); }
62 int flavour() const { return first.abspid(); }
63 double pTl1() const { return first.pT(); }
64 double pTl2() const { return second.pT(); }
65 };
66
67
68 struct Quadruplet {
69 Quadruplet (Dilepton z1, Dilepton z2): _z1(z1), _z2(z2) { }
70 enum class FlavCombi { mm=0, ee, me, em, undefined };
71 FourMomentum mom() const { return _z1.mom() + _z2.mom(); }
72 Dilepton getZ1() const { return _z1; }
73 Dilepton getZ2() const { return _z2; }
74 Dilepton _z1, _z2;
75 FlavCombi type() const {
76 if ( _z1.flavour() == 13 && _z2.flavour() == 13) { return FlavCombi::mm; }
77 else if (_z1.flavour() == 11 && _z2.flavour() == 11) { return FlavCombi::ee; }
78 else if (_z1.flavour() == 13 && _z2.flavour() == 11) { return FlavCombi::me; }
79 else if (_z1.flavour() == 11 && _z2.flavour() == 13) { return FlavCombi::em; }
80 else return FlavCombi::undefined;
81 }
82 };
83
84
85 vector<Quadruplet> getBestQuads(Particles& particles) {
86 // H->ZZ->4l pairing
87 // - Two same flavor opposite charged leptons
88 // - Ambiguities in pairing are resolved by choosing the combination
89 // that results in the smaller value of |mll - mZ| for each pair successively
90 vector<Quadruplet> quads {};
91
92 size_t n_parts = particles.size();
93 if (n_parts < 4) return quads;
94
95 // STEP 1: find SFOS pairs
96 vector<Dilepton> SFOS;
97 for (size_t i = 0; i < n_parts; ++i) {
98 for (size_t j = 0; j < i; ++j) {
99 if (particles[i].pid() == -particles[j].pid()) {
100 // sort such that the negative lepton is listed first
101 if (particles[i].pid() > 0) SFOS.push_back(Dilepton(make_pair(particles[i], particles[j])));
102 else SFOS.push_back(Dilepton(make_pair(particles[j], particles[i])));
103 }
104 }
105 }
106 if (SFOS.size() < 2) return quads;
107
108 // now we sort the SFOS pairs
109 std::sort(SFOS.begin(), SFOS.end(), [](const Dilepton& p1, const Dilepton& p2) {
110 return fabs(p1.mom().mass() - Z_mass) < fabs(p2.mom().mass() - Z_mass);
111 });
112
113 // Form all possible quadruplets passing the pT cuts
114 for (size_t k = 0; k < SFOS.size(); ++k) {
115 for (size_t l = k+1; l < SFOS.size(); ++l) {
116 if (deltaR(SFOS[k].first.mom(), SFOS[l].first.mom()) < 1e-13) continue;
117 if (deltaR(SFOS[k].first.mom(), SFOS[l].second.mom()) < 1e-13) continue;
118 if (deltaR(SFOS[k].second.mom(), SFOS[l].first.mom()) < 1e-13) continue;
119 if (deltaR(SFOS[k].second.mom(), SFOS[l].second.mom()) < 1e-13) continue;
120
121 vector<double> lep_pt { SFOS[k].pTl1(), SFOS[k].pTl2(), SFOS[l].pTl1(), SFOS[l].pTl2() };
122 std::sort(lep_pt.begin(), lep_pt.end(), std::greater<double>());
123 if (!(lep_pt[0] > 20*GeV && lep_pt[1] > 15*GeV && lep_pt[2] > 10*GeV)) continue;
124 quads.push_back( Quadruplet(SFOS[k], SFOS[l]) );
125 }
126 }
127 return quads;
128 }
129
130
131 bool passMassCuts(const Quadruplet& theQuad){
132 const vector<double> cuts_m34{ 5*GeV, 5*GeV, 12*GeV, 12*GeV, 50*GeV };
133 const vector<double> cuts_m4l{ 0, 100*GeV, 110*GeV, 140*GeV, 190*GeV };
134
135 double m4l = theQuad.mom().mass();
136 double mZ1 = theQuad.getZ1().mom().mass();
137 double mZ2 = theQuad.getZ2().mom().mass();
138
139 // Invariant-mass requirements
140 double cutval = cuts_m34.back();
141 for (size_t k = 0; k < cuts_m34.size(); ++k) {
142 if (cuts_m4l[k] > m4l) {
143 cutval = cuts_m34[k-1] + (cuts_m34[k] - cuts_m34[k-1])/(cuts_m4l[k] - cuts_m4l[k-1]) * (m4l - cuts_m4l[k-1]);
144 break;
145 }
146 }
147 return inRange(mZ1, 50*GeV, 106*GeV) && inRange(mZ2, cutval, 115*GeV);
148 }
149
150
151 bool pass_dRll(const Quadruplet& theQuad) {
152 const double dR_min_same = 0.1;
153 const double dR_min_opp = 0.2;
154 double dr_min_cross = dR_min_opp;
155 if (theQuad.getZ1().flavour() == theQuad.getZ2().flavour()) {
156 dr_min_cross = dR_min_same;
157 }
158 return !((deltaR(theQuad.getZ1().first, theQuad.getZ1().second) < dR_min_same) ||
159 (deltaR(theQuad.getZ2().first, theQuad.getZ2().second) < dR_min_same) ||
160 (deltaR(theQuad.getZ1().first, theQuad.getZ2().first) < dr_min_cross) ||
161 (deltaR(theQuad.getZ1().first, theQuad.getZ2().second) < dr_min_cross) ||
162 (deltaR(theQuad.getZ1().second, theQuad.getZ2().first) < dr_min_cross) ||
163 (deltaR(theQuad.getZ1().second, theQuad.getZ2().second) < dr_min_cross));
164 }
165
166
167 bool pass_Jpsi(const Quadruplet& theQuad){
168 Particles all_leps { theQuad.getZ1().first, theQuad.getZ1().second, theQuad.getZ2().first, theQuad.getZ2().second };
169 for (const Particle& lep1 : all_leps) {
170 for (const Particle& lep2 : all_leps) {
171 if (lep1.pid() == -lep2.pid() && (lep1.mom() + lep2.mom()).mass() < 5*GeV) return false;
172 }
173 }
174 return true;
175 }
176
177
178 // Handle 3 further CF stages - m12/34, dRmin, jpsi veto
179 bool passSelection (const Quadruplet& theQuad){
180 return passMassCuts(theQuad) && pass_dRll(theQuad) && pass_Jpsi(theQuad);
181 }
182
183
184 // Do the analysis
185 void analyze(const Event& event) {
186
187 //preselection of leptons for ZZ-> llll final state
188 Particles dressed_leptons;
189 for (auto lep : apply<LeptonFinder>(event, "muons").dressedLeptons()) { dressed_leptons.push_back(lep); }
190 for (auto lep : apply<LeptonFinder>(event, "elecs").dressedLeptons()) { dressed_leptons.push_back(lep); }
191
192 auto foundDressed = getBestQuads(dressed_leptons);
193 // if we don't find any quad, we can stop here
194 if (foundDressed.empty()) vetoEvent;
195
196 bool pass = passSelection(foundDressed[0]);
197 if (pass) {
198 double m4l = foundDressed[0].mom().mass()/GeV;
199 double pt4l = foundDressed[0].mom().pT()/GeV;
200 double y4l = foundDressed[0].mom().absrap();
201 Quadruplet::FlavCombi flavour = foundDressed[0].type();
202 _h["m4l_inclusive"]->fill(m4l);
203 if ( pt4l < 20.) _h["m4l_ptslice1"]->fill(m4l);
204 else if (pt4l < 50.) _h["m4l_ptslice2"]->fill(m4l);
205 else if (pt4l < 100.) _h["m4l_ptslice3"]->fill(m4l);
206 else if (pt4l < 600.) _h["m4l_ptslice4"]->fill(m4l);
207
208 if ( y4l < 0.4) _h["m4l_rapidityslice1"]->fill(m4l);
209 else if (y4l < 0.8) _h["m4l_rapidityslice2"]->fill(m4l);
210 else if (y4l < 1.2) _h["m4l_rapidityslice3"]->fill(m4l);
211 else if (y4l < 2.5) _h["m4l_rapidityslice4"]->fill(m4l);
212
213 if ( flavour == Quadruplet::FlavCombi::mm) _h["m4l_4mu"]->fill(m4l);
214 else if (flavour == Quadruplet::FlavCombi::ee) _h["m4l_4e"]->fill(m4l);
215 else if (flavour == Quadruplet::FlavCombi::me || flavour == Quadruplet::FlavCombi::em) {
216 _h["m4l_2e2mu"]->fill(m4l);
217 }
218 }
219
220 }
221
222
223 /// Finalize
224 void finalize() {
225 const double sf = crossSection() / femtobarn / sumOfWeights();
226 scale(_h, sf);
227 }
228
229
230 private:
231
232 map<string, Histo1DPtr> _h;
233 static constexpr double Z_mass = 91.1876;
234
235 };
236
237
238 RIVET_DECLARE_PLUGIN(ATLAS_2019_I1720442);
239
240}
|