Rivet analyses

title: CMS_2019_I1764472

Measurement of the differential ttbar production cross section as a function of the jet mass and top quark mass in boosted hadronic top quark decays.

Experiment: CMS (LHC)

Inspire ID: 1764472

Status: VALIDATED

Authors: - Dennis Schwarz - Roman Kogler - Johannes Haller

References: - TOP-19-005 - Phys.Rev.Lett. 124 (2020) 20, 202001 - arXiv: 1911.03800

Beams: p+ p+

Beam energies: (6500.0, 6500.0)GeV

Run details: - ttbar events at sqrt(s) = 13 TeV, lepton+jets selection at particle level

A measurement of the ttbar production cross section as a function of the jet mass of hadronic decays of boosted top quarks is presented. The measurement is carried out in the lepton+jets channel. As leptons, muons and electrons are defined with originate from the W boson decay. Jets are clustered from all stable particles excluding neutrinos in a two-step procedure using the XCone algorithm. At first, two large jets (R=1.2) are found aiming at a reconstruction of the two top quarks. Using the constituents of those jets, XCone is run again finding three subjets with R=0.4. A jet representing the hadronic and leptonic top quark decay is found by the angular distance to the leading lepton in the event. The final jets are constructed as the sum of the subjet four-momenta. The jet identified as the one containing the hadronic top quark decay is required to have \pt > 400 GeV. Furthermore, the jet mass of this jet has to be larger than the invariant mass of the combined system of second jet and lepton. Both, the differential and normalized differential ttbar production cross sections are measured as a function of the jet mass.

Source code:CMS_2019_I1764472.cc

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

include "Rivet/Analysis.hh"

include "Rivet/Projections/ChargedLeptons.hh"

include "Rivet/Projections/FastJets.hh"

include "Rivet/Projections/FinalState.hh"

include "Rivet/Projections/LeptonFinder.hh"

include "Rivet/Projections/PartonicTops.hh"

include "Rivet/Projections/PromptFinalState.hh"

include "Rivet/Projections/VetoedFinalState.hh"

include "fastjet/contrib/Nsubjettiness.hh"

include "fastjet/contrib/XConePlugin.hh"

namespace Rivet {

/// @brief Measurement of the jet mass for boosted top quarks at 13 TeV class CMS_2019_I1764472 : public Analysis { public:

/// Constructor
RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2019_I1764472);


/// @name Analysis methods
/// @{

