Rivet analyses referenceATLAS_2021_I1849535Inclusive 4-lepton cross sections at 13 TeVExperiment: ATLAS (LHC) Inspire ID: 1849535 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
Measurements of four-lepton differential and integrated fiducial cross-sections in events with two same-flavour, opposite-charge electron or muon pairs are presented. The data correspond to 139 fb$^{-1}$ of $\sqrt{s}=13$ TeV proton-proton collisions, collected by the ATLAS detector during Run 2 of the Large Hadron Collider (2015--2018). The final state has contributions from a number of interesting Standard Model processes that dominate in different four-lepton invariant mass regions, including single $Z$ boson production, Higgs boson production and on-shell $ZZ$ production, with a complex mix of interference terms, and possible contributions from physics beyond the Standard Model. The differential cross-sections include the four-lepton invariant mass inclusively, in slices of other kinematic variables, and in different lepton flavour categories. Also measured are dilepton invariant masses, transverse momenta, and angular correlation variables, in four regions of four-lepton invariant mass, each dominated by different processes. The measurements are corrected for detector effects and are compared with state-of-the-art Standard Model calculations, which are found to be consistent with the data. The $Z\to 4\ell$ branching fraction is extracted, giving a value of $\left(4.41 \pm 0.30\right) \times 10^{-6}$. Constraints on effective field theory parameters and a model based on a spontaneously broken $B-L$ gauge symmetry are also evaluated. Further reinterpretations can be performed with the provided information. Source code: ATLAS_2021_I1849535.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/PromptFinalState.hh"
5#include "Rivet/Projections/VetoedFinalState.hh"
6#include "Rivet/Projections/ChargedFinalState.hh"
7#include "Rivet/Projections/LeptonFinder.hh"
8
9namespace Rivet {
10
11
12 /// M4l lineshape analysis
13 class ATLAS_2021_I1849535 : public Analysis {
14 public:
15
16 /// Constructor
17 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2021_I1849535);
18
19 void init() {
20
21 // Selection
22 Cut el_fid_sel = (Cuts::abseta < 2.47) && (Cuts::pT > 7*GeV);
23 Cut mu_fid_sel = (Cuts::abseta < 2.7) && (Cuts::pT > 5*GeV);
24
25 PromptFinalState photons(Cuts::abspid == PID::PHOTON);
26 PromptFinalState elecs(Cuts::abspid == PID::ELECTRON);
27 PromptFinalState muons(Cuts::abspid == PID::MUON && mu_fid_sel);
28 elecs.acceptTauDecays(true);
29 muons.acceptTauDecays(true);
30
31
32 // Final state including all charged particles
33 declare(ChargedFinalState(), "CFS");
34
35 LeptonFinder dressed_elecs(elecs, photons, 0.1, el_fid_sel);
36 declare(dressed_elecs, "elecs");
37 declare(muons, "muons");
38
39 // Book histos
40 book(_h["m4l_paper"], 1,1,1);
41 book(_h["m4l_4mu_paper"], 2,1,1);
42 book(_h["m4l_4e_paper"], 3,1,1);
43 book(_h["m4l_2e2mu_paper"],4,1,1);
44
45 book(_h["mZ1_Z_paper"],5,1,1);
46 book(_h["mZ1_H_paper"],6,1,1);
47 book(_h["mZ1_offshell_paper"],7,1,1);
48 book(_h["mZ1_ZZ_paper"],8,1,1);
49
50 book(_h["mZ2_Z_paper"],9,1,1);
51 book(_h["mZ2_H_paper"],10,1,1);
52 book(_h["mZ2_offshell_paper"],11,1,1);
53 book(_h["mZ2_ZZ_paper"],12,1,1);
54
55 book(_h["ptZ1_Z_paper"],13,1,1);
56 book(_h["ptZ1_H_paper"],14,1,1);
57 book(_h["ptZ1_offshell_paper"],15,1,1);
58 book(_h["ptZ1_ZZ_paper"],16,1,1);
59
60 book(_h["ptZ2_Z_paper"],17,1,1);
61 book(_h["ptZ2_H_paper"],18,1,1);
62 book(_h["ptZ2_offshell_paper"],19,1,1);
63 book(_h["ptZ2_ZZ_paper"],20,1,1);
64
65 book(_h["costhetastar1_Z_paper"],21,1,1);
66 book(_h["costhetastar1_H_paper"],22,1,1);
67 book(_h["costhetastar1_offshell_paper"],23,1,1);
68 book(_h["costhetastar1_ZZ_paper"],24,1,1);
69
70 book(_h["costhetastar2_Z_paper"],25,1,1);
71 book(_h["costhetastar2_H_paper"],26,1,1);
72 book(_h["costhetastar2_offshell_paper"],27,1,1);
73 book(_h["costhetastar2_ZZ_paper"],28,1,1);
74
75 book(_h["dy_Z1Z2_Z_paper"] ,29,1,1);
76 book(_h["dy_Z1Z2_H_paper"] ,30,1,1);
77 book(_h["dy_Z1Z2_offshell_paper"],31,1,1);
78 book(_h["dy_Z1Z2_ZZ_paper"] ,32,1,1);
79
80 book(_h["dphi_Z1Z2_Z_paper"] ,33,1,1);
81 book(_h["dphi_Z1Z2_H_paper"] ,34,1,1);
82 book(_h["dphi_Z1Z2_offshell_paper"],35,1,1);
83 book(_h["dphi_Z1Z2_ZZ_paper"] ,36,1,1);
84
85 book(_h["dphi_l1l2_Z_paper"] ,37,1,1);
86 book(_h["dphi_l1l2_H_paper"] ,38,1,1);
87 book(_h["dphi_l1l2_offshell_paper"],39,1,1);
88 book(_h["dphi_l1l2_ZZ_paper"] ,40,1,1);
89
90 book(_h["m4l_ptslice1_paper"],41,1,1);
91 book(_h["m4l_ptslice2_paper"],42,1,1);
92 book(_h["m4l_ptslice3_paper"],43,1,1);
93 book(_h["m4l_ptslice4_paper"],44,1,1);
94 book(_h["m4l_ptslice5_paper"],45,1,1);
95
96 book(_h["m4l_yslice1_paper"],46,1,1);
97 book(_h["m4l_yslice2_paper"],47,1,1);
98 book(_h["m4l_yslice3_paper"],48,1,1);
99 book(_h["m4l_yslice4_paper"],49,1,1);
100 book(_h["m4l_yslice5_paper"],50,1,1);
101
102 }
103
104
105 /// Generic dilepton candidate
106 struct Dilepton : public ParticlePair {
107 Dilepton() { }
108 Dilepton(ParticlePair _particlepair) : ParticlePair(_particlepair) {
109 assert(first.abspid() == second.abspid());
110 }
111 FourMomentum mom() const { return first.momentum() + second.momentum(); }
112 operator FourMomentum() const { return mom(); }
113 static bool cmppT(const Dilepton& lx, const Dilepton& rx) { return lx.mom().pT() < rx.mom().pT(); }
114 int flavour() const { return first.abspid(); }
115 double pTl1() const { return first.pT(); }
116 double pTl2() const { return second.pT(); }
117 };
118
119
120 struct Quadruplet {
121 Quadruplet (Dilepton z1, Dilepton z2): _z1(z1), _z2(z2) { }
122 enum class FlavCombi { mm=0, ee, me, em, undefined };
123 FourMomentum mom() const { return _z1.mom() + _z2.mom(); }
124 Dilepton getZ1() const { return _z1; }
125 Dilepton getZ2() const { return _z2; }
126 Dilepton _z1, _z2;
127 FlavCombi type() const {
128 if ( _z1.flavour() == 13 && _z2.flavour() == 13) { return FlavCombi::mm; }
129 else if (_z1.flavour() == 11 && _z2.flavour() == 11) { return FlavCombi::ee; }
130 else if (_z1.flavour() == 13 && _z2.flavour() == 11) { return FlavCombi::me; }
131 else if (_z1.flavour() == 11 && _z2.flavour() == 13) { return FlavCombi::em; }
132 else return FlavCombi::undefined;
133 }
134 };
135
136
137 bool passesTruthIsolation(Quadruplet quad, const Particles charged_tracks, Particles& truthLeptons ){
138 bool pass =true;
139 Particles leps;
140 leps.push_back(quad._z1.first);
141 leps.push_back(quad._z2.first);
142 leps.push_back(quad._z1.second);
143 leps.push_back(quad._z2.second);
144 for (auto &lep : leps){
145 double pTinCone = -lep.pT();
146 for (const Particle& track : charged_tracks) {
147 if (deltaR(lep.momentum(), track.momentum()) < 0.3)
148 pTinCone += track.pT();
149 }
150 for (const Particle& tlep: truthLeptons) {
151 float dR= deltaR(lep.momentum(), tlep.momentum());
152 if ( dR>0 && dR < 0.3)
153 pTinCone -= tlep.pT();
154 }
155 if (pTinCone > 0.16* lep.pT()){
156 pass=false;
157 }
158 }
159 return pass;
160 }
161
162 vector<Quadruplet> getBestQuads(Particles& particles, bool drcut = true) {
163 // H->ZZ->4l pairing
164 // - Two same flavor opposite charged leptons
165 // - Ambiguities in pairing are resolved by choosing the combination
166 // that results in the smaller value of |mll - mZ| for each pair successively
167 vector<Quadruplet> quads {};
168
169 size_t n_parts = particles.size();
170 if (n_parts < 4) return quads;
171
172 // STEP 1: find SFOS pairs
173 vector<Dilepton> SFOS;
174 for (size_t i = 0; i < n_parts; ++i) {
175 for (size_t j = 0; j < i; ++j) {
176 if (particles[i].pid() == -particles[j].pid()) {
177 // sort such that the negative lepton is listed first
178 Dilepton sfos;
179 if (particles[i].pid() > 0) sfos = make_pair(particles[i], particles[j]);
180 else sfos = make_pair(particles[j], particles[i]);
181
182 if (sfos.mom().mass() > _ll_mass*GeV && (!drcut || deltaR(particles[i],particles[j]) > _dRll)) SFOS.push_back(sfos);
183 }
184 }
185 }
186 if (SFOS.size() < 2) return quads;
187
188 // now we sort the SFOS pairs
189 std::sort(SFOS.begin(), SFOS.end(), [](const Dilepton& p1, const Dilepton& p2) {
190 return fabs(p1.mom().mass() - Z_mass) < fabs(p2.mom().mass() - Z_mass);
191 });
192
193 //form all possible quadruplets, passing the pt cuts, the dR cuts and the mll cuts
194 for (size_t k = 0; k < SFOS.size(); ++k) {
195 for (size_t l = k+1; l < SFOS.size(); ++l) {
196 if(drcut) {
197 if (deltaR(SFOS[k].first.mom(), SFOS[l].first.mom()) < _dRll) continue;
198 if (deltaR(SFOS[k].first.mom(), SFOS[l].second.mom()) < _dRll) continue;
199 if (deltaR(SFOS[k].second.mom(), SFOS[l].first.mom()) < _dRll) continue;
200 if (deltaR(SFOS[k].second.mom(), SFOS[l].second.mom()) < _dRll) continue;
201 }
202 if ( (SFOS[k].first.pid() == -SFOS[l].first.pid()) && ((SFOS[k].first.mom() + SFOS[l].first.mom()).mass() < _ll_mass*GeV)) continue;
203 if ( (SFOS[k].first.pid() == -SFOS[l].second.pid()) && ((SFOS[k].first.mom() + SFOS[l].second.mom()).mass() < _ll_mass*GeV)) continue;
204 if ( (SFOS[k].second.pid() == -SFOS[l].first.pid()) && ((SFOS[k].second.mom() + SFOS[l].first.mom()).mass() < _ll_mass*GeV)) continue;
205 if ( (SFOS[k].second.pid() == -SFOS[l].second.pid()) && ((SFOS[k].second.mom() + SFOS[l].second.mom()).mass() < _ll_mass*GeV)) continue;
206
207 //think technically this should happen before quad formation now and with all leptons not just those in quad so commenting out
208 //std::vector<double> lep_pt { SFOS[k].pTl1(), SFOS[k].pTl2(), SFOS[l].pTl1(), SFOS[l].pTl2() };
209 //std::sort(lep_pt.begin(), lep_pt.end(), std::greater<double>());
210 //if (!(lep_pt[0] > _pt_lep1*GeV && lep_pt[1] > _pt_lep2*GeV && lep_pt[2] > _pt_lep3*GeV)) continue;
211 quads.push_back( Quadruplet(SFOS[k], SFOS[l]) );
212 }
213 }
214 return quads;
215 }
216
217
218 bool passPtLeptons(const Particles& particles) {
219 size_t n_parts = particles.size();
220 if (n_parts < 4) return false;
221 // cut on pT of leptons
222 return ( particles[0].mom().pt() > _pt_lep1*GeV && particles[1].mom().pt() > _pt_lep2*GeV && particles[2].mom().pt() > _pt_lep3*GeV) ;
223 }
224
225
226 // Do the analysis
227 void analyze(const Event& event) {
228
229 const Particles charged_tracks = apply<ChargedFinalState>(event, "CFS").particles();
230
231 // Preselection of leptons for ZZ-> llll final state
232 Particles dressed_leptons = apply<ParticleFinder>(event, "muons").particlesByPt() +
233 apply<ParticleFinder>(event, "elecs").particlesByPt();
234 isortByPt(dressed_leptons);
235
236 auto foundDressedNoDrll = getBestQuads(dressed_leptons,false);
237
238 //now doing pt cut before quad formation so also apply this here
239 if (!passPtLeptons(dressed_leptons)) vetoEvent;
240
241 auto foundDressed = getBestQuads(dressed_leptons);
242 // if we don't find any quad, we can stop here
243 if (foundDressed.empty()) vetoEvent;
244 if (!passesTruthIsolation(foundDressed[0], charged_tracks, dressed_leptons)) vetoEvent;
245
246 double m4l = foundDressed[0].mom().mass()/GeV;
247 double pt4l = foundDressed[0].mom().pT()/GeV;
248 double y4l = foundDressed[0].mom().absrap();
249 double mZ1 = foundDressed[0].getZ1().mom().mass()/GeV;
250 double mZ2 = foundDressed[0].getZ2().mom().mass()/GeV;
251 double ptZ1 = foundDressed[0].getZ1().mom().pT()/GeV;
252 double ptZ2 = foundDressed[0].getZ2().mom().pT()/GeV;
253 double dy_Z1Z2 = fabs(foundDressed[0].getZ1().mom().rapidity() - foundDressed[0].getZ2().mom().rapidity());
254 double dphi_Z1Z2 = deltaPhi(foundDressed[0].getZ1().mom(), foundDressed[0].getZ2().mom());
255 double dphi_l1l2 = deltaPhi(dressed_leptons[0].mom(),dressed_leptons[1].mom());
256
257 _h["m4l_paper"]->fill(m4l);
258 if ( pt4l < 10.) _h["m4l_ptslice1_paper"]->fill(m4l);
259 else if (pt4l < 20.) _h["m4l_ptslice2_paper"]->fill(m4l);
260 else if (pt4l < 50.) _h["m4l_ptslice3_paper"]->fill(m4l);
261 else if (pt4l < 100.) _h["m4l_ptslice4_paper"]->fill(m4l);
262 else if (pt4l < 600.) _h["m4l_ptslice5_paper"]->fill(m4l);
263
264 if (y4l < 0.3) _h["m4l_yslice1_paper"]->fill(m4l);
265 else if (y4l < 0.6) _h["m4l_yslice2_paper"]->fill(m4l);
266 else if (y4l < 0.9) _h["m4l_yslice3_paper"]->fill(m4l);
267 else if (y4l < 1.2) _h["m4l_yslice4_paper"]->fill(m4l);
268 else if (y4l < 2.5) _h["m4l_yslice5_paper"]->fill(m4l);
269
270
271 Quadruplet::FlavCombi flavour = foundDressed[0].type();
272 if ( flavour == Quadruplet::FlavCombi::mm) { _h["m4l_4mu_paper"]->fill(m4l); }
273 else if (flavour == Quadruplet::FlavCombi::ee) { _h["m4l_4e_paper"]->fill(m4l); }
274 else if (flavour == Quadruplet::FlavCombi::me || flavour == Quadruplet::FlavCombi::em) {
275 _h["m4l_2e2mu_paper"]->fill(m4l);
276 }
277
278 // polarization variables
279 // Get four-momentum of the first lepton pair
280 const FourMomentum pcom = foundDressed.at(0).getZ1().mom();
281 const Vector3 betacom = pcom.betaVec();
282 const Vector3 unitboostvec = betacom.unit();
283 const LorentzTransform comboost = LorentzTransform::mkFrameTransformFromBeta(betacom);
284 // Get four-momentum of the negative lepton w.r.t. the first lepton pair
285 const FourMomentum p1com = comboost.transform(foundDressed.at(0).getZ1().first.mom());
286 double costhetastar1 = cos(p1com.p3().angle(unitboostvec));
287
288 // Get four-momentum of the second lepton pair
289 const FourMomentum pcom2 = foundDressed.at(0).getZ2().mom();
290 const Vector3 betacom2 = pcom2.betaVec();
291 const Vector3 unitboostvec2 = betacom2.unit();
292 const LorentzTransform comboost2 = LorentzTransform::mkFrameTransformFromBeta(betacom2);
293 // Get four-momentum of the negative lepton w.r.t. the second lepton pair
294 const FourMomentum p2com = comboost2.transform(foundDressed.at(0).getZ2().first.mom());
295 double costhetastar2 = cos(p2com.p3().angle(unitboostvec2));
296
297 //fill m4l binned variables
298 if (60 < m4l && m4l < 100.) {
299 _h["mZ1_Z_paper"]->fill(mZ1);
300 _h["mZ2_Z_paper"]->fill(mZ2);
301 _h["ptZ1_Z_paper"]->fill(ptZ1);
302 _h["ptZ2_Z_paper"]->fill(ptZ2);
303 _h["dy_Z1Z2_Z_paper"]->fill(dy_Z1Z2);
304 _h["dphi_Z1Z2_Z_paper"]->fill(dphi_Z1Z2);
305 _h["dphi_l1l2_Z_paper"]->fill(dphi_l1l2);
306 _h["costhetastar1_Z_paper"]->fill(costhetastar1 );
307 _h["costhetastar2_Z_paper"]->fill(costhetastar2 );
308 }
309 else if(120 < m4l && m4l < 130 ){
310 _h["mZ1_H_paper"]->fill(mZ1);
311 _h["mZ2_H_paper"]->fill(mZ2);
312 _h["ptZ1_H_paper"]->fill(ptZ1);
313 _h["ptZ2_H_paper"]->fill(ptZ2);
314 _h["dy_Z1Z2_H_paper"]->fill(dy_Z1Z2);
315 _h["dphi_Z1Z2_H_paper"]->fill(dphi_Z1Z2);
316 _h["dphi_l1l2_H_paper"]->fill(dphi_l1l2);
317 _h["costhetastar1_H_paper"]->fill(costhetastar1 );
318 _h["costhetastar2_H_paper"]->fill(costhetastar2 );
319 }
320 else if(180 < m4l && m4l < 2000){
321 _h["mZ1_ZZ_paper"]->fill(mZ1);
322 _h["mZ2_ZZ_paper"]->fill(mZ2);
323 _h["ptZ1_ZZ_paper"]->fill(ptZ1);
324 _h["ptZ2_ZZ_paper"]->fill(ptZ2);
325 _h["dy_Z1Z2_ZZ_paper"]->fill(dy_Z1Z2);
326 _h["dphi_Z1Z2_ZZ_paper"]->fill(dphi_Z1Z2);
327 _h["dphi_l1l2_ZZ_paper"]->fill(dphi_l1l2);
328 _h["costhetastar1_ZZ_paper"]->fill(costhetastar1 );
329 _h["costhetastar2_ZZ_paper"]->fill(costhetastar2 );
330 }
331 else{
332 _h["mZ1_offshell_paper"]->fill(mZ1);
333 _h["mZ2_offshell_paper"]->fill(mZ2);
334 _h["ptZ1_offshell_paper"]->fill(ptZ1);
335 _h["ptZ2_offshell_paper"]->fill(ptZ2);
336 _h["dy_Z1Z2_offshell_paper"]->fill(dy_Z1Z2);
337 _h["dphi_Z1Z2_offshell_paper"]->fill(dphi_Z1Z2);
338 _h["dphi_l1l2_offshell_paper"]->fill(dphi_l1l2);
339 _h["costhetastar1_offshell_paper"]->fill(costhetastar1);
340 _h["costhetastar2_offshell_paper"]->fill(costhetastar2);
341 }
342
343 }//end analysis
344
345 /// Finalize
346 void finalize() {
347 const double sf = crossSection() / femtobarn / sumOfWeights();
348 scale(_h, sf);
349 }
350
351 private:
352
353 map<string, Histo1DPtr> _h;
354 static constexpr double Z_mass = 91.1876;
355 static constexpr float _pt_lep1 = 20.;
356 static constexpr float _pt_lep2 = 10.;
357 static constexpr float _pt_lep3 = 0.;
358 static constexpr float _ll_mass = 5.;
359 static constexpr float _dRll = 0.05;
360
361 }; // end class ATLAS_2021_I1849535
362
363 RIVET_DECLARE_PLUGIN(ATLAS_2021_I1849535);
364} // end namespace rivet
|