Rivet analyses referenceCMS_2022_I2129461tW single top-quark differential cross sections at 13 TeVExperiment: CMS (LHC) Inspire ID: 2129461 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
Measurements of the inclusive and normalised differential cross sections are presented for the production of single top quarks in association with a W boson in proton-proton collisions at a centre-of-mass energy of 13 TeV. The data used were recorded with the CMS detector at the LHC during 2016-2018, and correspond to an integrated luminosity of $138 fb^{-1}$. Events containing one electron and one muon in the final state are analysed. For the inclusive measurement, a multivariate discriminant, exploiting the kinematic properties of the events is used to separate the signal from the dominant $t\bar{t}$ background. A cross section of $79.2 \pm 0.9 (stat) ^{7.7}_{-8.0} (syst) \pm 1.2 (lumi) pb$ is obtained, consistent with the predictions of the standard model. For the differential measurements, a fiducial region is defined according to the detector acceptance, and the requirement of exactly one jet coming from the fragmentation of a bottom quark. The resulting distributions are unfolded to particle level and agree with the predictions at next-to-leading order in perturbative quantum chromodynamics. Source code: CMS_2022_I2129461.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Projections/ChargedLeptons.hh"
6#include "Rivet/Projections/LeptonFinder.hh"
7#include "Rivet/Projections/MissingMomentum.hh"
8#include "Rivet/Projections/PromptFinalState.hh"
9#include "Rivet/Projections/VetoedFinalState.hh"
10
11namespace Rivet {
12
13
14 /// @brief Normalised differential cross sections for 13 TeV tW-channel single top-quark production
15 class CMS_2022_I2129461 : public Analysis {
16 public:
17
18 /// Constructor
19 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2022_I2129461);
20
21
22 /// Book histograms and initialise projections before the run
23 void init() {
24
25 // final state of all stable particles
26 Cut particle_cut = (Cuts::abseta < 5.0) and (Cuts::pT > 0.*MeV);
27 FinalState fs(particle_cut);
28
29 // select charged leptons
30 ChargedLeptons charged_leptons(fs);
31
32 // select final state photons for dressed lepton clustering
33 IdentifiedFinalState photons(fs);
34 photons.acceptIdPair(PID::PHOTON);
35
36 // select final state prompt charged leptons
37 PromptFinalState prompt_leptons(charged_leptons);
38 prompt_leptons.acceptMuonDecays(true);
39 prompt_leptons.acceptTauDecays(true);
40
41 // select final state prompt photons
42 PromptFinalState prompt_photons(photons);
43 prompt_photons.acceptMuonDecays(true);
44 prompt_photons.acceptTauDecays(true);
45
46 // Dressed leptons from selected prompt charged leptons and photons
47 Cut lepton_cut = ((Cuts::abseta < 2.4) &&
48 (Cuts::pT > 20.*GeV) &&
49 (((Cuts::abspid == PID::ELECTRON) &&
50 ((Cuts::abseta < 1.4442) or (Cuts::abseta > 1.566))) ||
51 (Cuts::abspid == PID::MUON)));
52
53 LeptonFinder dressed_leptons(prompt_photons, prompt_leptons, 0.1, lepton_cut);
54 declare(dressed_leptons, "LeptonFinder");
55
56 // Jets
57 VetoedFinalState fsForJets(fs);
58 fsForJets.addVetoOnThisFinalState(dressed_leptons);
59 declare(
60 // excludes all neutrinos by default
61 FastJets(fsForJets, JetAlg::ANTIKT, 0.4),
62 "Jets");
63
64 // pTmiss
65 declare(MissingMomentum(fs), "MET");
66
67 // book histograms
68 book(_hist_norm_lep1_pt, "d06-x01-y08");
69 book(_hist_norm_lep1lep2jet1_pz, "d07-x01-y08");
70 book(_hist_norm_jet1_pt, "d08-x01-y08");
71 book(_hist_norm_lep1lep2jet1_m, "d09-x01-y08");
72 book(_hist_norm_lep1lep2_dphi, "d10-x01-y08");
73 book(_hist_norm_lep1lep2jet1met_mt, "d11-x01-y08");
74 }
75
76
77 /// @brief Perform the per-event analysis
78 void analyze(const Event& event) {
79 DressedLeptons dressedLeptons = apply<LeptonFinder>(event, "LeptonFinder").dressedLeptons();
80
81 // Require at least two dressed leptons
82 if (dressedLeptons.size() < 2) vetoEvent;
83
84 // Require that the leading ones are an electron and a muon of opposite charge
85 if (abs(dressedLeptons.at(0).pid() + dressedLeptons.at(1).pid()) != 2) vetoEvent;
86
87 // Require that the leading lepton has at least 25 GeV of pT
88 if (dressedLeptons.at(0).pt() <= 25.*GeV) vetoEvent;
89
90 // Require that all identified lepton pairs have at least 20 GeV of invariant mass
91 for (int iL = 0; iL < int(dressedLeptons.size()); iL++) {
92 for (int jL = iL + 1; jL < int(dressedLeptons.size()); jL++) {
93 if ((dressedLeptons.at(iL).momentum() + dressedLeptons.at(jL).momentum()).mass() <= 20.*GeV) vetoEvent;
94 }
95 }
96
97 // Jet object ID
98 Cut jet_cut((Cuts::abseta < 2.4) and (Cuts::pT > 20.*GeV));
99 vector<Jet> jets = apply<FastJets>(
100 event,
101 "Jets"
102 ).jets(jet_cut);
103
104 // ignore jets that overlap with dressed leptons within dR < 0.4 and collect them in loose jets and jets categories
105 Jets cleanedLooseJets;
106 Jets cleanedJets;
107 idiscardIfAnyDeltaRLess(jets, dressedLeptons, 0.4);
108 for (const Jet& jet: jets) {
109 if (jet.pt() > 30.*GeV) cleanedJets.push_back(jet);
110 else cleanedLooseJets.push_back(jet);
111 }
112
113 if (cleanedJets.size() != 1) vetoEvent; // select events with exactly one jet...
114 if (not cleanedJets.at(0).bTagged()) vetoEvent; // ...that must be b-tagged...
115 if (cleanedLooseJets.size() > 0) vetoEvent; // ...and no loose jets
116
117 // fill the histograms
118 _hist_norm_lep1_pt->fill(min(max(dressedLeptons.at(0).pt(), 26.*GeV), 149.*GeV) / GeV);
119 FourMomentum llj = dressedLeptons.at(0).momentum() + dressedLeptons.at(1).momentum() + cleanedJets.at(0).momentum();
120 _hist_norm_lep1lep2jet1_pz->fill(min(max(abs((llj).pz()), 1.*GeV), 449.*GeV) / GeV);
121 _hist_norm_jet1_pt->fill(min(max(cleanedJets.at(0).pt(), 31.*GeV), 149.*GeV) / GeV);
122 _hist_norm_lep1lep2jet1_m->fill(min(max((llj).mass(), 51.*GeV), 399.*GeV) / GeV);
123 _hist_norm_lep1lep2_dphi->fill( abs(deltaPhi(dressedLeptons.at(0).phi(), dressedLeptons.at(1).phi())) / PI );
124
125 Vector3 tmpPtmiss3 = apply<MissingMomentum>(event, "MET").vectorPt();
126 tmpPtmiss3.setZ(0.);
127 FourMomentum ptmiss = FourMomentum(tmpPtmiss3.mod(), -tmpPtmiss3.x(), -tmpPtmiss3.y(), 0.);
128 FourMomentum tmpM = llj + ptmiss;
129 _hist_norm_lep1lep2jet1met_mt->fill(min(max(sqrt(tmpM.E2() - tmpM.pz2()), 101.*GeV), 499.*GeV) / GeV);
130 }
131
132 /// @brief Normalise histograms after the run
133 void finalize() {
134 normalize(_hist_norm_lep1_pt);
135 normalize(_hist_norm_lep1lep2jet1_pz);
136 normalize(_hist_norm_jet1_pt);
137 normalize(_hist_norm_lep1lep2jet1_m);
138 normalize(_hist_norm_lep1lep2_dphi);
139 normalize(_hist_norm_lep1lep2jet1met_mt);
140 }
141
142 private:
143 // Declaration of histograms
144 Histo1DPtr _hist_norm_lep1_pt;
145 Histo1DPtr _hist_norm_lep1lep2jet1_pz;
146 Histo1DPtr _hist_norm_jet1_pt;
147 Histo1DPtr _hist_norm_lep1lep2jet1_m;
148 Histo1DPtr _hist_norm_lep1lep2_dphi;
149 Histo1DPtr _hist_norm_lep1lep2jet1met_mt;
150 };
151
152 RIVET_DECLARE_PLUGIN(CMS_2022_I2129461);
153}
|