fastjet

name

fastjet

description

Use when the user needs to run jet clustering in Python — anti-kt, kt, Cambridge/Aachen algorithms — from arrays of particle constituents (calorimeter towers, particle-flow objects, truth particles). Backed by the `fastjet` Python bindings. Does NOT cover reading pre-clustered jets from a file (use `physlite-basics` for PHYSLITE small-R jets), 4-vector arithmetic on existing jets (use `vector`), or jet energy calibration. Disambiguator phrase: fastjet Python jet clustering.

fastjet

Overview

fastjet is the standard jet-finding library in HEP. The Python bindings (fastjet PyPI package, also known as pyjet in older versions) expose the anti-kt, Cambridge/Aachen, and kt clustering algorithms. In ATLAS physics analyses the main use is truth-level jet clustering for generator studies or jet substructure calculations on particle-level events.

The library provides two interfaces: the array-oriented (columnar) interface, which accepts awkward arrays directly and handles many events at once without Python loops, and the classic (object-oriented) interface, which mirrors the C++ API and operates one event at a time on lists of PseudoJet objects. Prefer the columnar interface whenever you are iterating over events.

When to Use

• Clustering particle-level or parton-level events from a Monte Carlo generator
• Computing jet substructure variables (N-subjettiness, soft drop) for truth jets
• Validating detector-level jets against truth clustering
• Grooming jets for boosted object analyses

Key Concepts

Concept	Notes
fastjet.ClusterSequence	Dispatch point: pass ak.Array → columnar; pass list[PseudoJet] → classic
AwkwardClusterSequence	Internal type created for columnar input; inclusive_jets() returns ak.Array
fastjet.PseudoJet	Four-vector container (px, py, pz, E) — classic interface only
fastjet.JetDefinition	Algorithm + R parameter, e.g. antikt_algorithm with R=0.4
fastjet.ClusterSequenceArea	Adds jet area computation (needed for pileup subtraction) — classic only
fastjet.Selector	Filter jets by pT, eta, etc. — classic interface
Jet constituents (columnar)	cluster.constituents() returns jagged ak.Array of particles per jet
Jet constituents (classic)	jet.constituents() returns list of PseudoJets
User index	pseudo_jet.set_user_index(i) links back to the original particle

Canonical Patterns

Cluster jets from a list of particles

import fastjet as fj
import numpy as np

# Four-vectors: (px, py, pz, E) per particle
px = np.array([...])
py = np.array([...])
pz = np.array([...])
E  = np.array([...])

pseudojets = [fj.PseudoJet(px[i], py[i], pz[i], E[i]) for i in range(len(px))]

jet_def = fj.JetDefinition(fj.antikt_algorithm, 0.4)
cs = fj.ClusterSequence(pseudojets, jet_def)

jets = fj.sorted_by_pt(cs.inclusive_jets(ptmin=25.0))   # 25 GeV (use 25000.0 MeV for ATLAS)

for j in jets:
    print(f"pT={j.pt():.1f} GeV  eta={j.eta():.2f}  phi={j.phi():.2f}")
    print(f"  {len(j.constituents())} constituents")

Columnar interface — single event (px/py/pz/E fields)

Pass an ak.Array directly; no PseudoJet construction loop needed. Output is also an ak.Array (MomentumArray4D), not a list.

import fastjet
import awkward as ak

particles = ak.Array([
    {"px": 1.2,  "py": 3.2,  "pz":   5.4, "E":  23.5},
    {"px": 32.2, "py": 64.21, "pz": 543.34, "E": 755.12},
    {"px": 32.45, "py": 63.21, "pz": 543.14, "E": 835.56},
])

jetdef = fastjet.JetDefinition(fastjet.antikt_algorithm, 0.4)
cluster = fastjet.ClusterSequence(particles, jetdef)

jets = cluster.inclusive_jets(min_pt=25.0)  # MomentumArray4D; 25 GeV (use 25000.0 MeV for ATLAS)
print(jets.pt, jets.eta, jets.phi)

Columnar interface — multi-event (jagged awkward array)

The columnar interface is designed for many events at once. Pass a jagged array (outer dim = events, inner dim = particles per event); output is also jagged.

import fastjet
import awkward as ak

# Each sub-array is one event; lengths may differ (jagged)
events = ak.Array([
    [{"px": 1.2,  "py": 3.2,  "pz":   5.4, "E":  23.5},
     {"px": 32.2, "py": 64.21, "pz": 543.34, "E": 755.12}],
    [{"px": 32.45, "py": 63.21, "pz": 543.14, "E": 835.56},
     {"px": -1.2, "py": -3.2, "pz": -5.4, "E": 23.5}],
])

