rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ATLAS_2022_I2077570

Z + high transverse momentum jets at ATLAS
Experiment: ATLAS (LHC)
Inspire ID: 2077570
Status: VALIDATED
Authors:
  • Alexandre Laurier
References: Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • pp -> Z + jets at 13 TeV with pT(jet) > 100 GeV

Cross-section measurements for a $Z$ boson produced in association with high-transverse-momentum jets ($p_\text{T}\geq 100$ GeV) and decaying into a charged-lepton pair ($e^+e^-$, $\mu^+\mu^-$) are presented. The measurements are performed using proton-proton collisions at $\sqrt{s} = 13$ TeV corresponding to an integrated luminosity of 139 fb${}^{-1}$ collected by the ATLAS experiment at the LHC. Measurements of angular correlations between the $Z$ boson and the closest jet are performed in events with at least one jet with $p_\text{T}\geq 500$ GeV. Event topologies of particular interest are the collinear emission of a $Z$ boson in dijet events and a boosted $Z$ boson recoiling against a jet. Fiducial cross sections are compared with state-of-the-art theoretical predictions. The data are found to agree with next-to-next-to-leading-order predictions by NNLOjet and with the next-to-leading-order multi-leg generators MadGraph5_aMC@NLO and Sherpa.

Source code: ATLAS_2022_I2077570.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FastJets.hh"
  4#include "Rivet/Projections/FinalState.hh"
  5#include "Rivet/Projections/IdentifiedFinalState.hh"
  6#include "Rivet/Projections/VetoedFinalState.hh"
  7#include "Rivet/Projections/PromptFinalState.hh"
  8#include "Rivet/Projections/DressedLeptons.hh"
  9
 10namespace Rivet {
 11
 12  /// Colinear Z + Jets in pp at 13 TeV
 13  class ATLAS_2022_I2077570 : public Analysis {
 14  public:
 15
 16    /// Constructor
 17    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2022_I2077570);
 18
 19    /// Book histograms and initialise projections before the run
 20    void init() {
 21
 22      // Get options
 23      _mode = 0;
 24      if ( getOption("LMODE") == "ELEL") _mode = 1;
 25      if ( getOption("LMODE") == "MUMU") _mode = 2;
 26
 27      // AntiKt4TruthWZJets prescription
 28      // Photons
 29      FinalState all_photons(Cuts::abspid == PID::PHOTON);
 30      
 31      // Muons
 32      PromptFinalState bare_mu(Cuts::abspid == PID::MUON, true); // true = use muons from prompt tau decays
 33      DressedLeptons all_dressed_mu(all_photons, bare_mu, 0.1, Cuts::abseta < 2.5, true);
 34      
 35      // Electrons
 36      PromptFinalState bare_el(Cuts::abspid == PID::ELECTRON, true); // true = use electrons from prompt tau decays
 37      DressedLeptons all_dressed_el(all_photons, bare_el, 0.1, Cuts::abseta < 2.5, true);
 38      
 39      //Jet forming
 40      VetoedFinalState vfs(FinalState(Cuts::abseta < 4.5));
 41      vfs.addVetoOnThisFinalState(all_dressed_el);
 42      vfs.addVetoOnThisFinalState(all_dressed_mu);
 43            
 44      FastJets jet(vfs, FastJets::ANTIKT, 0.4, JetAlg::Muons::ALL, JetAlg::Invisibles::NONE);
 45      declare(jet, "Jets");
 46
 47
 48      // Current definition of leptons + jets
 49      PromptFinalState photons(Cuts::abspid == PID::PHOTON);
 50      PromptFinalState electrons(Cuts::abspid == PID::ELECTRON);
 51      PromptFinalState muons(Cuts::abspid == PID::MUON);
 52
 53      // Kinematic cuts for leptons
 54      const Cut cuts_lep = Cuts::pT > 25*GeV && Cuts::abseta < 2.5;
 55
 56      DressedLeptons dressed_electrons(photons, electrons, 0.1, cuts_lep);
 57      declare(dressed_electrons, "DressedElectrons");
 58
 59      DressedLeptons dressed_muons(photons, muons, 0.1, cuts_lep);
 60      declare(dressed_muons, "DressedMuons");
 61
 62
 63      book(_h["ZpT"],        1, 1, 1);
 64      book(_h["jetpT"],      2, 1, 1);
 65      book(_h["NJets"],      3, 1, 1);
 66      book(_h["NJets500"],   4, 1, 1);
 67      book(_h["minDR"],      5, 1, 1);
 68      book(_h["rZJ"],        6, 1, 1);
 69      book(_h["rZJ_coll"],   7, 1, 1);
 70      book(_h["rZJ_b2b"],    8, 1, 1);
 71      book(_h["NJets_coll"], 9, 1, 1);
 72      book(_h["NJets_b2b"], 10, 1, 1);
 73      book(_h["HT"],        11, 1, 1);
 74      book(_h["minDR600"],  12, 1, 1);
 75      book(_h["NJets600"],  13, 1, 1);
 76
 77    }
 78
 79    /// Perform the per-event analysis
 80    void analyze(const Event& event) {
 81
 82
 83      // access fiducial electrons and muons
 84      const Particle *l1 = nullptr, *l2 = nullptr;
 85      auto muons = apply<DressedLeptons>(event, "DressedMuons").dressedLeptons();
 86      auto elecs = apply<DressedLeptons>(event, "DressedElectrons").dressedLeptons();
 87
 88      // lepton-jet overlap removal (note: muons are not included in jet finding)
 89      // Jets eta < 2.5, pT > 30GeV for overlap removal
 90      Jets jets = apply<FastJets>(event, "Jets").jetsByPt(Cuts::pT > 30*GeV && Cuts::absrap < 2.5);
 91      // Remove all jets within dR < 0.2 of a dressed lepton
 92      idiscardIfAnyDeltaRLess(jets, muons, 0.2);
 93      idiscardIfAnyDeltaRLess(jets, elecs, 0.2);
 94      
 95      // Remove leptons within dR < 0.4 of a jet
 96      idiscardIfAnyDeltaRLess(muons, jets, 0.4);
 97      idiscardIfAnyDeltaRLess(elecs, jets, 0.4);
 98
 99      int nElecs = elecs.size();
100      int nMuons = muons.size();
101
102      if (nElecs + nMuons != 2) vetoEvent; // dilepton cut
103      if (_mode == 2 && nMuons !=2) vetoEvent;
104      if (_mode == 1 && nElecs !=2) vetoEvent;
105
106      if (elecs.size()==2){
107        l1=&elecs[0];
108        l2=&elecs[1];
109      }
110      else if (muons.size()==2){
111        l1=&muons[0];
112        l2=&muons[1];
113      }
114      else vetoEvent;
115
116      // if not e+e- or mu+mu- pair, veto
117      if (l1->pid() + l2->pid() !=0) vetoEvent;
118
119      // Dilepton selection :: Z mass peak.
120      FourMomentum ll = (l1->mom() + l2->mom());
121      double Zm  = ll.mass(); 
122      if ( !inRange(Zm/GeV, 71.0, 111.0) ) vetoEvent;
123      double ZpT = ll.pT();
124      
125      // Calculate the observables
126      double jet0pT = 0.;
127      double cljetpT = 0.;
128      double minDR = 99.;
129
130      // Require jets to be above 100GeV
131      ifilter_select(jets, Cuts::pT > 100*GeV);
132      double HTjet = sum(jets, Kin::pT, 0.);
133      for (const Jet& j : jets) {
134        // find minDR and closest jet to Z boson, only with 100GeV+ jets
135        double dr = deltaR(j, ll ,RAPIDITY);
136        if (dr < minDR) {
137          minDR = dr;
138          cljetpT = j.pT();
139        }
140      }
141      const size_t Njets = jets.size();
142
143      // NJets >= 1      
144      if (Njets < 1) vetoEvent;
145      // Exclusive NJets, jet pT > 100 GeV
146
147      _h["NJets"]->fill(Njets);
148
149      // Z pT
150      _h["ZpT"]->fill(ZpT/GeV);
151      
152      //Leading jet pT 
153      jet0pT = jets[0].pT();
154      _h["jetpT"]->fill(jet0pT/GeV);
155
156      // HT
157      double HT = HTjet + l1->pT() + l2->pT();
158      _h["HT"]->fill(HT/GeV);
159
160      // HTJet > 600 GeV selection
161      if (HTjet >= 600.) {
162        double minDR_HTjet = 99.;
163        // closest jet to Z
164        for (const Jet& j : jets) {
165          const double dr = deltaR(j, ll, RAPIDITY);
166          if (dr < minDR_HTjet)  minDR_HTjet = dr;
167        }
168        // Fill histograms of HTjet > 600 GeV
169        _h["NJets600"]->fill(Njets);
170        _h["minDR600"]->fill(minDR_HTjet);
171      } // end of HTjet > 600 GeV
172
173      // Our high pT phase-space
174      if (jet0pT/GeV < 500.0) vetoEvent;
175      
176      // Exclusive NJets, jet pT > 500 GeV
177      _h["NJets500"]->fill(Njets);
178      
179      // Z/J pT ratio
180      _h["rZJ"]->fill(ZpT/cljetpT);
181      _h["minDR"]->fill(minDR);
182      
183      // Phase space with DR<1.4       
184      if (minDR < 1.4) {
185        _h["NJets_coll"]->fill(Njets);
186        _h["rZJ_coll"]->fill(ZpT/cljetpT);
187      // Phase space with DR>2.0
188      } else if (minDR >2.0) {
189        _h["NJets_b2b"]->fill(Njets);
190        _h["rZJ_b2b"]->fill(ZpT/cljetpT);
191      }
192
193    }
194
195    void finalize() {
196      double xsec = crossSectionPerEvent()/picobarn;
197      // Analysis measures Z->ll(ee or mumu)
198      // If file contains both Z->ee and Z->mumu, divide xs by 2
199      if (_mode == 0) xsec *= 0.5;
200      // Normalize to cross section.
201      scale(_h, xsec);
202    } // end of finalize
203
204    //@}
205
206    // define histograms
207    size_t _mode;
208    map<string, Histo1DPtr> _h;
209 
210  };
211
212  RIVET_DECLARE_PLUGIN(ATLAS_2022_I2077570);
213}