Rivet analyses referenceCMS_2021_I1978840Measurement of $\mathrm{W}^{\pm}\gamma$ differential cross sections in proton-proton collisions at $\sqrt{s}=13\,\mathrm{TeV}$ and effective field theory constraintsExperiment: CMS (LHC) Inspire ID: 1978840 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
Differential cross section measurements of $\mathrm{W}^{\pm}\gamma$ production in proton-proton collisions at $\sqrt{s} = 13\,\mathrm{TeV}$ are presented. The data set used in this study was collected with the CMS detector at the CERN LHC in 2016--2018 with an integrated luminosity of $138\,\mathrm{fb}^{-1}$. Candidate events containing an electron or muon, a photon, and missing transverse momentum are selected. The measurements are compared with standard model predictions computed at next-to-leading and next-to-next-to-leading orders in perturbative quantum chromodynamics. Constraints on the presence of TeV-scale new physics affecting the $\mathrm{WW}\gamma$ vertex are determined within an effective field theory framework, focusing on the $\mathcal{O}_{3W}$ operator. A simultaneous measurement of the photon transverse momentum and the azimuthal angle of the charged lepton in a special reference frame is performed. Source code: CMS_2021_I1978840.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Event.hh"
4#include "Rivet/Math/LorentzTrans.hh"
5#include "Rivet/Particle.hh"
6#include "Rivet/Projections/ChargedLeptons.hh"
7#include "Rivet/Projections/DressedLeptons.hh"
8#include "Rivet/Projections/FastJets.hh"
9#include "Rivet/Projections/FinalState.hh"
10#include "Rivet/Projections/IdentifiedFinalState.hh"
11#include "Rivet/Projections/MissingMomentum.hh"
12#include "Rivet/Projections/PromptFinalState.hh"
13#include "Rivet/Projections/VetoedFinalState.hh"
14#include "Rivet/Projections/InvisibleFinalState.hh"
15
16namespace Rivet {
17
18
19 /// @brief Measurement of W/gamma differential cross sections at sqrt(s) = 13 TeV
20 class CMS_2021_I1978840 : public Analysis {
21 public:
22
23 struct WGammaRivetVariables {
24 bool is_wg_gen;
25
26 double l0_pt;
27 double l0_eta;
28 double l0_phi;
29 double l0_M;
30 int l0_q;
31 unsigned l0_abs_pdgid;
32
33 double p0_pt;
34 double p0_eta;
35 double p0_phi;
36 double p0_M;
37 bool p0_frixione;
38 double p0_frixione_sum;
39
40 double l0p0_dr;
41 double mt_cluster;
42
43 double n0_pt;
44 double n0_eta;
45 double n0_phi;
46 double n0_M;
47
48 double met_pt;
49 double met_phi;
50
51 double true_phi;
52 double true_phi_f;
53
54 int n_jets;
55
56 WGammaRivetVariables() { resetVars(); }
57
58 void resetVars() {
59 is_wg_gen = false;
60 l0_pt = 0.;
61 l0_eta = 0.;
62 l0_phi = 0.;
63 l0_M = 0.;
64 l0_q = 0;
65 l0_abs_pdgid = 0;
66 p0_pt = 0.;
67 p0_eta = 0.;
68 p0_phi = 0.;
69 p0_M = 0.;
70 p0_frixione = false;
71 n0_pt = 0.;
72 n0_eta = 0.;
73 n0_phi = 0.;
74 n0_M = 0.;
75 met_pt = 0.;
76 met_phi = 0.;
77 true_phi = 0.;
78 true_phi_f = 0.;
79 l0p0_dr = 0.;
80 mt_cluster = 0.;
81 n_jets = 0;
82 }
83 };
84
85 struct WGSystem {
86 int lepton_charge;
87
88 FourMomentum wg_system;
89
90 FourMomentum c_w_boson;
91 FourMomentum c_charged_lepton;
92 FourMomentum c_neutrino;
93 FourMomentum c_photon;
94
95 FourMomentum r_w_boson;
96 FourMomentum r_charged_lepton;
97 FourMomentum r_neutrino;
98 FourMomentum r_photon;
99
100 WGSystem(Particle const& lep, FourMomentum const& neu, Particle const& pho, bool verbose);
101
102 double phi();
103 double symphi();
104 };
105
106 double photon_iso_dr_ = 0.4;
107 double lepton_pt_cut_ = 30.;
108 double lepton_abs_eta_cut_ = 2.5;
109 double photon_pt_cut_ = 30.;
110 double photon_abs_eta_cut_ = 2.5;
111 double missing_pt_cut_ = 40.;
112 double lepton_photon_dr_cut_ = 0.7;
113 double dressed_lepton_cone_ = 0.1;
114
115 double eft_lepton_pt_cut_ = 80.;
116 double eft_photon_pt_cut_ = 150.;
117
118 double jet_pt_cut_ = 30.;
119 double jet_abs_eta_cut_ = 2.5;
120 double jet_dr_cut_ = 0.4;
121
122 WGammaRivetVariables vars_;
123 map<string, Histo1DPtr> _h;
124
125
126 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2021_I1978840);
127
128
129 void init() {
130 vars_.resetVars();
131
132 FinalState fs;
133 declare(fs, "FinalState");
134
135 // Jets - all final state particles excluding neutrinos
136 VetoedFinalState vfs;
137 vfs.vetoNeutrinos();
138 FastJets fastjets(vfs, FastJets::ANTIKT, 0.4);
139 declare(fastjets, "Jets");
140
141 // Dressed leptons
142 ChargedLeptons charged_leptons(fs);
143 PromptFinalState prompt_leptons(charged_leptons);
144 declare(prompt_leptons, "PromptLeptons");
145
146 IdentifiedFinalState photons(fs);
147 photons.acceptIdPair(PID::PHOTON);
148 PromptFinalState prompt_photons(photons);
149 prompt_photons.acceptMuonDecays(true);
150 prompt_photons.acceptTauDecays(true);
151
152 DressedLeptons dressed_leptons(prompt_photons, prompt_leptons, dressed_lepton_cone_,
153 Cuts::open(),
154 /*useDecayPhotons*/ false);
155 declare(dressed_leptons, "DressedLeptons");
156
157 // Photons
158 VetoedFinalState vetoed_prompt_photons(prompt_photons);
159 vetoed_prompt_photons.addVetoOnThisFinalState(dressed_leptons);
160 declare(vetoed_prompt_photons, "Photons");
161
162 // Invisibles
163 InvisibleFinalState invisibles(true, true, true);
164 declare(invisibles, "Invisibles");
165
166 // MET
167 declare(MissingMomentum(fs), "MET");
168
169 // Booking of histograms
170 book(_h["baseline_photon_pt"], 1, 1, 1);
171 book(_h["baseline_photon_eta"], 8, 1, 1);
172 book(_h["baseline_leppho_dr"], 15, 1, 1);
173 book(_h["baseline_leppho_deta"], 22, 1, 1);
174 book(_h["baseline_mt_cluster"], 29, 1, 1);
175 book(_h["baseline_njet"], 36, 1, 1);
176 book(_h["raz_leppho_deta"], 40, 1, 1);
177 book(_h["eft_photon_pt_phi_0"], 54, 1, 1);
178 book(_h["eft_photon_pt_phi_1"], 55, 1, 1);
179 book(_h["eft_photon_pt_phi_2"], 56, 1, 1);
180 }
181
182 /// Perform the per-event analysis
183 void analyze(const Event& event) {
184 vars_.resetVars();
185
186 const Particles leptons = apply<FinalState>(event, "DressedLeptons").particlesByPt();
187
188 if (leptons.size() == 0) {
189 vetoEvent;
190 }
191 auto l0 = leptons.at(0);
192
193 const Particles photons = apply<FinalState>(event, "Photons")
194 .particlesByPt(DeltaRGtr(l0, lepton_photon_dr_cut_));
195 if (photons.size() == 0) {
196 vetoEvent;
197 }
198 auto p0 = photons.at(0);
199
200
201 const Particles invisibles = apply<FinalState>(event, "Invisibles").particlesByPt();
202 FourMomentum n0;
203 for (auto const& inv : invisibles) {
204 n0 += inv.momentum();
205 }
206
207 if (invisibles.size() == 0) {
208 vetoEvent;
209 }
210
211 FourMomentum met = apply<MissingMomentum>(event, "MET").missingMomentum();
212 // Redefine the MET
213 met = FourMomentum(met.pt(), met.px(), met.py(), 0.);
214
215 // Filter jets on pT, eta and DR with lepton and photon
216 const Jets jets = apply<FastJets>(event, "Jets").jetsByPt([&](Jet const& j) {
217 return j.pt() > jet_pt_cut_ && j.abseta() < jet_abs_eta_cut_ && deltaR(j, l0) > jet_dr_cut_ && deltaR(j, p0) > jet_dr_cut_;
218 });
219
220 if (leptons.size() >= 1 && photons.size() >= 1 && invisibles.size() >= 1) {
221 // Populate variables
222 vars_.is_wg_gen = true;
223 vars_.l0_pt = l0.pt()/GeV;
224 vars_.l0_eta = l0.eta();
225 vars_.l0_phi = l0.phi(PhiMapping::MINUSPI_PLUSPI);
226 vars_.l0_M = l0.mass()/GeV;
227 vars_.l0_q = l0.charge();
228 vars_.l0_abs_pdgid = l0.abspid();
229 vars_.p0_pt = p0.pt()/GeV;
230 vars_.p0_eta = p0.eta();
231 vars_.p0_phi = p0.phi(PhiMapping::MINUSPI_PLUSPI);
232 vars_.p0_M = p0.mass()/GeV;
233 vars_.n_jets = jets.size();
234 vars_.l0p0_dr = deltaR(l0, p0);
235 vars_.n0_pt = n0.pt()/GeV;
236 vars_.n0_eta = n0.eta();
237 vars_.n0_phi = n0.phi(PhiMapping::MINUSPI_PLUSPI);
238 vars_.n0_M = n0.mass()/GeV;
239 vars_.met_pt = met.pt()/GeV;
240 vars_.met_phi = met.phi(PhiMapping::MINUSPI_PLUSPI);
241
242 // Here we build a list of particles to cluster jets, to be used in the photon isolation.
243 // The selection of particles that we want to veto from this isolation sum is non-trivial:
244 // the leading pT lepton, the leading pT photon that has DeltaR > 0.7 from the leading pT
245 // lepton, the invisibles and any tau decay products. Therefore, the selection is done
246 // here instead of in the initialise method.
247 Particles finalparts_iso = apply<FinalState>(event, "FinalState").particles();
248 Particles filtered_iso;
249 for (Particle const& p : finalparts_iso) {
250 bool veto_particle = false;
251 if (p.genParticle() == l0.genParticle() || p.genParticle() == p0.genParticle() || p.fromTau()) {
252 veto_particle = true;
253 }
254 for (auto const& inv : invisibles) {
255 if (p.genParticle() == inv.genParticle()) {
256 veto_particle = true;
257 }
258 }
259 if (!veto_particle) {
260 filtered_iso.push_back(p);
261 }
262 }
263 auto proj = getProjection<FastJets>("Jets");
264 proj.reset();
265 proj.calc(filtered_iso);
266 auto jets_iso = proj.jets();
267
268 vars_.p0_frixione = true;
269 double frixione_sum = 0.;
270
271 // Apply Frixione isolation to the photon:
272 auto jparts = sortBy(jets_iso, [&](Jet const& part1, Jet const& part2) {
273 return deltaR(part1, p0) < deltaR(part2, p0);
274 });
275
276 for (auto const& ip : jparts) {
277 double dr = deltaR(ip, p0);
278 if (dr >= photon_iso_dr_) {
279 break;
280 }
281 frixione_sum += ip.pt();
282 if (frixione_sum > (p0.pt() * ((1. - cos(dr)) / (1. - cos(photon_iso_dr_))))) {
283 vars_.p0_frixione = false;
284 }
285 }
286
287 // Now calculate EFT phi observables
288 auto wg_system = WGSystem(l0, n0, p0, false);
289
290 vars_.true_phi = wg_system.phi();
291 vars_.true_phi_f = wg_system.symphi();
292
293 // Calculate mTcluster
294 auto cand1 = l0.momentum() + p0.momentum();
295 auto full_system = cand1 + met;
296 double mTcluster2 = sqr(sqrt(cand1.mass2() + cand1.pt2()) + met.pt()) - full_system.pt2();
297 if (mTcluster2 > 0) {
298 vars_.mt_cluster = sqrt(mTcluster2);
299 } else {
300 vars_.mt_cluster = 0.;
301 }
302
303 if (vars_.l0_pt > lepton_pt_cut_ && fabs(vars_.l0_eta) < lepton_abs_eta_cut_ &&
304 vars_.p0_pt > photon_pt_cut_ && fabs(vars_.p0_eta) < photon_abs_eta_cut_ &&
305 vars_.p0_frixione && vars_.l0p0_dr > lepton_photon_dr_cut_ &&
306 vars_.met_pt > missing_pt_cut_) {
307 _h["baseline_photon_pt"]->fill(vars_.p0_pt / GeV);
308 _h["baseline_photon_eta"]->fill(vars_.p0_eta);
309 _h["baseline_leppho_dr"]->fill(vars_.l0p0_dr);
310 _h["baseline_leppho_deta"]->fill(vars_.l0_eta - vars_.p0_eta);
311 _h["baseline_mt_cluster"]->fill(vars_.mt_cluster / GeV);
312 double fill_n_jets = vars_.n_jets >= 2 ? 2. : double(vars_.n_jets);
313 _h["baseline_njet"]->fill(fill_n_jets);
314 if (vars_.n_jets == 0 && vars_.mt_cluster > 150.) {
315 _h["raz_leppho_deta"]->fill(vars_.l0_eta - vars_.p0_eta);
316 }
317 }
318
319 if (vars_.l0_pt > eft_lepton_pt_cut_ && fabs(vars_.l0_eta) < lepton_abs_eta_cut_ &&
320 vars_.p0_pt > eft_photon_pt_cut_ && fabs(vars_.p0_eta) < photon_abs_eta_cut_ &&
321 vars_.p0_frixione && vars_.l0p0_dr > lepton_photon_dr_cut_ &&
322 vars_.met_pt > missing_pt_cut_ && vars_.n_jets == 0) {
323 double absphi = fabs(vars_.true_phi_f);
324 if (absphi > 0. && absphi <= (PI / 6.)) {
325 _h["eft_photon_pt_phi_0"]->fill(vars_.p0_pt / GeV);
326 } else if (absphi > (PI / 6.) && absphi <= (PI / 3.)) {
327 _h["eft_photon_pt_phi_1"]->fill(vars_.p0_pt / GeV);
328 } else if (absphi > (PI / 3.) && absphi <= (PI / 2.)) {
329 _h["eft_photon_pt_phi_2"]->fill(vars_.p0_pt / GeV);
330 }
331 }
332 }
333 }
334
335
336 void finalize() {
337 double flavor_factor = 3. / 2.; // account for the fact that tau events are vetoed
338 // Scale according to cross section
339 for (string x :
340 {"baseline_photon_pt", "baseline_photon_eta", "baseline_leppho_dr", "baseline_leppho_deta",
341 "baseline_mt_cluster", "baseline_njet", "raz_leppho_deta", "eft_photon_pt_phi_0",
342 "eft_photon_pt_phi_1", "eft_photon_pt_phi_2"}) {
343 scale(_h[x], flavor_factor * crossSection() / femtobarn / sumOfWeights());
344 }
345
346 // Since these are really 2D, we need to divide by the y bin width:
347 for (string x : {"eft_photon_pt_phi_0", "eft_photon_pt_phi_1", "eft_photon_pt_phi_2"}) {
348 scale(_h[x], 1. / (PI / 6.));
349 }
350 }
351
352 };
353
354
355 CMS_2021_I1978840::WGSystem::WGSystem(Particle const& lep, FourMomentum const& neu,
356 Particle const& pho, bool verbose) {
357 lepton_charge = lep.charge3() / 3;
358 wg_system += lep.momentum();
359 wg_system += neu;
360 wg_system += pho.momentum();
361
362 if (verbose) {
363 cout << "> charge: " << lepton_charge << "\n";
364 cout << "> wg_system: " << wg_system << "\n";
365 cout << "> lepton : " << lep.momentum() << "\n";
366 cout << "> neutrino : " << neu << "\n";
367 cout << "> photon : " << pho.momentum() << "\n";
368 }
369
370 auto boost = LorentzTransform::mkFrameTransformFromBeta(wg_system.betaVec());
371
372 c_charged_lepton = boost.transform(lep);
373 c_neutrino = boost.transform(neu);
374 c_photon = boost.transform(pho);
375
376 if (verbose) {
377 cout << "> c_lepton : " << c_charged_lepton << "\n";
378 cout << "> c_neutrino : " << c_neutrino << "\n";
379 cout << "> c_photon : " << c_photon << "\n";
380 }
381
382 FourMomentum c_w_boson;
383 c_w_boson += c_charged_lepton;
384 c_w_boson += c_neutrino;
385
386 Vector3 r_uvec = wg_system.vector3().unit();
387 if (verbose) {
388 cout << "> c_w_boson : " << c_w_boson << endl;
389 cout << "> r_uvec : " << r_uvec << endl;
390 }
391
392 Vector3 z_uvec = c_w_boson.vector3().unit();
393 Vector3 y_uvec = z_uvec.cross(r_uvec).unit();
394 Vector3 x_uvec = y_uvec.cross(z_uvec).unit();
395 if (verbose) {
396 cout << "> x_uvec : " << x_uvec << endl;
397 cout << "> y_uvec : " << y_uvec << endl;
398 cout << "> z_uvec : " << z_uvec << endl;
399 }
400
401 Matrix3 rot_matrix;
402 rot_matrix.setRow(0, x_uvec).setRow(1, y_uvec).setRow(2, z_uvec);
403 auto rotator = LorentzTransform();
404 rotator = rotator.postMult(rot_matrix);
405 if (verbose) {
406 cout << "> rotator : " << rotator << endl;
407 }
408 r_w_boson = rotator.transform(c_w_boson);
409 r_charged_lepton = rotator.transform(c_charged_lepton);
410 r_neutrino = rotator.transform(c_neutrino);
411 r_photon = rotator.transform(c_photon);
412
413 if (verbose) {
414 cout << "> r_lepton : " << r_charged_lepton << endl;
415 cout << "> r_neutrino : " << r_neutrino << endl;
416 cout << "> r_photon : " << r_photon << endl;
417 cout << "> r_w_boson : " << r_w_boson << endl;
418 }
419 }
420
421 double CMS_2021_I1978840::WGSystem::phi() {
422 double lep_phi = r_charged_lepton.phi(PhiMapping::MINUSPI_PLUSPI);
423 return mapAngleMPiToPi(lepton_charge > 0 ? (lep_phi) : (lep_phi + PI));
424 }
425
426 double CMS_2021_I1978840::WGSystem::symphi() {
427 double lep_phi = phi();
428 if (lep_phi > PI / 2.) {
429 return PI - lep_phi;
430 } else if (lep_phi < -1. * (PI / 2.)) {
431 return -1. * (PI + lep_phi);
432 } else {
433 return lep_phi;
434 }
435 }
436
437
438 RIVET_DECLARE_PLUGIN(CMS_2021_I1978840);
439
440}
|