Rivet analyses

title: LHCB_2025_I2922449

Measurement of b-jet mass with and without grooming

Experiment: LHCB (LHC)

Inspire ID: 2922449

Status: VALIDATED

Authors: - Ezra D. Lesser

References: - Phys. Lett. B 869 (2025) 139854 - DOI:10.1016/J.PhysLetB.2025.139854 - arXiv: 2505.11955 - Expt page: LHCb-PAPER-2025-009 - CERN-EP-2025-097

Beams: p+ p+

Beam energies: (6500.0, 6500.0)GeV

Run details: - Hard QCD events from proton-proton interactions at 13 TeV centre-of-mass energy in which $B^\pm$ meson decay is suppressed since they are fully reconstructed from the $J/\psi K^\pm$ decay channel before the b-jet mass measurement is performed.

The LHCb Collaboration presents measurements of the $b$-jet mass both with and without grooming using proton-proton data collected at a center-of-mass energy of $\sqrt{s} = 13$ TeV. $B^\pm$ mesons are fully reconstructed from the $J/\psi K^\pm$ decay channel. Jets with transverse momentum $p_\mathrm{T,jet} > 10$ GeV are reconstructed within the rapidity range $2.5 < y_\mathrm{jet} < 4.0$ using the anti-$k_\mathrm{T}$ algorithm with resolution paparmeter $R = 0.5$. Groomed jets are produced using Soft Drop with $z_\mathrm{cut} = 0.1$ and $\beta = 0$. Both groomed and ungroomed jets are further evaluated using the Winner-Take-All flavor algorithm, and figures are produced both with and without WTA flavor tagging.

Source code:LHCB_2025_I2922449.cc

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

include "Rivet/Analysis.hh"

include "Rivet/Analyses/LHCbCommon.hh"

include "Rivet/Projections/FastJets.hh"

include "fastjet/ClusterSequence.hh"

include "fastjet/JetDefinition.hh"

include "fastjet/contrib/LundGenerator.hh" // For Soft Drop grooming

namespace Rivet {

/// @brief Measurements of the $b$-jet mass both with and without grooming in proton-proton collision events at 13 TeV class LHCB_2025_I2922449 : public Analysis { ///////////////////////////////////////////////////////////////////////////////////////////////// // Analysis parameters

// Particle reconstruction
const double ETA_MIN_PARTICLES = 0.; // Min pseudorapidity of final-state particles
const double ETA_MAX_PARTICLES = 7.; // Max "

// Jet reconstruction
const double JET_R = 0.5;             // Jet resolution parameter for anti-kT
const double RAPIDITY_MIN_JETS = 2.5; // Min rapidity of constructed B-jets
const double RAPIDITY_MAX_JETS = 4.0; // Max "
const double PT_MIN_JETS = 10.;       // Min jet transverse momentum
const double PT_MAX_JETS = 100.;      // Max "

// HF reconstruction
const vector<PdgId> CH_B_MESON_PID = {PID::BPLUS, PID::BMINUS};

// WTA reclustering
const double SMALL_NUMBER = 1e-4;
const double& MAX_WTA_DISTANCE = SMALL_NUMBER;
const fastjet::JetDefinition WTA_JET_DEF = fastjet::JetDefinition(fastjet::cambridge_algorithm,
                                                                  fastjet::JetDefinition::max_allowable_R,
                                                                  fastjet::WTA_pt_scheme);

// Soft Drop grooming
const double SD_ZCUT = 0.1;
const double SD_BETA = 0.;
const fastjet::contrib::LundGenerator B_JET_LG = fastjet::contrib::LundGenerator(
    fastjet::cambridge_algorithm);

// Jet pT bin edges
const vector<double> PT_BIN_EDGES = {10., 12., 15., 20., 30., 50., 100.};

public:

/// Constructor
RIVET_DEFAULT_ANALYSIS_CTOR(LHCB_2025_I2922449);


/// @name Analysis methods
/// @{
/////////////////////////////////////////////////////////////////////////////////////////////////
/// Book histograms and initialise projections before the run
void init() {

  // Initialise and register projections
  // HVRecoFinalState projection selects also targeted charged B mesons
  Cut particle_selector = Cuts::etaIn(ETA_MIN_PARTICLES, ETA_MAX_PARTICLES);
  const LHCb::HVRecoFinalState hvfs(particle_selector, CH_B_MESON_PID, particle_selector);
  declare(hvfs, "hvfs");

  // The final-state particles declared above are clustered using FastJet with
  // the anti-kT algorithm and a jet-radius parameter 0.5
  // muons are included and neutrinos are excluded from the clustering
  FastJets jetfs(hvfs, JetAlg::ANTIKT, JET_R, JetMuons::ALL, JetInvisibles::NONE);
  declare(jetfs, "jets");

  // Book histograms

  book(_hh, PT_BIN_EDGES);        // HEPData Table 12-17
  book(_hh_gr, PT_BIN_EDGES);     // HEPData Table 0-5
  book(_hh_WTA, PT_BIN_EDGES);    // HEPData Table 36-41
  book(_hh_gr_WTA, PT_BIN_EDGES); // HEPData Table 18-23
  for (size_t i = 0; i < _hh->numBins(); ++i) {
    book(_hh->bin(i + 1), 13 + i, 1, 1);
    book(_hh_gr->bin(i + 1), 1 + i, 1, 1);
    book(_hh_WTA->bin(i + 1), 37 + i, 1, 1);
    book(_hh_gr_WTA->bin(i + 1), 19 + i, 1, 1);
  }

  return;
}


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

