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