Rivet analyses


title: ATLAS_2017_I1627873

EW Zjj using early Run-2 data

Experiment: ATLAS (LHC)

Inspire ID: 1627873

Status: VALIDATED

Authors: - Christian Johnson - Deepak Kar - Christian Gutschow

References: - Expt page: ATLAS-STDM-2016-09 - JHEP 1404 (2014) 031 - arXiv: 1709.10264

Beams: p+ p+

Beam energies: (6500.0, 6500.0)GeV

Run details: - Inclusive Z+jets events at 13 TeV, with Z decaying to muons or electrons

The cross-section for the production of two jets in association with a leptonically decaying Z boson ($Zjj$) is measured in proton-proton collisions at a centre-of-mass energy of 13 TeV, using data recorded with the ATLAS detector at the Large Hadron Collider, corresponding to an integrated luminosity of 3.2 fb$^{-1}$. The electroweak $Zjj$ cross-section is extracted in a fiducial region chosen to enhance the electroweak contribution relative to the dominant Drell-Yan $Zjj$ process, which is constrained using a data-driven approach. The measured fiducial electroweak cross-section is $\sigma_{EWZjj}$ = 119$\pm$16(stat.)$\pm$20(syst.)$\pm$2(lumi.) fb for dijet invariant mass greater than 250 GeV, and 34.2$\pm$5.8(stat.)$\pm$5.5(syst.)$\pm$0.7(lumi.) fb for dijet invariant mass greater than 1 TeV. Standard Model predictions are in agreement with the measurements. The inclusive $Zjj$ cross-section is also measured in six different fiducial regions with varying contributions from electroweak and Drell-Yan $Zjj$ production.

Source code:ATLAS_2017_I1627873.cc

