rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2021_I1852328

WW + $\geq$1 jet production at 13 TeV
Experiment: ATLAS (LHC)
Inspire ID: 1852328
Status: VALIDATED
Authors:
  • Hannes Mildner
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • pp -> WW + jet production at 13 TeV

Fiducial and differential measurements of $W^+W^-$ production in events with at least one hadronic jet are presented. These measurements are sensitive to the properties of gauge boson self interactions and provide a test of perturbative quantum chromodynamics and the electroweak theory. The analysis is performed using proton--proton collision data collected at $\sqrt{s}=13$ TeV with the ATLAS experiment, corresponding to an integrated luminosity of 139 fb$^{-1}$. Events with exactly one oppositely charged electron-muon pair and at least one hadronic jet with a transverse momentum of $p_{\mathrm{T}}>30$ GeV and a pseudorapidity of $|\eta|<4.5$ are selected. After subtracting the background contributions and correcting for detector effects, the jet-inclusive $W^+W^-+\ge 1$ fiducial cross-section and $W^+W^-$+jets differential cross-sections with respect to several kinematic variables are measured, for the first time at the LHC. These measurements include leptonic quantities, like the lepton transverse momenta as well as the transverse mass of the $W^+W^-$ system, and also jet related observables like the leading jet transverse momentum and the jet multiplicity.

