Rivet analyses referenceCMS_2022_I2080534EW W+W- production at 13 TeVExperiment: CMS (LHC) Inspire ID: 2080534 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
An observation is reported of the electroweak production of a $W^+W^-$ pair in association with two jets, with both $W$ bosons decaying leptonically. The data sample corresponds to an integrated luminosity of 138fb$^{-1}$ of proton-proton collisions at $\sqrt{s} = 13$ TeV, collected by the CMS detector at the CERN LHC. Events are selected by requiring exactly two opposite-sign leptons (electrons or muons) and two jets with large pseudorapidity separation and high dijet invariant mass. Events are categorized based on the flavor of the final-state leptons. A signal is observed with a significance of 5.6 standard deviations (5.2 expected) with respect to the background-only hypothesis. The measured fiducial cross section is 10.2 $\pm$ 2.0 fb and this value is consistent with the standard model prediction of 9.1 $\pm$ 0.6 fb. Source code: CMS_2022_I2080534.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Projections/LeptonFinder.hh"
6#include "Rivet/Projections/MissingMomentum.hh"
7#include "Rivet/Projections/PromptFinalState.hh"
8
9namespace Rivet {
10
11 /// @brief EW W+W- production at 13 TeV
12 class CMS_2022_I2080534 : public Analysis {
13 public:
14
15 /// Constructor
16 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2022_I2080534);
17
18 /// @name Analysis methods
19 /// @{
20 void init() {
21
22 // The basic final-state projection:
23 // all final-state particles within
24 // the given eta acceptance
25 const FinalState fs(Cuts::abseta < 4.9 && Cuts::pT > 100 * MeV);
26
27 // The final-state particles declared above are clustered using FastJet with
28 // the anti-kT algorithm and a jet-radius parameter 0.4
29 // muons and neutrinos are excluded from the clustering
30 FastJets jetfs(fs, JetAlg::ANTIKT, 0.4);
31 declare(jetfs, "jets");
32
33 // FinalState of prompt photons and bare muons and electrons in the event
34 PromptFinalState photons(Cuts::abspid == PID::PHOTON);
35 PromptFinalState bare_leps(Cuts::abspid == PID::MUON || Cuts::abspid == PID::ELECTRON);
36
37 // Dress the prompt bare leptons with prompt photons within dR < 0.1,
38 // and apply some fiducial cuts on the dressed leptons
39 Cut lepton_cuts = Cuts::abseta < 2.5 && Cuts::pT > 10 * GeV;
40 LeptonFinder dressed_leps(bare_leps, photons, 0.1, lepton_cuts);
41 declare(dressed_leps, "leptons");
42
43 // Missing momentum
44 declare(MissingMomentum(fs), "MET");
45
46 // Book histograms
47 book(_h_xsec_exp, 1, 1, 1);
48 }
49
50 /// Perform the per-event analysis
51 void analyze(const Event& event) {
52 // Retrieve dressed leptons, sorted by pT
53 const Particles& leptons = apply<LeptonFinder>(event, "leptons").particlesByPt();
54
55 // Retrieve clustered jets, sorted by pT, with a minimum pT cut
56 Jets jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::pT > 30 * GeV && Cuts::abseta < 4.7);
57
58 // Retrieve Etmiss variable
59 const double Etmiss = apply<MissingMomentum>(event, "MET").missingEt();
60
61 // Remove jets with dressed leptons within dR < 0.4
62 idiscardIfAnyDeltaRLess(jets, leptons, 0.4);
63 if (jets.size() < 2) vetoEvent;
64
65 // Veto event if there is at least one b-jet
66 for (const Jet& j : jets) {
67 if (j.bTagged()) vetoEvent;
68 }
69
70 // Apply lepton cuts
71 if (leptons.size() < 2) vetoEvent;
72 if (leptons[0].pT() < 25 * GeV || leptons[1].pT() < 13 * GeV) vetoEvent;
73 if (leptons.size() > 2 && leptons[2].pT() > 10 * GeV) vetoEvent;
74 if (leptons[0].pid() * leptons[1].pid() != -11 * 13 &&
75 leptons[0].pid() * leptons[1].pid() != -11 * 11 &&
76 leptons[0].pid() * leptons[1].pid() != -13 * 13) vetoEvent;
77
78 FourMomentum l1 = leptons[0].momentum();
79 FourMomentum l2 = leptons[1].momentum();
80 FourMomentum ll = (l1 + l2);
81 const double mll = ll.mass();
82
83 if (ll.pT() < 30*GeV) vetoEvent;
84 if (mll < 50*GeV) vetoEvent;
85
86 // Apply a missing-momentum cut
87 if (Etmiss < 20*GeV) vetoEvent;
88
89 FourMomentum j1 = jets[0].momentum();
90 FourMomentum j2 = jets[1].momentum();
91 const double mjj = (j1 + j2).mass();
92 const double detajj = abs(j1.eta() - j2.eta());
93
94 // Apply VBS cuts
95 if (mjj < 300 * GeV || detajj < 2.5) vetoEvent;
96
97 _h_xsec_exp->fill("Exclusive"s);
98 }
99
100 /// Normalise histograms etc., after the run
101 void finalize() {
102 const double norm = crossSection() / femtobarn / sumOfWeights();
103
104 scale(_h_xsec_exp, norm);
105 }
106
107 //@}
108
109 /// @name Histograms
110 //@{
111 BinnedHistoPtr<string> _h_xsec_exp;
112 //@}
113 };
114
115 RIVET_DECLARE_PLUGIN(CMS_2022_I2080534);
116
117} // namespace Rivet
|