jetdef = fastjet.JetDefinition(fastjet.antikt_algorithm, 0.4)
cluster = fastjet.ClusterSequence(events, jetdef)

jets = cluster.inclusive_jets(min_pt=25.0)  # jagged: [[jets_ev0], [jets_ev1]]; 25 GeV (use 25000.0 MeV for ATLAS)
print(ak.num(jets, axis=1))   # number of jets per event
print(jets.pt)                # pT of all jets, jagged

Columnar interface — pt/eta/phi/M coordinates with vector

import fastjet
import awkward as ak
import vector
vector.register_awkward()

events = ak.Array(
    [
        [{"pt": 71.8, "eta": 2.72, "phi": 1.11, "M": 0.14},
         {"pt": 3.42, "eta": 1.24, "phi": 1.21, "M": 0.0}],
        [{"pt": 55.0, "eta": -1.5, "phi": 0.5,  "M": 0.14}],
    ],
    with_name="Momentum4D",
)

jetdef = fastjet.JetDefinition(fastjet.antikt_algorithm, 0.4)
cluster = fastjet.ClusterSequence(events, jetdef)
jets = cluster.inclusive_jets(min_pt=25.0)  # 25 GeV (use 25000.0 MeV for ATLAS)

Columnar interface — accessing constituents without a Python loop

# constituents() returns a jagged array: [event][jet] → ak.Array of particles
constituents = cluster.constituents(min_pt=25.0)  # 25 GeV (use 25000.0 MeV for ATLAS)
# constituents[event_i][jet_j] is an ak.Array of the contributing particles
# e.g. sum constituent pT:
print(ak.sum(constituents.pt, axis=-1))

Jet substructure: N-subjettiness

import fastjet as fj

# requires fastjet contrib (fjcontrib)
from fastjet import contrib

tau_calc = contrib.Nsubjettiness(1, contrib.OnePass_KT_Axes(),
                                  contrib.UnnormalizedMeasure(beta=1.0))

for jet in jets:
    tau1 = tau_calc(jet)

Soft drop grooming

from fastjet import contrib

sd = contrib.SoftDrop(beta=0.0, zcut=0.1, R=1.0)
groomed_jets = [sd(j) for j in jets]

Jet area (for pileup)

ghost_area = fj.AreaDefinition(fj.active_area, fj.GhostedAreaSpec(5.0))
cs_area = fj.ClusterSequenceArea(pseudojets, jet_def, ghost_area)
for j in cs_area.inclusive_jets(ptmin=25.0):  # 25 GeV (use 25000.0 MeV for ATLAS)
    print(f"area = {j.area():.3f}")

Gotchas

• Prefer the columnar interface over per-particle loops: Building PseudoJet objects one-by-one in a Python loop ([fj.PseudoJet(p.px, ...) for p in particles]) crosses the Python interpreter on every particle and every event. Pass an ak.Array to ClusterSequence directly — this pushes the conversion into C++ and processes all events without Python overhead.
• Classic interface is required for ClusterSequenceArea and fjcontrib: The columnar interface supports only ClusterSequence. Jet area computation and substructure tools (Nsubjettiness, SoftDrop) require the classic PseudoJet-based path.
• Units: fastjet has no built-in unit system — ensure all four-vectors use the same units (usually GeV).
• phi() range: fastjet returns phi in [0, 2π); some downstream code expects (-π, π]. Use fj.PseudoJet.phi_std() or subtract 2π if needed.
• pyjet vs fastjet: the older pyjet package uses numpy structured arrays; the newer fastjet package (from Scikit-HEP) uses PseudoJet objects and is the recommended path.
• fjcontrib availability: substructure tools (Nsubjettiness, SoftDrop) require the fastjet package built with contrib support, or a separate fjcontrib install.
• with_name="Momentum4D" required for pt/eta/phi/M input: Without vector.register_awkward() and with_name="Momentum4D", the columnar interface only auto-detects (px, py, pz, E) field names.

Interop

• pyhepmc: Read HepMC3 events, extract final-state particles into an ak.Array with (px, py, pz, E) fields, then pass directly to ClusterSequence using the columnar interface.
• pylhe: Read LHE parton-level events for parton-jet matching studies.
• vector / awkward: Use vector.register_awkward() + with_name="Momentum4D" to pass (pt, eta, phi, M) arrays to the columnar interface. Classic interface still requires converting to PseudoJet objects first.

Docs

https://fastjet.fr/ Python bindings: https://scikit-hep.org/fastjet/