Source code: ATLAS_2021_I1852328.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FastJets.hh"
  4#include "Rivet/Projections/FinalState.hh"
  5#include "Rivet/Projections/MissingMomentum.hh"
  6#include "Rivet/Projections/PromptFinalState.hh"
  7#include "Rivet/Projections/LeptonFinder.hh"
  8#include "Rivet/Projections/VetoedFinalState.hh"
  9
 10namespace Rivet {
 11
 12
 13  /// @brief WW production at 13 TeV
 14  class ATLAS_2021_I1852328 : public Analysis {
 15  public:
 16
 17    /// Constructor
 18    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2021_I1852328);
 19
 20
 21    /// Book histograms and initialise projections before the run
 22    void init() {
 23
 24      const FinalState fs(Cuts::abseta < 5.);
 25
 26      // Project photons for dressing
 27      FinalState photons(Cuts::abspid == PID::PHOTON);
 28
 29      // Cuts for leptons
 30      Cut lepton_cuts = (Cuts::abseta < 2.5) && (Cuts::pT > 27*GeV);
 31
 32      // Project dressed leptons (e/mu not from tau) with pT > 27 GeV and |eta| < 2.5
 33      PromptFinalState lep_bare(Cuts::abspid == PID::MUON || Cuts::abspid == PID::ELECTRON);
 34      lep_bare.acceptTauDecays(false);
 35      LeptonFinder lep_dressed(lep_bare, photons, 0.1, lepton_cuts);
 36      declare(lep_dressed,"lep_dressed");
 37
 38      // Leptons (+ dressed photons) to be removed from jets
 39      PromptFinalState bare_mu(Cuts::abspid == PID::MUON, TauDecaysAs::PROMPT);
 40      LeptonFinder all_dressed_mu(bare_mu, photons, 0.1, Cuts::abseta < 2.5);
 41      PromptFinalState bare_el(Cuts::abspid == PID::ELECTRON, TauDecaysAs::PROMPT);
 42      LeptonFinder all_dressed_el(bare_el, photons, 0.1, Cuts::abseta < 2.5);
 43
 44      // Define hadrons as everything but dressed leptons (for jet clustering)
 45      VetoedFinalState hadrons(fs);
 46      hadrons.addVetoOnThisFinalState(all_dressed_el);
 47      hadrons.addVetoOnThisFinalState(all_dressed_mu);
 48      declare(hadrons, "hadrons");
 49
 50      // Get MET from generic invisibles
 51      MissingMomentum mm(fs);
 52      declare(mm, "met");
 53
 54      // Project jets
 55      FastJets jets(hadrons, JetAlg::ANTIKT, 0.4, JetMuons::ALL, JetInvisibles::DECAY);
 56      declare(jets, "jets");
 57
 58      book(_h["xs_inf"], 1, 1, 1);
 59      book(_h["xs_bveto_inf"], 1, 1, 2);
 60      book(_h["lep0pt"], 2, 1, 1);
 61      book(_h["lep0pt_bveto"], 2, 1, 2);
 62      book(_h["lep1pt"], 5, 1, 1);
 63      book(_h["lep1pt_bveto"], 5, 1, 2);
 64      book(_h["jet0pt"], 8, 1, 1);
 65      book(_h["jet0pt_bveto"], 8, 1, 2);
 66      book(_h["htjet"], 11, 1, 1);
 67      book(_h["htjet_bveto"], 11, 1, 2);
 68      book(_h["st"], 14, 1, 1);
 69      book(_h["st_bveto"], 14, 1, 2);
 70      book(_h["mt"], 17, 1, 1);
 71      book(_h["mt_bveto"], 17, 1, 2);
 72      book(_h["mll_inf"], 20, 1, 1);
 73      book(_h["mll_bveto_inf"], 20, 1, 2);
 74      book(_h["ptll"], 23, 1, 1);
 75      book(_h["ptll_bveto"], 23, 1, 2);
 76      book(_h["dphill_inf"], 26, 1, 1);
 77      book(_h["dphill_bveto_inf"], 26, 1, 2);
 78      book(_h["yll_inf"], 29, 1, 1);
 79      book(_h["yll_bveto_inf"], 29, 1, 2);
 80      book(_h["costhetastar_inf"], 32, 1, 1);
 81      book(_h["costhetastar_bveto_inf"], 32, 1, 2);
 82      book(_h["njet"], 35, 1, 1);
 83      book(_h["njet_bveto"], 35, 1, 2);
 84
 85      book(_h["mll_jet0pt200_inf"], 38, 1, 1);
 86      book(_h["mll_jet0pt200_bveto_inf"], 38, 1, 2);
 87      book(_h["dphill_jet0pt200_inf"], 41, 1, 1);
 88      book(_h["dphill_jet0pt200_bveto_inf"], 41, 1, 2);
 89
 90      book(_h["dphil1j0_lep0pt200_inf"], 44, 1, 1);
 91      book(_h["dphil1j0_lep0pt200_bveto_inf"], 44, 1, 2);
 92      book(_h["drl1j0_lep0pt200_inf"], 47, 1, 1);
 93      book(_h["drl1j0_lep0pt200_bveto_inf"], 47, 1, 2);
 94      book(_h["rl1l0_lep0pt200_inf"], 50, 1, 1);
 95      book(_h["rl1l0_lep0pt200_bveto_inf"], 50, 1, 2);
 96      book(_h["rl1j0_lep0pt200"], 53, 1, 1);
 97      book(_h["rl1j0_lep0pt200_bveto"], 53, 1, 2);
 98
 99    }
100
101
102    /// Perform the per-event analysis
103    void analyze(const Event& event) {
104
105      // Get invisible
106      const MissingMomentum& met = apply<MissingMomentum>(event, "met");
107
108      // Find leptons and sort by pT
109      const Particles& leptons = apply<ParticleFinder>(event, "lep_dressed").particlesByPt();
110
111      const Jets& jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::absrap < 4.5 && Cuts::pT > 30*GeV);
112      const Jets& btagging_jets = apply<FastJets>(event, "jets").jetsByPt(Cuts::absrap < 2.5 && Cuts::pT > 20*GeV);
113
114      // Define observables
115
116      const FourMomentum dilep = leptons.size()>1 ? leptons[0].momentum() + leptons[1].momentum() : FourMomentum(0,0,0,0);
117      const double lep0pt = leptons.size()>0 ? leptons[0].pT()/GeV : -1;
118      const double lep1pt = leptons.size()>1 ? leptons[1].pT()/GeV : -1;
119      const double mll = leptons.size()>1 ? dilep.mass()/GeV : -1;
120      const double ptll = leptons.size()>1 ? dilep.pT()/GeV : -1;
121      const double yll = leptons.size()>1 ? dilep.absrap(): -1;
122      const double dphill = leptons.size()>1 ? fabs(deltaPhi(leptons[0], leptons[1])) : -1.;
123      const double costhetastar = leptons.size()>1 ? fabs(tanh((leptons[0].eta() - leptons[1].eta()) / 2)) : -1;
124      const double jet0pt = jets.size()>0 ? jets[0].pT()/GeV : -1;
125      const double dphil1j0 = leptons.size()>1 and jets.size()>0 ? fabs(deltaPhi(leptons[1], jets[0])) : -1.;
126      const double drl1j0 = leptons.size()>1 and jets.size()>0 ? fabs(deltaR(leptons[1], jets[0])) : -1.;
127      const double rl1l0 = leptons.size()>1 ? leptons[1].pT()/leptons[0].pT() : -1;
128      const double rl1j0 = leptons.size()>1 and jets.size()>0 ? leptons[1].pT()/jets[0].pT() : -1;
129      Vector3 dilep_tvec(dilep.px(), dilep.py(),0);
130      Vector3 met_tvec = met.vectorMissingPt();
131      double et_ll = std::sqrt(dilep_tvec.mod2() + dilep.mass()*dilep.mass());
132      const double mt =sqrt(std::pow(et_ll + met.met() ,2) - (dilep_tvec + met_tvec).mod2());
133      const int njet = jets.size();
134
135      size_t nbtags = count(btagging_jets, hasBTag());
136      const double htjet = sum(jets, pT, 0.0);
137      const double st = sum(leptons, pT, htjet);
138
139      if (leptons.size()!=2 ) vetoEvent;
140      if (leptons[0].abspid() == leptons[1].abspid()) vetoEvent;
141      if (leptons[0].pid()*leptons[1].pid()>0) vetoEvent;
142      if (dilep.mass() <= 85*GeV)  vetoEvent;
143      if (jets.empty()) vetoEvent;
144
145      _h["xs_inf"]->fill(0.5);
146      _h["lep0pt"]->fill(lep0pt);
147      _h["lep1pt"]->fill(lep1pt);
148      _h["mll_inf"]->fill(mll);
149      _h["ptll"]->fill(ptll);
150      _h["yll_inf"]->fill(yll);
151      _h["dphill_inf"]->fill(dphill);
152      _h["costhetastar_inf"]->fill(costhetastar);
153      _h["jet0pt"]->fill(jet0pt);
154      _h["mt"]->fill(mt);
155      _h["njet"]->fill(njet);
156      _h["htjet"]->fill(htjet/GeV);
157      _h["st"]->fill(st/GeV);
158
159      if (lep0pt>200){
160	      _h["dphil1j0_lep0pt200_inf"]->fill(dphil1j0);
161	      _h["drl1j0_lep0pt200_inf"]->fill(drl1j0);
162	      _h["rl1l0_lep0pt200_inf"]->fill(rl1l0);
163	      _h["rl1j0_lep0pt200"]->fill(rl1j0);
164      }
165      if(jet0pt>200){
166	     _h["mll_jet0pt200_inf"]->fill(mll);
167	     _h["dphill_jet0pt200_inf"]->fill(dphill);
168      }
169
170
171      if(nbtags == 0) {
172        _h["xs_bveto_inf"]->fill(0.5);
173        _h["lep0pt_bveto"]->fill(lep0pt);
174        _h["lep1pt_bveto"]->fill(lep1pt);
175        _h["mll_bveto_inf"]->fill(mll);
176        _h["ptll_bveto"]->fill(ptll);
177        _h["yll_bveto_inf"]->fill(yll);
178        _h["dphill_bveto_inf"]->fill(dphill);
179        _h["costhetastar_bveto_inf"]->fill(costhetastar);
180        _h["jet0pt_bveto"]->fill(jet0pt);
181        _h["mt_bveto"]->fill(mt);
182        _h["njet_bveto"]->fill(njet);
183        _h["htjet_bveto"]->fill(htjet);
184        _h["st_bveto"]->fill(st);
185
186         if (lep0pt>200){
187        	  _h["dphil1j0_lep0pt200_bveto_inf"]->fill(dphil1j0);
188 	          _h["drl1j0_lep0pt200_bveto_inf"]->fill(drl1j0);
189	          _h["rl1l0_lep0pt200_bveto_inf"]->fill(rl1l0);
190    	      _h["rl1j0_lep0pt200_bveto"]->fill(rl1j0);
191           }
192
193         if(jet0pt>200){
194	          _h["mll_jet0pt200_bveto_inf"]->fill(mll);
195	          _h["dphill_jet0pt200_bveto_inf"]->fill(dphill);
196          }
197      }
198    }
199
200
201    /// Normalise histograms etc., after the run
202    void finalize() {
203
204       const double sf = crossSectionPerEvent()/femtobarn;
205
206       scale(_h, sf);
207       for (auto& hist : _h) {
208          if (hist.first.find("_inf") != string::npos) {
209            const size_t nBins = hist.second->numBins();
210            auto& overflow = hist.second->bin(nBins+1);
211            hist.second->fill(hist.second->bin(nBins).xMid(), overflow.sumW());
212            overflow.reset();
213          }
214       }
215    }
216
217
218  private:
219
220    /// @name Histograms
221    map<string, Histo1DPtr> _h;
222
223  };
224
225
226  RIVET_DECLARE_PLUGIN(ATLAS_2021_I1852328);
227
228}