Rivet analyses referenceATLAS_2023_I2690799Measurement of electroweak 4 lepton + 2 jet productionExperiment: ATLAS (LHC) Inspire ID: 2690799 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
Differential cross-sections are measured for the production of four charged leptons in association with two jets. These measurements are sensitive to final states in which the jets are produced via the strong interaction as well as to the purely-electroweak vector boson scattering process. The analysis is performed using proton-proton collision data collected by ATLAS at $\sqrt{s}$=13 TeV and with an integrated luminosity of 140 fb${}^{-1}$. The data are corrected for the effects of detector inefficiency and resolution and are compared to state-of-the-art Monte Carlo event generator predictions. The differential cross-sections are used to search for anomalous weak-boson self-interactions that are induced by dimension-six and dimension-eight operators in Standard Model effective field theory. Source code: ATLAS_2023_I2690799.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#include "Rivet/Projections/PromptFinalState.hh"
8
9namespace Rivet {
10
11
12 /// @brief Measurement of electroweak 4 lepton + 2 jet production
13 class ATLAS_2023_I2690799 : public Analysis {
14 public:
15
16 /// Constructor
17 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2023_I2690799);
18
19 /// A Z Boson structrue:
20 struct Z_can {
21 int lep_1, lep_2;
22 double dZ;
23 FourMomentum p;
24 };
25
26 /// Book histograms and initialise projections before the run
27 void init() {
28
29 // Initialise and register projections
30
31
32 // FinalState of prompt photons and bare muons and electrons in the event
33 PromptFinalState photons(Cuts::abspid == PID::PHOTON);
34 PromptFinalState electrons(Cuts::abspid == PID::ELECTRON);
35 PromptFinalState muons(Cuts::abspid == PID::MUON);
36
37 // Dress the prompt bare leptons with prompt photons within dR < 0.05,
38 // and apply some fiducial cuts on the dressed leptons
39 Cut cuts_el = (Cuts::pT > 7*GeV) && ( Cuts::abseta < 2.47);
40 Cut cuts_mu = (Cuts::pT > 5*GeV) && ( Cuts::abseta < 2.7);
41
42 LeptonFinder dressed_electrons(electrons, photons, 0.1, cuts_el);
43 declare(dressed_electrons, "DressedElectrons");
44
45 LeptonFinder dressed_muons(muons, photons, 0.1, cuts_mu);
46 declare(dressed_muons, "DressedMuons");
47
48 // We define AntiKt4TruthWZJets here, the leptons here are not used for next step selection.
49 // Muons
50 PromptFinalState bare_mu(Cuts::abspid == PID::MUON, TauDecaysAs::PROMPT); // true = use muons from prompt tau decays
51 LeptonFinder all_dressed_mu(bare_mu, photons, 0.1, Cuts::abseta < 2.5);
52 //
53 // Electrons
54 PromptFinalState bare_el(Cuts::abspid == PID::ELECTRON, TauDecaysAs::PROMPT); // true = use electrons from prompt tau decays
55 LeptonFinder all_dressed_el(bare_el, photons, 0.1, Cuts::abseta < 2.5);
56
57 VetoedFinalState vfs(FinalState(Cuts::abseta < 4.5));
58 // The final-state particles declared above are clustered using FastJet with
59 // the anti-kT algorithm and a jet-radius parameter 0.4
60 // AntiKt4TruthWZJets
61 vfs.addVetoOnThisFinalState(all_dressed_el);
62 vfs.addVetoOnThisFinalState(all_dressed_mu);
63 FastJets jetfs(vfs, JetAlg::ANTIKT, 0.4, JetMuons::ALL, JetInvisibles::DECAY);
64 declare(jetfs, "jets");
65
66 // Book histograms
67 book( _h["m_jj_SR"], 7, 1, 1);
68 book( _h["m_jj_CR"], 8, 1, 1);
69 book( _h["m_4l_SR"], 9, 1, 1);
70 book( _h["m_4l_CR"], 10, 1, 1);
71 book( _h["pt_4l_SR"], 11, 1, 1);
72 book( _h["pt_4l_CR"], 12, 1, 1);
73 book( _h["dphi_SR"], 13, 1, 1);
74 book( _h["dphi_CR"], 14, 1, 1);
75 book( _h["dY_SR"], 15, 1, 1);
76 book( _h["dY_CR"], 16, 1, 1);
77 book( _h["costs1_SR"], 17, 1, 1);
78 book( _h["costs1_CR"], 18, 1, 1);
79 book( _h["costs3_SR"], 19, 1, 1);
80 book( _h["costs3_CR"], 20, 1, 1);
81 book( _h["pt_jj_SR"], 21, 1, 1);
82 book( _h["pt_jj_CR"], 22, 1, 1);
83 book( _h["pt_4ljj_SR"], 23, 1, 1);
84 book( _h["pt_4ljj_CR"], 24, 1, 1);
85 book( _h["st_4ljj_SR"], 25, 1, 1);
86 book( _h["st_4ljj_CR"], 26, 1, 1);
87
88 }
89
90
91 /// Perform the per-event analysis
92 void analyze(const Event& event) {
93
94 // Retrieve dressed leptons, sorted by pT
95 DressedLeptons muons = apply<LeptonFinder>(event, "DressedMuons").dressedLeptons();
96 DressedLeptons elecs = apply<LeptonFinder>(event, "DressedElectrons").dressedLeptons();
97 DressedLeptons leps;
98 leps.insert( leps.end(), muons.begin(), muons.end() );
99 leps.insert( leps.end(), elecs.begin(), elecs.end() );
100
101 size_t n_leps = leps.size();
102 if (n_leps < 4) vetoEvent;
103
104 //leading and sub-leading leps pt
105 sort(leps.begin(), leps.end(), cmpMomByPt);
106 if (leps[0].pT() < 20*GeV || leps[1].pT() < 20*GeV) vetoEvent;
107 //lepton seperation
108 double dR_min = 9999;
109 for (size_t i = 0; i < n_leps - 1; ++i) {
110 for (size_t j = i + 1; j < n_leps; ++j) {
111 double dR = deltaR(leps[i], leps[j]);
112 if (dR < dR_min) dR_min = dR;
113 }
114 }
115 if (dR_min < 0.05) vetoEvent;
116
117
118 //ZZ pairing
119 size_t sfos_pair = 0;
120 double minm2l = 9999.;
121 vector<Z_can> z;
122 Particle *l1, *l3;
123 for (size_t i = 0; i < n_leps - 1; ++i) {
124 for (size_t j = i + 1; j < n_leps; ++j) {
125
126 if (leps[i].pid() != -1*leps[j].pid()) continue;
127 ++sfos_pair;
128 FourMomentum twoL = leps[i].momentum() + leps[j].momentum();
129 if (twoL.mass() < minm2l) minm2l = twoL.mass();
130
131 double deltaM = fabs( twoL.mass() - Zmass_PDG);
132 Z_can tmp_z;
133 tmp_z.dZ = deltaM;
134 tmp_z.p = twoL;
135 if (leps[i].pid() < 0) {
136 tmp_z.lep_1 = i;
137 tmp_z.lep_2 = j;
138 } else{
139 tmp_z.lep_1 = j;
140 tmp_z.lep_2 = i;
141 }
142 z.push_back(tmp_z);
143 }
144 }
145 if (minm2l < 5*GeV) vetoEvent; //All SFOS ll pairs should have m2l > 5 GeV.
146 if (sfos_pair < 2) vetoEvent;
147
148
149 std::sort(z.begin(), z.end(), [](const Z_can& z1, const Z_can& z2) { return z1.dZ < z2.dZ;});
150
151 int z1_id = -1;
152 int z2_id = -1;
153 FourMomentum fourL;
154 bool found_4l = false;
155 for (size_t i = 0; i < z.size() - 1; ++i) {
156 for (size_t j = i + 1; j < z.size(); ++j) {
157 if (z[i].lep_1 == z[j].lep_1 ||
158 z[i].lep_1 == z[j].lep_2 ||
159 z[i].lep_2 == z[j].lep_1 ||
160 z[i].lep_2 == z[j].lep_2) continue;
161 fourL = z[i].p + z[j].p;
162 if (fourL.mass() > 130*GeV) {
163 z1_id = i;
164 z2_id = j;
165 found_4l =true;
166 break;
167 }
168 }
169 if (found_4l) break;
170 }
171 if (!found_4l) vetoEvent;
172
173 // To judge which is the leading pair
174
175 FourMomentum twoL_1, twoL_2;
176
177 if (abs(z[z1_id].p.rap()) > abs(z[z2_id].p.rap())){
178 twoL_1 = z[z1_id].p;
179 twoL_2 = z[z2_id].p;
180 int i = z[z1_id].lep_1;
181 int j = z[z2_id].lep_1;
182 l1 = &leps[i];
183 l3 = &leps[j];
184 } else {
185 twoL_1 = z[z2_id].p;
186 twoL_2 = z[z1_id].p;
187 int i = z[z2_id].lep_1;
188 int j = z[z1_id].lep_1;
189 l1 = &leps[i];
190 l3 = &leps[j];
191 }
192
193 // Retrieve clustered jets, sorted by pT, with a minimum pT cut
194 Jets jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 30*GeV && Cuts::absrap < 4.5);
195 if (jets.size() < 2) vetoEvent;
196 if (jets[0].pt() < 40*GeV) vetoEvent;
197 FourMomentum twoJ;
198 double dY = -1;
199 int jet_id1 = -1, jet_id2 = -1;
200 bool found_2j = false;
201 for (size_t i = 0; i < jets.size() - 1; ++i) {
202 for (size_t j = i + 1; j < jets.size(); ++j) {
203 if (jets[i].rap()*jets[j].rap() > 0) continue; // First pair in A/C
204 dY = fabs(jets[i].rap() - jets[j].rap());
205 twoJ = jets[i].momentum() + jets[j].momentum();
206 jet_id1 = i;
207 jet_id2 = j;
208 found_2j = true;
209 break;
210 }
211 if (found_2j) break;
212 }
213 if (!found_2j) vetoEvent;
214 if (dY < 2) vetoEvent;
215 if (twoJ.mass() < 300*GeV) vetoEvent;
216
217 // Calculate Observables!
218
219 double phi1, phi2;
220 if (jets[jet_id1].rap() > jets[jet_id2].rap()) {
221 phi1 = jets[jet_id1].phi(); phi2 = jets[jet_id2].phi();
222 }
223 else {
224 phi1 = jets[jet_id2].phi(); phi2 = jets[jet_id1].phi();
225 }
226 const double dphi = mapAngleMPiToPi(phi1 - phi2);
227 const double mjj = twoJ.mass();
228 const double m4l = fourL.mass();
229 const double pt4l = fourL.pT();
230 const double ptjj = twoJ.pT();
231 const double pt4ljj = (fourL + twoJ).pT();
232 const double st4ljj = twoL_1.pT() + twoL_2.pT() + jets[jet_id1].pT() + jets[jet_id2].pT();
233 const double centrality = abs(fourL.rap() - 0.5*(jets[jet_id1].rap() + jets[jet_id2].rap())) \
234 / abs(jets[jet_id1].rap() - jets[jet_id2].rap());
235
236 //negetive leps1 in Z rest frame work
237 Vector3 beta2lcom_1 = twoL_1.betaVec();
238 LorentzTransform com2lboost_1 = LorentzTransform::mkFrameTransformFromBeta(beta2lcom_1);
239 FourMomentum l1_2lcom = com2lboost_1.transform(l1->momentum());
240 const double costs1 = cos(twoL_1.angle(l1_2lcom));
241
242
243 //negetive leps3 in Z rest frame work
244 Vector3 beta2lcom_2 = twoL_2.betaVec();
245 LorentzTransform com2lboost_2 = LorentzTransform::mkFrameTransformFromBeta(beta2lcom_2);
246 FourMomentum l3_2lcom = com2lboost_2.transform(l3->momentum());
247 const double costs3 = cos(twoL_2.angle(l3_2lcom));
248
249 const string suff = (centrality < 0.4)? "_SR" : "_CR";
250 _h["dphi"+suff]->fill(dphi);
251 _h["m_jj"+suff]->fill(mjj/GeV);
252 _h["dY"+suff]->fill(dY);
253 _h["m_4l"+suff]->fill(m4l/GeV);
254 _h["pt_4l"+suff]->fill(pt4l/GeV);
255 _h["pt_jj"+suff]->fill(ptjj/GeV);
256 _h["pt_4ljj"+suff]->fill(pt4ljj/GeV);
257 _h["st_4ljj"+suff]->fill(st4ljj/GeV);
258 _h["costs1"+suff]->fill(costs1/GeV);
259 _h["costs3"+suff]->fill(costs3/GeV);
260 }
261
262
263 /// Normalise histograms etc., after the run
264 void finalize() {
265
266 const double norm = crossSection()/sumOfWeights()/femtobarn;
267 scale(_h, norm);
268
269 }
270
271 private:
272
273 map<string, Histo1DPtr> _h;
274 const double Zmass_PDG = 91.1876;
275 };
276
277
278 RIVET_DECLARE_PLUGIN(ATLAS_2023_I2690799);
279
280}
|