void init() {

  // Prompt leptons
  ChargedLeptons charged_leptons;
  PromptFinalState prompt_leptons(charged_leptons);
  declare(prompt_leptons, "PromptLeptons");

  // Final state particles for jet clustering
  VetoedFinalState fs_jets;
  fs_jets.vetoNeutrinos();

  // First XCone jet clustering step
  /// @todo Clean this up with something in-place
  fastjet::contrib::PseudoXConePlugin* plugin_xcone = new fastjet::contrib::PseudoXConePlugin(2, 1.2,
                                                                                              2.0);
  declare(FastJets(fs_jets, plugin_xcone), "FatJets");

  // Partonic tops for decay channel definition
  declare(PartonicTops(TopDecay::E_MU, PromptEMuFromTau::NO), "LeptonicTops");
  declare(PartonicTops(TopDecay::HADRONIC), "HadronicTops");

  // Book histograms
  book(_hist_mass, "d01-x01-y01");
  book(_hist_mass_norm, "d02-x01-y01");
}


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

  // Decay mode check
  const Particles& leptonicTops = apply<PartonicTops>(event, "LeptonicTops").particlesByPt();
  const Particles& hadronicTops = apply<PartonicTops>(event, "HadronicTops").particlesByPt();
  if (leptonicTops.size() != 1 || hadronicTops.size() != 1) vetoEvent;


  // Get prompt leptons
  const PromptFinalState& prompt_leptons = apply<PromptFinalState>(event, "PromptLeptons");
  const Particles& leptons = prompt_leptons.particles();
  if (leptons.empty()) vetoEvent;

  // Select leading lepton
  Particle lepton;
  for (const Particle& l : leptons) {
    if (l.pT() > lepton.pT()) lepton = l;
  }
  if (lepton.pT() < 60 * GeV) vetoEvent;

  // Get the fat jets
  const Jets& fatjets = apply<FastJets>(event, "FatJets").jets();

  // Get index of hadronic jet by distance to lepton
  int ihad = 0;
  int ilep = 1;

  double dR0 = deltaR(lepton, fatjets.at(0));
  double dR1 = deltaR(lepton, fatjets.at(1));

  if (dR0 < dR1) {
    ihad = 1;
    ilep = 0;
  }

  // Get jet constituents
  const Particles& phad = fatjets.at(ihad).particles();
  const Particles& plep = fatjets.at(ilep).particles();

  // Cluster subjets
  FinalState fs_dummy;
  /// @todo Avoid rebuilding the plugin for every event
  fastjet::JetDefinition::Plugin* plugin_subhad = new fastjet::contrib::PseudoXConePlugin(3, 0.4, 2.0);
  fastjet::contrib::PseudoXConePlugin* plugin_sublep = new fastjet::contrib::PseudoXConePlugin(3, 0.4,
                                                                                               2.0);
  FastJets hadsubcluster(fs_dummy, plugin_subhad);
  FastJets lepsubcluster(fs_dummy, plugin_sublep);
  hadsubcluster.calc(phad);
  lepsubcluster.calc(plep);

  Jets subjets_had = hadsubcluster.jets();
  Jets subjets_lep = lepsubcluster.jets();

  // Subtract the lepton four vector from closest subjet if dR<0.4
  Jets subjets_had_clean;
  double dRmin_had = 0.4;
  unsigned int i_dRmin_had = 0;
  bool found_match_had = false;
  for (unsigned int i = 0; i < subjets_had.size(); i++) {
    double dR = deltaR(subjets_had[i], lepton);
    if (dR < dRmin_had) {
      dRmin_had = dR;
      i_dRmin_had = i;
      found_match_had = true;
    }
  }
  for (unsigned int i = 0; i < subjets_had.size(); i++) {
    Jet subjet = subjets_had[i];
    if (found_match_had && i == i_dRmin_had) {
      subjet = Jet(subjets_had[i].momentum() - lepton.momentum(), subjets_had[i].particles(),
                   subjets_had[i].tags());
    }
    subjets_had_clean.push_back(subjet);
  }
  std::sort(subjets_had_clean.begin(), subjets_had_clean.end(), cmpMomByPt);

  // do the same for lep jets
  Jets subjets_lep_clean;
  double dRmin_lep = 0.4;
  unsigned int i_dRmin_lep = 0;
  bool found_match_lep = false;
  for (unsigned int i = 0; i < subjets_lep.size(); i++) {
    double dR = deltaR(subjets_lep[i], lepton);
    if (dR < dRmin_lep) {
      dRmin_lep = dR;
      i_dRmin_lep = i;
      found_match_lep = true;
    }
  }
  for (unsigned int i = 0; i < subjets_lep.size(); i++) {
    Jet subjet = subjets_lep[i];
    if (found_match_lep && i == i_dRmin_lep) {
      subjet = Jet(subjets_lep[i].momentum() - lepton.momentum(), subjets_lep[i].particles(),
                   subjets_lep[i].tags());
    }
    subjets_lep_clean.push_back(subjet);
  }
  std::sort(subjets_lep_clean.begin(), subjets_lep_clean.end(), cmpMomByPt);

  // Subjet cuts
  if (subjets_had_clean.size() != 3) vetoEvent;
  if (subjets_lep_clean.size() != 3) vetoEvent;
  for (Jet jet : subjets_had_clean) {
    if (jet.pT() < 30 * GeV) vetoEvent;
    if (jet.abseta() > 2.5) vetoEvent;
  }

  // Combine subjets to final jets
  FourMomentum hadjet;
  for (Jet subjet : subjets_had_clean) {
    if (subjet.abseta() < 2.5) hadjet += subjet.momentum();
  }
  FourMomentum lepjet;
  for (Jet subjet : subjets_lep_clean) {
    if (subjet.abseta() < 2.5) lepjet += subjet.momentum();
  }

  // Jet pT cuts
  if (hadjet.pT() < 400 * GeV) vetoEvent;
  if (lepjet.pT() < 10 * GeV) vetoEvent;

  // m(hadjet) > m(lepjet+lepton)
  FourMomentum secondJetLepton = lepjet + lepton.momentum();
  if (hadjet.mass() < secondJetLepton.mass()) vetoEvent;

  // Fill histograms
  _hist_mass->fill(hadjet.mass() / GeV);
  _hist_mass_norm->fill(hadjet.mass() / GeV);
}


/// Normalise and scale histograms
void finalize() {
  const double sf = crossSection() / femtobarn / sumOfWeights();
  scale(_hist_mass, sf);
  normalize(_hist_mass_norm, 1.0, false);
}

/// @}

private:

// Histograms
Histo1DPtr _hist_mass, _hist_mass_norm;

};

RIVET_DECLARE_PLUGIN(CMS_2019_I1764472);

} ```