Rivet analyses

title: ATLAS_2016_I1424838

Event shapes in leptonic $Z$-events

Experiment: ATLAS (LHC)

Inspire ID: 1424838

Status: VALIDATED

Authors: - Holger Schulz

References: - Expt page: ATLAS-STDM-2014-07 - arXiv: 1602.08980 - Eur. Phys. J. C 76(7), 1-40, (2016)

Beams: p+ p+

Beam energies: (3500.0, 3500.0)GeV

Run details: - $Z \to \ell\ell$ with $p_\perp(\ell)>20$~GeV

Measurement of transverse event-shape observables ($N_\text{ch}$, $\sum \pT$, thrust, beam-thrust, $F$-parameter, and spherocity) in $Z \to \ell\ell$ events at 7~TeV for different $\pT(Z)$ regions.

Source code:ATLAS_2016_I1424838.cc

```c++

include "Rivet/Analysis.hh"

include "Rivet/Projections/Beam.hh"

include "Rivet/Projections/ChargedFinalState.hh"

include "Rivet/Projections/DileptonFinder.hh"

include "Rivet/Projections/FParameter.hh"

include "Rivet/Projections/Spherocity.hh"

include "Rivet/Projections/Thrust.hh"

include "Rivet/Projections/VetoedFinalState.hh"

namespace Rivet {

/// @brief Event shapes in leptonic $Z$-events class ATLAS_2016_I1424838 : public Analysis { public:

/// Constructor
RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2016_I1424838);


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

  // Charged particles inside acceptance region
  const ChargedFinalState cfs(Cuts::abseta < 2.5 && Cuts::pT > 500 * MeV);
  declare(cfs, "CFS");

  // Dilepton finders
  DileptonFinder zfinder(91.2 * GeV, 0.1,
                         Cuts::abseta < 2.4 && Cuts::pT > 20.0 * GeV && Cuts::abspid == PID::ELECTRON,
                         Cuts::massIn(66 * GeV, 116 * GeV));
  declare(zfinder, "DileptonFinder");
  DileptonFinder zfinder_mu(91.2 * GeV, 0.1,
                            Cuts::abseta < 2.4 && Cuts::pT > 20.0 * GeV && Cuts::abspid == PID::MUON,
                            Cuts::massIn(66 * GeV, 116 * GeV));
  declare(zfinder_mu, "DileptonFinderMu");

  // This CFS only contains charged particles inside the acceptance excluding the leptons
  VetoedFinalState remfs(cfs);
  remfs.addVetoOnThisFinalState(zfinder);
  remfs.addVetoOnThisFinalState(zfinder_mu);
  declare(remfs, "REMFS");

  const FParameter fparam(remfs);
  declare(fparam, "FParameter_");

  const Spherocity sphero(remfs);
  declare(sphero, "Spherocity_");


  // Booking of ES histos
  for (size_t alg = 0; alg < 5; ++alg) {
    // Book the inclusive histograms
    size_t offset = alg * 6;
    book(_h_Elec_Ntrk[alg], offset + 1, 1, 1);
    book(_h_Elec_SumPt[alg], offset + 2, 1, 1);
    book(_h_Elec_Beamthrust[alg], offset + 3, 1, 1);
    book(_h_Elec_Thrust[alg], offset + 4, 1, 1);
    book(_h_Elec_FParam[alg], offset + 5, 1, 1);
    book(_h_Elec_Spherocity[alg], offset + 6, 1, 1);
    book(_h_Muon_Ntrk[alg], offset + 31, 1, 1);
    book(_h_Muon_SumPt[alg], offset + 32, 1, 1);
    book(_h_Muon_Beamthrust[alg], offset + 33, 1, 1);
    book(_h_Muon_Thrust[alg], offset + 34, 1, 1);
    book(_h_Muon_FParam[alg], offset + 35, 1, 1);
    book(_h_Muon_Spherocity[alg], offset + 36, 1, 1);
  }
}


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

  // Check for Z boson in event
  const DileptonFinder& zfinder = apply<DileptonFinder>(event, "DileptonFinder");
  MSG_DEBUG("Num e+ e- pairs found = " << zfinder.bosons().size());
  const bool isElec = zfinder.bosons().size() == 1;

  const DileptonFinder& zfinder_mu = apply<DileptonFinder>(event, "DileptonFinderMu");
  MSG_DEBUG("Num mu+ mu- pairs found = " << zfinder_mu.bosons().size());
  const bool isMuon = zfinder_mu.bosons().size() == 1;

  // Only accept events with exactly two electrons or exactly two muons
  if (isElec && isMuon) vetoEvent;
  if (!(isElec || isMuon)) vetoEvent;

  // This determines the Zpt phase-space
  double zpT = -1000;
  if (isElec) zpT = zfinder.bosons()[0].pT();
  if (isMuon) zpT = zfinder_mu.bosons()[0].pT();

  unsigned int alg = 4; //< for > 25 GeV
  if (zpT < 6 * GeV)
    alg = 1;
  else if (inRange(zpT / GeV, 6, 12))
    alg = 2;
  else if (inRange(zpT / GeV, 12, 25))
    alg = 3;
  assert(alg < 5);
  assert(alg > 0);

  // All charged particles within |eta|<2.5 except the leptons from Z-decay
  const VetoedFinalState& remfs = apply<VetoedFinalState>(event, "REMFS");
  // sumPt and Beamthrust (the latter will only be filled if the min Nch criterion is met)
  // and Thrust preparation
  double sumPt = 0.0, beamThrust = 0.0;
  vector<Vector3> momenta;
  for (const Particle& p : remfs.particles()) {
    const double pT = p.pT();
    sumPt += pT;
    beamThrust += pT * exp(-p.abseta());
    const Vector3 mom = p.mom().pTvec();
    momenta.push_back(mom);
  }

  // Fill inclusive histos
  if (isElec) {
    _h_Elec_Ntrk[alg]->fill(remfs.size());
    _h_Elec_Ntrk[0]->fill(remfs.size());
    _h_Elec_SumPt[alg]->fill(sumPt);
    _h_Elec_SumPt[0]->fill(sumPt);
  }
  if (isMuon) {
    _h_Muon_Ntrk[alg]->fill(remfs.size());
    _h_Muon_Ntrk[0]->fill(remfs.size());
    _h_Muon_SumPt[alg]->fill(sumPt);
    _h_Muon_SumPt[0]->fill(sumPt);
  }

  // Skip event shape calculation if we don't match the minimum Nch criterion
  if (remfs.size() >= 2) {

    // Eventshape calculations

    // Calculate transverse Thrust using all charged FS particles except the lepton
    // This is copied/inspired from the CMS_6000011_S8957746 analysis
    if (momenta.size() == 2) {
      // We need to use a ghost so that Thrust.calc() doesn't return 1.
      momenta.push_back(Vector3(1e-10 * MeV, 0., 0.));
    }
    Thrust thrustC;
    thrustC.calc(momenta);

    double thrust = thrustC.thrust();

    // F-Parameter
    const FParameter& fparam = apply<FParameter>(event, "FParameter_");
    // Spherocity
    const Spherocity& sphero = apply<Spherocity>(event, "Spherocity_");

    // Histos differential in NMPI

    // Fill inclusive histos
    if (isElec) {
      _h_Elec_Thrust[alg]->fill(thrust);
      _h_Elec_Thrust[0]->fill(thrust);
      _h_Elec_FParam[alg]->fill(fparam.F());
      _h_Elec_FParam[0]->fill(fparam.F());
      _h_Elec_Spherocity[alg]->fill(sphero.spherocity());
      _h_Elec_Spherocity[0]->fill(sphero.spherocity());
      _h_Elec_Beamthrust[alg]->fill(beamThrust / GeV);
      _h_Elec_Beamthrust[0]->fill(beamThrust / GeV);
    }
    if (isMuon) {
      _h_Muon_Thrust[alg]->fill(thrust);
      _h_Muon_Thrust[0]->fill(thrust);
      _h_Muon_FParam[alg]->fill(fparam.F());
      _h_Muon_FParam[0]->fill(fparam.F());
      _h_Muon_Spherocity[alg]->fill(sphero.spherocity());
      _h_Muon_Spherocity[0]->fill(sphero.spherocity());
      _h_Muon_Beamthrust[alg]->fill(beamThrust / GeV);
      _h_Muon_Beamthrust[0]->fill(beamThrust / GeV);
    }
  }
}


/// Normalise histograms etc., after the run
void finalize() {
  for (size_t alg = 0; alg < 5; ++alg) {
    normalize(_h_Elec_Ntrk[alg], 100.);
    normalize(_h_Elec_SumPt[alg], 100.);
    normalize(_h_Elec_Beamthrust[alg], 100.);
    normalize(_h_Elec_Thrust[alg]);
    normalize(_h_Elec_FParam[alg]);
    normalize(_h_Elec_Spherocity[alg]);
    normalize(_h_Muon_Ntrk[alg], 100.);
    normalize(_h_Muon_SumPt[alg], 100.);
    normalize(_h_Muon_Beamthrust[alg], 100.);
    normalize(_h_Muon_Thrust[alg]);
    normalize(_h_Muon_FParam[alg]);
    normalize(_h_Muon_Spherocity[alg]);
  }
}

private:

Histo1DPtr _h_Elec_Ntrk[5];
Histo1DPtr _h_Elec_SumPt[5];
Histo1DPtr _h_Elec_Beamthrust[5];
Histo1DPtr _h_Elec_Thrust[5];
Histo1DPtr _h_Elec_FParam[5];
Histo1DPtr _h_Elec_Spherocity[5];

Histo1DPtr _h_Muon_Ntrk[5];
Histo1DPtr _h_Muon_SumPt[5];
Histo1DPtr _h_Muon_Beamthrust[5];
Histo1DPtr _h_Muon_Thrust[5];
Histo1DPtr _h_Muon_FParam[5];
Histo1DPtr _h_Muon_Spherocity[5];

};

RIVET_DECLARE_PLUGIN(ATLAS_2016_I1424838);

} ```