```c++ // -- C++ --

include "Rivet/Analysis.hh"

include "Rivet/Projections/FastJets.hh"

include "Rivet/Projections/FinalState.hh"

include "Rivet/Projections/LeptonFinder.hh"

include "Rivet/Projections/PromptFinalState.hh"

include "Rivet/Projections/VetoedFinalState.hh"

namespace Rivet {

/// @brief EW Zjj using early Run-2 data class ATLAS_2017_I1627873 : public Analysis { public:

/// Constructor
RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2017_I1627873);


/// Book histograms and initialise projections before the run
void init() {

  _mode = 0;
  if (getOption("TYPE") == "EW_ONLY") _mode = 1;

  FinalState fs(Cuts::abseta < 5.0);

  FinalState photon_fs(Cuts::abspid == PID::PHOTON);

  PromptFinalState electron_fs(Cuts::abspid == PID::ELECTRON);
  PromptFinalState muon_fs(Cuts::abspid == PID::MUON);

  LeptonFinder dressed_electrons(electron_fs, photon_fs, 0.1, Cuts::abseta < 2.47 && Cuts::pT > 25 * GeV);
  declare(dressed_electrons, "DressedElectrons");

  LeptonFinder dressed_muons(muon_fs, photon_fs, 0.1, Cuts::abseta < 2.47 && Cuts::pT > 25 * GeV);
  declare(dressed_muons, "DressedMuons");

  VetoedFinalState remfs(fs);
  remfs.addVetoOnThisFinalState(dressed_electrons);
  remfs.addVetoOnThisFinalState(dressed_muons);

  FastJets jets(remfs, JetAlg::ANTIKT, 0.4, JetMuons::ALL, JetInvisibles::ALL);
  declare(jets, "Jets");

  if (_mode)
    book(_h, 3, 1, 1);
  else
    book(_h, 2, 1, 1);
}


/// Perform the per-event analysis
void analyze(const Event& event) {

  const Jets& jets = apply<FastJets>(event, "Jets").jetsByPt(Cuts::pT > 25 * GeV && Cuts::abseta < 4.4);
  DressedLeptons electrons = apply<LeptonFinder>(event, "DressedElectrons").dressedLeptons();
  DressedLeptons muons = apply<LeptonFinder>(event, "DressedMuons").dressedLeptons();

  // Overlap Removal
  idiscardIfAnyDeltaRLess(electrons, jets, 0.4);
  idiscardIfAnyDeltaRLess(muons, jets, 0.4);

  Particle lep1, lep2;
  if (electrons.size() == 2 && muons.empty()) {
    lep1 = electrons[0];
    lep2 = electrons[1];
    if (lep1.charge3() == lep2.charge3()) vetoEvent;
  }
  else if (electrons.empty() && muons.size() == 2) {
    lep1 = muons[0];
    lep2 = muons[1];
    if (lep1.charge3() == lep2.charge3()) vetoEvent;
  }
  else
    vetoEvent;

  if (jets.size() < 2) vetoEvent;

  const FourMomentum dilepton = lep1.mom() + lep2.mom();
  if (!inRange(dilepton.mass(), 81.0 * GeV, 101.0 * GeV)) vetoEvent;

  const double jet1pt = jets[0].pT();
  const double jet2pt = jets[1].pT();
  const double mjj = (jets[0].mom() + jets[1].mom()).mass();
  const double zpt = (lep1.mom() + lep2.mom()).pT();

  size_t ngapjets = 0;
  Jet thirdjet;
  for (size_t i = 2; i < jets.size(); ++i) {
    const Jet j = jets[i];
    if (_isBetween(j, jets[0], jets[1])) {
      if (!ngapjets) thirdjet = j;
      ++ngapjets;
    }
  }

  const double ptbal_vec = (jets[0].mom() + jets[1].mom() + lep1.mom() + lep2.mom()).pT();
  const double ptbal_sc = jets[0].pT() + jets[1].pT() + lep1.pT() + lep2.pT();
  const double ptbalance2 = ptbal_vec / ptbal_sc;

  const double ptbal3_vec = (jets[0].mom() + jets[1].mom() + thirdjet.mom() + lep1.mom() + lep2.mom())
                                .pT();
  const double ptbal3_sc = jets[0].pT() + jets[1].pT() + thirdjet.pT() + lep1.pT() + lep2.pT();
  const double ptbalance3 = ptbal3_vec / ptbal3_sc;


  //categories: baseline, high-PT, EW-enriched, QCD-enriched, high-mass, EW-enriched and high-mass
  if (!(jet1pt > 55 * GeV && jet2pt > 45 * GeV)) vetoEvent;

  if (_mode) {
    if (zpt > 20.0 * GeV && !ngapjets && ptbalance2 < 0.15 && mjj > 250.0 * GeV) {
      _h->fill(_h->bin(1).xEdge());
    }
    if (zpt > 20.0 * GeV && !ngapjets && ptbalance2 < 0.15 && mjj > 1000.0 * GeV) {
      _h->fill(_h->bin(2).xEdge());
    }
  }
  else {
    _h->fill(_h->bin(1).xEdge());
    if (jet1pt > 85.0 * GeV && jet2pt > 75.0 * GeV) _h->fill(_h->bin(2).xEdge());
    if (zpt > 20.0 * GeV && ngapjets == 0 && ptbalance2 < 0.15 && mjj > 250.0 * GeV)
      _h->fill(_h->bin(3).xEdge());
    if (zpt > 20.0 * GeV && ngapjets && ptbalance3 < 0.15 && mjj > 250.0 * GeV)
      _h->fill(_h->bin(4).xEdge());
    if (mjj > 1000.0 * GeV) _h->fill(_h->bin(5).xEdge());
    if (zpt > 20.0 * GeV && !ngapjets && ptbalance2 < 0.15 && mjj > 1000.0 * GeV)
      _h->fill(_h->bin(3).xEdge());
  }
}


/// Normalise histograms etc., after the run
void finalize() {

  const double factor = crossSection() / (_mode ? femtobarn : picobarn) / sumOfWeights();
  scale(_h, factor);
}

bool _isBetween(const Jet probe, const Jet boundary1, const Jet boundary2) {
  double y_p = probe.rapidity();
  double y_b1 = boundary1.rapidity();
  double y_b2 = boundary2.rapidity();

  double y_min = std::min(y_b1, y_b2);
  double y_max = std::max(y_b1, y_b2);

  if (y_p > y_min && y_p < y_max)
    return true;
  else
    return false;
}

///@}

private:

size_t _mode;

/// @name Histograms
///@{
BinnedHistoPtr<string> _h;
///@}

};

RIVET_DECLARE_PLUGIN(ATLAS_2017_I1627873);

} ```