  // retrieve reconstructed charged B mesons from the HVRecoFinalState projection
  Particles B_mesons = apply<LHCb::HVRecoFinalState>(event, "hvfs").hfHadrons();

  if (B_mesons.empty()) {
    MSG_DEBUG("No charged B mesons found. Veto event...");
    vetoEvent;
  }

  // listing B mesons in DEBUG mode
  if (getLog().isActive(Log::Level::DEBUG)) {
    MSG_DEBUG("B meson(s) found with ");
    bool add_comma = false;
    for (const Particle& B_meson : B_mesons) {
      if (add_comma) MSG_DEBUG(", ");
      MSG_DEBUG("rapidity=" << B_meson.rap() << " and pT=" << B_meson.pt());
      add_comma = true;
    }
    MSG_DEBUG(std::endl);
  }

  // Retrieve clustered jets, sorted by pT, with applied rapidity and pT cuts
  Cut jet_selector = Cuts::rapIn(RAPIDITY_MIN_JETS, RAPIDITY_MAX_JETS)
      & Cuts::ptIn(PT_MIN_JETS, PT_MAX_JETS);
  Jets jets = apply<FastJets>(event, "jets").jetsByPt(jet_selector);
  if (jets.empty()) {
    MSG_DEBUG("Could not reconstruct jets in LHCb detector fiducial phase space" << std::endl);
    vetoEvent;
  }

  Jets B_jets;
  for (const Jet& jet : jets) {
    // Jet contains undecayed B+ or B-
    if (jet.containsParticleId(CH_B_MESON_PID)) B_jets.push_back(jet);
  }

  MSG_DEBUG("#" << B_jets.size() << " B jets found" << endl);
  if (B_jets.empty()) {
    vetoEvent;
  }
  ////////////////////////////////////////////////////////////////////////////////////////

  // Loop over jets and apply operations & fill histograms
  Jets B_jets_WTA;
  for (const Jet& B_jet : B_jets) {
    // Find the WTA axis of the B jet using C/A algorithm
    fastjet::ClusterSequence B_jet_WTA_cs(B_jet.pseudojet().constituents(), WTA_JET_DEF);
    Jet B_jet_WTA = fastjet::sorted_by_pt(B_jet_WTA_cs.inclusive_jets())[0];
    assert(B_jet.pseudojet().constituents().size()
           == B_jet_WTA.pseudojet()
                  .constituents()
                  .size()); // Make sure that reclustering preserves all particles

    // Check if charged B meson is on WTA axis (arXiv:2205.01117)
    bool WTA_tagged = false;
    for (const Particle& B_meson : B_mesons) {
      double WTA_distance = B_jet_WTA.pseudojet().delta_R(B_meson);
      if (WTA_distance < MAX_WTA_DISTANCE) {
        WTA_tagged = true;
        break;
      }
    }

    // Groom the jets with Soft Drop
    Jet B_jet_gr;
    for (fastjet::contrib::LundDeclustering& ld : B_JET_LG.result(B_jet)) {
      if (ld.z() > SD_ZCUT * std::pow(ld.Delta() / JET_R, SD_BETA)) { // SD condition
        B_jet_gr = ld.pair();
        break;
      }
    }
    // Check that charged B meson survives
    bool gr_tagged = false;
    for (const Particle& B_meson : B_mesons) {
      if (gr_tagged || !B_jet_gr.pseudojet().has_constituents()) break;
      for (const fastjet::PseudoJet& particle : B_jet_gr.pseudojet().constituents()) {
        if ((particle.delta_R(B_meson) < SMALL_NUMBER)
            && fuzzyEquals(particle.pt(), B_meson.pt(), SMALL_NUMBER)) {
          // Should be the same particle (could always compare barcodes?!)
          gr_tagged = true;
          break;
        }
      }
    }

    if (gr_tagged) {
      MSG_DEBUG("** B jet survived grooming" << endl);
    }
    else {
      MSG_DEBUG("** B jet did NOT survive grooming" << endl);
    }

    // Fill histograms
    MSG_DEBUG("Filling histograms..." << endl);
    double jpt = B_jet.pT() / GeV;
    double mpt = B_jet.pseudojet().m() / jpt;
    _hh->fill(jpt, mpt);
    if (WTA_tagged) _hh_WTA->fill(jpt, mpt);
    if (gr_tagged) {
      // Note: Groomed observable uses ungroomed jet pT
      double grmpt = B_jet_gr.pseudojet().m() / jpt;
      _hh_gr->fill(jpt, grmpt);
      if (WTA_tagged) _hh_gr_WTA->fill(jpt, grmpt);
    }
  }

  return;
}


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

  normalize(_hh_WTA);
  normalize(_hh_gr);
  normalize(_hh);
  normalize(_hh_gr_WTA);

  return;
}

/// @}

/////////////////////////////////////////////////////////////////////////////////////////////////
/// @name Histograms
/// @{
Histo1DGroupPtr _hh, _hh_gr, _hh_WTA, _hh_gr_WTA;
/// @}

};

RIVET_DECLARE_PLUGIN(LHCB_2025_I2922449);

} ```