Rivet analyses referenceCMS_2020_I1814328W$^+$W$^-$ boson pair production in proton-proton collisions at $\sqrt{s} =$ 13 TeVExperiment: CMS (LHC) Inspire ID: 1814328 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
A measurement of the W W boson pair production cross section in proton-proton collisions at sqrt(s) = 13 TeV is presented. The data used in this study were collected with the CMS detector at the LHC and correspond to an integrated luminosity of 35.9 fb-1. The W+W- candidate events are selected by first requiring two oppositely-charged leptons (electrons or muons). Two methods for reducing background contributions are employed. In the first one, a sequence of requirements on kinematic quantities is applied allowing a measurement of the total cross section: 117.6 +/- 6.8 pb which agrees well with the theoretical cross section. Fiducial cross sections are also reported for events with zero jets or one jet, and the change in the zero-jet fiducial cross section with the jet pT threshold is measured. Normalized differential cross sections are reported within the fiducial region; these are compared to theoretical predictions based on quantum chromodynamics at next-to-next-to-leading-order accuracy. A second method for suppressing background contributions employs two random forest classifiers. The analysis based on this method includes a measurement of the total cross section and also a measurement of the normalized jet multiplicity distribution in W+W- events. Finally, a dilepton invariant mass distribution is used to probe for physics beyond the standard model in the context of an effective field theory and limits on the presence of dimension-6 operators are derived. Source code: CMS_2020_I1814328.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Projections/DressedLeptons.hh"
6#include "Rivet/Projections/MissingMomentum.hh"
7#include "Rivet/Projections/PromptFinalState.hh"
8
9namespace Rivet {
10
11
12 /// @brief Measurements of W+W- boson pair production in pp collisions at 13 TeV
13 class CMS_2020_I1814328 : public Analysis {
14 public:
15
16 /// Constructor
17 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2020_I1814328);
18
19
20 /// @name Analysis methods
21 /// @{
22
23 /// Book histograms and initialise projections before the run
24 void init() {
25
26 // Initialise and register projections
27
28 // The basic final-state projection:
29 // all final-state particles within
30 // the given eta acceptance
31 const FinalState fs(Cuts::abseta < 4.9);
32 const FinalState fsjet4p5(Cuts::abseta < 4.5);
33 const FinalState fsjet2p5(Cuts::abseta < 2.5);
34
35 // The final-state particles declared above are clustered using FastJet with
36 // the anti-kT algorithm and a jet-radius parameter 0.4
37 // muons and neutrinos are excluded from the clustering
38 FastJets jet4p5fs(fsjet4p5, FastJets::ANTIKT, 0.4);
39 declare(jet4p5fs, "jets4p5");
40 FastJets jet2p5fs(fsjet2p5, FastJets::ANTIKT, 0.4);
41 declare(jet2p5fs, "jets2p5");
42
43 // FinalState of prompt photons and bare muons and electrons in the event
44 PromptFinalState photons(Cuts::abspid == PID::PHOTON);
45 PromptFinalState bare_leps(Cuts::abspid == PID::MUON || Cuts::abspid == PID::ELECTRON);
46 bare_leps.acceptTauDecays(false);
47
48 // Dress the prompt bare leptons with prompt photons within dR < 0.1,
49 // and apply some fiducial cuts on the dressed leptons
50 Cut lepton_cuts = Cuts::abseta < 2.5 && Cuts::pT > 25*GeV;
51 DressedLeptons dressed_leps(photons, bare_leps, 0.1, lepton_cuts);
52 declare(dressed_leps, "leptons");
53
54 // Missing momentum
55 declare(MissingMomentum(fs), "MET");
56
57 // Book histograms
58 book(_h_WW_njets_norm , 2, 1, 1);
59 book(_h_WW_mll_norm , 4, 1, 1);
60 book(_h_WW_ptlmax_norm, 5, 1, 1);
61 book(_h_WW_ptlmin_norm, 6, 1, 1);
62 book(_h_WW_dphill_norm, 7, 1, 1);
63 book(_h_WW_njet0 , 8, 1, 1);
64
65 }
66
67 /// Perform the per-event analysis
68 void analyze(const Event& event) {
69
70 // Apply a missing-momentum cut
71 if (apply<MissingMomentum>(event, "MET").missingPt() < 20*GeV) return;
72
73 // Retrieve dressed leptons, sorted by pT
74 vector<DressedLepton> leptons = apply<DressedLeptons>(event, "leptons").dressedLeptons();
75
76 // Retrieve clustered jets, sorted by pT, with a minimum pT cut
77 Jets jets25 = apply<FastJets>(event, "jets4p5").jetsByPt(Cuts::pT > 25*GeV);
78 Jets jets30 = apply<FastJets>(event, "jets4p5").jetsByPt(Cuts::pT > 30*GeV);
79 Jets jets35 = apply<FastJets>(event, "jets4p5").jetsByPt(Cuts::pT > 35*GeV);
80 Jets jets45 = apply<FastJets>(event, "jets4p5").jetsByPt(Cuts::pT > 45*GeV);
81 Jets jets60 = apply<FastJets>(event, "jets4p5").jetsByPt(Cuts::pT > 60*GeV);
82 Jets jetsNj = apply<FastJets>(event, "jets2p5").jetsByPt(Cuts::pT > 30*GeV);
83
84 // Remove all jets within dR < 0.4 of a dressed lepton
85 idiscardIfAnyDeltaRLess(jets25, leptons, 0.4);
86 idiscardIfAnyDeltaRLess(jets30, leptons, 0.4);
87 idiscardIfAnyDeltaRLess(jets35, leptons, 0.4);
88 idiscardIfAnyDeltaRLess(jets45, leptons, 0.4);
89 idiscardIfAnyDeltaRLess(jets60, leptons, 0.4);
90 idiscardIfAnyDeltaRLess(jetsNj, leptons, 0.4);
91
92 if (leptons.size() == 2 && leptons[0].pid() * leptons[1].pid() < 0) {
93 FourMomentum dilCand = leptons[0].momentum() + leptons[1].momentum();
94 if (dilCand.mass() > 20*GeV && dilCand.pT() > 30*GeV){
95 double ptlmax = leptons[0].pT(); double ptlmin = leptons[1].pT();
96 if (ptlmax < ptlmin) {
97 ptlmax = leptons[1].pT(); ptlmin = leptons[0].pT();
98 }
99
100 if (leptons[0].abspid() != leptons[1].abspid()) {
101 _h_WW_njets_norm ->fill(min((double)jetsNj.size()+1, 2.999));
102 _h_WW_mll_norm ->fill(min(dilCand.mass()/GeV, 1499.999));
103 _h_WW_ptlmax_norm->fill(min(ptlmax/GeV, 399.999));
104 _h_WW_ptlmin_norm->fill(min(ptlmin/GeV, 149.999));
105 _h_WW_dphill_norm->fill(deltaPhi(leptons[0], leptons[1]));
106 }
107
108 if (jets25.size() == 0) _h_WW_njet0->fill(1.0);
109 if (jets30.size() == 0) _h_WW_njet0->fill(2.0);
110 if (jets35.size() == 0) _h_WW_njet0->fill(3.0);
111 if (jets45.size() == 0) _h_WW_njet0->fill(4.0);
112 if (jets60.size() == 0) _h_WW_njet0->fill(5.0);
113
114 }
115 }
116
117 }
118
119
120 /// @todo Replace with barchart()
121 void normalizeToSum(Histo1DPtr hist) {
122 double sum = 0.;
123 for (size_t i = 0; i < hist->numBins(); ++i) {
124 sum += hist->bin(i).height();
125 }
126 scale(hist, 1./sum);
127 }
128
129
130 /// Normalise histograms etc., after the run
131 void finalize() {
132
133 double norm = (sumOfWeights() != 0) ? crossSection()/picobarn/sumOfWeights() : 1.0;
134
135 normalizeToSum(_h_WW_njets_norm );
136 normalizeToSum(_h_WW_mll_norm );
137 normalizeToSum(_h_WW_ptlmax_norm);
138 normalizeToSum(_h_WW_ptlmin_norm);
139 normalizeToSum(_h_WW_dphill_norm);
140
141 scale(_h_WW_njet0, norm);
142
143 }
144
145 /// @}
146
147
148 /// @name Histograms
149 /// @{
150 Histo1DPtr _h_WW_njets_norm;
151 Histo1DPtr _h_WW_mll_norm, _h_WW_ptlmax_norm, _h_WW_ptlmin_norm, _h_WW_dphill_norm;
152 Histo1DPtr _h_WW_njet0;
153 /// @}
154
155 };
156
157
158
159 RIVET_DECLARE_PLUGIN(CMS_2020_I1814328);
160
161}
|