Rivet analyses referenceATLAS_2019_I1738841Same-sign $WW$ Observation at 13 TeVExperiment: ATLAS (LHC) Inspire ID: 1738841 Status: VALIDATED Authors:
Beam energies: (6500.0, 6500.0) GeV Run details:
This Letter presents the observation and measurement of electroweak production of a same-sign $W$ boson pair in association with two jets using 36.1 fb$^{-1}$ of proton-proton collision data recorded at a center-of-mass energy of $\sqrt{s}=13$ TeV by the ATLAS detector at the Large Hadron Collider. The analysis is performed in the detector fiducial phase-space region, defined by the presence of two same-sign leptons, electron or muon, and at least two jets with a large invariant mass and rapidity difference. A total of 122 candidate events are observed for a background expectation of $69 \pm 7$ events, corresponding to an observed signal significance of 6.5 standard deviations. The measured fiducial signal cross section is $\sigma^{\textrm {fid.}}=2.89^{+0.51}_{-0.48} \textrm{(stat.)}\ ^{+0.29}_{-0.28} \textrm{(syst.)}\textrm{ fb}$. Source code: ATLAS_2019_I1738841.cc 1#include "Rivet/Analysis.hh"
2#include "Rivet/Projections/FinalState.hh"
3#include "Rivet/Projections/FastJets.hh"
4#include "Rivet/Projections/PromptFinalState.hh"
5#include "Rivet/Projections/LeptonFinder.hh"
6#include "Rivet/Projections/MissingMomentum.hh"
7#include "Rivet/Projections/VetoedFinalState.hh"
8
9namespace Rivet {
10
11
12 /// @brief same-sign WW at 13 TeV
13 class ATLAS_2019_I1738841 : public Analysis {
14 public:
15
16 /// Constructor
17 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2019_I1738841);
18
19
20 /// @name Analysis methods
21 /// @{
22
23 /// Book histograms and initialise projections before the run
24 void init() {
25
26 FinalState allLeps(Cuts::abspid == PID::ELECTRON || Cuts::abspid == PID::MUON);
27 PromptFinalState promptLeps(allLeps);
28 FinalState photons(Cuts::abspid == PID::PHOTON);
29
30 Cut dressedLep_cuts = (Cuts::abseta < 2.5) && (Cuts::pT > 27*GeV);
31 LeptonFinder dressedLeps(promptLeps, photons, 0.1, dressedLep_cuts);
32 declare(dressedLeps, "dressedLeptons");
33
34 // veto on leptons from prompt tau decays
35 VetoedFinalState lepsFromTaus(PromptFinalState(allLeps, TauDecaysAs::PROMPT));
36 lepsFromTaus.addVetoOnThisFinalState(promptLeps);
37 LeptonFinder vetoLeps(lepsFromTaus, photons, 0.1, dressedLep_cuts);
38 declare(vetoLeps, "vetoLeptons");
39
40 declare(MissingMomentum(), "eTmiss");
41
42 VetoedFinalState vfs(FinalState(Cuts::abseta < 4.5));
43 vfs.addVetoOnThisFinalState(dressedLeps);
44 FastJets jets(vfs, JetAlg::ANTIKT, 0.4, JetMuons::ALL, JetInvisibles::NONE);
45 declare(jets, "jets");
46
47 // Histograms
48 book(_c, 1, 1, 1);
49
50 }
51
52
53 /// Perform the per-event analysis
54 void analyze(const Event& event) {
55
56 if (apply<LeptonFinder>(event, "vetoLeptons").dressedLeptons().size()) vetoEvent;
57
58 const DressedLeptons &dressedLeptons = apply<LeptonFinder>(event, "dressedLeptons").dressedLeptons();
59 if (dressedLeptons.size() != 2) vetoEvent;
60 const Particle& lep1 = dressedLeptons[0];
61 const Particle& lep2 = dressedLeptons[1];
62
63 const Jets& jets = apply<FastJets>(event,"jets").jetsByPt(Cuts::pT > 35*GeV);
64 if (jets.size() < 2) vetoEvent;
65
66 // pT of the two leading jets
67 const Jet& tagJet1 = jets[0];
68 const Jet& tagJet2 = jets[1];
69 if (tagJet1.pT() < 65.) vetoEvent;
70 if (tagJet2.pT() < 35.) vetoEvent;
71
72 // leptons isolation from jets (dR >= 0.3)
73 float f_lj_dRmin = 999;
74 for (const Jet& j : jets) {
75 float dR = deltaR(j, lep1);
76 if(dR < f_lj_dRmin) f_lj_dRmin = dR;
77 dR = deltaR(j, lep2);
78 if(dR < f_lj_dRmin) f_lj_dRmin = dR;
79 }
80 if (f_lj_dRmin < 0.3) vetoEvent;
81
82 // same sign leptons
83 if (lep1.pid() * lep2.pid() < 0) vetoEvent;
84
85 // lepton isolation
86 const float f_ll_dR = deltaR(lep1.momentum(), lep2.momentum());
87 if (f_ll_dR < 0.3) vetoEvent;
88
89 // m_ll > 20 GeV
90 const float f_ll_m = (lep1.momentum() + lep2.momentum()).mass()/GeV;
91 if (f_ll_m < 20.) vetoEvent;
92
93 // MET > 30 GeV
94 const double ETmiss = apply<MissingMomentum>(event, "eTmiss").missingPt()/GeV;
95 if (ETmiss < 30.) vetoEvent;
96
97 // m_jj >= 500 GeV
98 const float f_tjets_m = (tagJet1.momentum() + tagJet2.momentum()).mass()/GeV;
99 if (f_tjets_m < 500.) vetoEvent;
100
101 // dY_jj >= 2
102 const float f_tjets_dY = deltaRap(tagJet1, tagJet2);
103 if (f_tjets_dY < 2) vetoEvent;
104
105 // --- fill histograms ---
106 _c->fill();
107 }
108
109
110 /// Normalise histograms etc., after the run
111 void finalize() {
112 scale(_c, crossSection() / femtobarn / sumOfWeights());
113 }
114
115 /// @}
116
117
118 /// @name Histograms
119 /// @{
120 CounterPtr _c;
121 /// @}
122
123 };
124
125 RIVET_DECLARE_PLUGIN(ATLAS_2019_I1738841);
126}
|