feat(analysis): per-step secondary multiplicity plots

Replace the event-level n_sec confusion matrix with two step-resolved
secondary-multiplicity comparisons:

- sec_count_per_step: overlay histogram of how many secondaries a single
  step emits, rollout series vs reference.
- sec_count_per_step_by_species: heatmap of per-step multiplicity of one
  species (zero row included) against species, drawn as one panel per
  rollout plus a reference panel, raw counts on a log color scale.

Both are backed by a new sources.secondaries_by_step view, which tags each
secondary with its emitting step — (event_id, parent_id, birth position)
on the rollout side, the row index on the reference side — so neither plot
needs a join against the step frame. Steps that emitted nothing are
recovered by subtraction from the chunk's step count, keeping both specs
sum-mergeable across condor chunks.

The rollout multiplicity is derived from the actual secondary birth rows
rather than the n_sec_pred column, which records the predicted count
before the per-event max-tracks cap.

_render_heatmap gained reference-panel and log-color support;
marginal_distance_summary sets neither key and is unchanged.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

tests/test_analysis_reduce.py

@@ -13,6 +13,7 @@ from giant.analysis.sources import (
     open_side,
     physical_steps,
     secondaries,
+    secondaries_by_step,
 )
 from giant.data.loader import EVENT_ID_FILE_STRIDE
 
@@ -159,18 +160,17 @@ def test_secondaries_rollout_vs_reference_align():
     assert t["pdg"].to_list() == [22, 22]
 
 
-def test_sec_count_by_event_zero_fills_events_with_no_secondaries():
-    r_phys = physical_steps(_rollout_frame(), Side.rollout)
-    r_sec = secondaries(_rollout_frame(), Side.rollout)
-    ev, n = R.sec_count_by_event(r_phys, r_sec)
-    # event 1 has one secondary track; event 2 has none and must still appear (as 0),
-    # not silently drop out of a plain group_by on the secondaries frame alone.
-    assert dict(zip(ev.tolist(), n.tolist())) == {1: 1, 2: 0}
+def test_secondaries_by_step_keys_each_secondary_to_its_emitting_step():
+    r = secondaries_by_step(_rollout_frame(), Side.rollout).collect()
+    assert r["pdg"].to_list() == [22]
+    # the rollout key is (event_id, parent_id, birth position) — the parent
+    # step's post_pos, copied verbatim onto the child's birth row.
+    assert r["step_key"][0] == {"event_id": 1, "parent_id": 0, "pre_x": 0.0, "pre_y": 0.0, "pre_z": 1.0}
 
-    t_all = _reference_frame()
-    t_sec = secondaries(t_all, Side.reference)
-    ev, n = R.sec_count_by_event(t_all, t_sec)
-    assert dict(zip(ev.tolist(), n.tolist())) == {1: 1, 2: 1}
+    t = secondaries_by_step(_reference_frame(), Side.reference).collect()
+    assert t["pdg"].to_list() == [22, 22]
+    # one row per emitting step; the empty-list step drops out entirely
+    assert [k["_row"] for k in t["step_key"]] == [0, 2]
 
 
 def test_leakage_fraction():

tests/test_catalog.py

@@ -10,10 +10,10 @@ from giant.analysis.catalog import (
     Bundle,
     PlotSpec,
     _containment_depths,
-    _integer_confusion,
     _ks_statistic,
 )
 from giant.analysis.context import Context, build_context
+from giant.analysis.grouping import pdg_label
 from giant.analysis.sources import RolloutSpec
 from tests.test_analysis_reduce import _reference_frame, _rollout_frame
 
@@ -160,8 +160,9 @@ def _validate_payload(r, names: list[str]) -> None:
 # data-dependent edges (event_total_edep), concat-then-mean/std (shower_
 # longitudinal), concat-then-max-edge (leakage_fraction), pdg-keyed sum with a
 # ratio (species_edep_share), a chunkable=False passthrough (router_gating),
-# nested sum-merge into a scorecard (marginal_distance_summary), concat-then-
-# event-id-join (n_sec_confusion), and concat-then-per-event-derived-quantity
+# sum-mergeable-with-a-zero-fill-denominator (sec_count_per_step{,_by_species}),
+# nested sum-merge into a scorecard (marginal_distance_summary), and
+# concat-then-per-event-derived-quantity
 # (shower_containment_depth_90, reusing the profile matrix's own merge shape).
 _CHUNK_EQUIVALENCE_IDS = [
     "marginal_edep",
@@ -170,9 +171,10 @@ _CHUNK_EQUIVALENCE_IDS = [
     "shower_longitudinal",
     "leakage_fraction",
     "sec_count_per_species",
+    "sec_count_per_step",
+    "sec_count_per_step_by_species",
     "router_gating",
     "marginal_distance_summary",
-    "n_sec_confusion",
     "shower_containment_depth_90",
 ]
 
@@ -217,7 +219,7 @@ def test_chunked_matches_unchunked(two_ctx: Context, spec_id: str):
 
 
 # ---------------------------------------------------------------------------
-# new (gitea #76) reductions: KS distance, confusion matrix, containment depth
+# new (gitea #76) reductions: KS distance and containment depth
 # ---------------------------------------------------------------------------
 
 
@@ -228,28 +230,6 @@ def test_ks_statistic():
     assert _ks_statistic([10, 0], [0, 0]) == 1.0  # one side empty, other isn't -> maximal mismatch
 
 
-def test_integer_confusion_matches_event_pairing():
-    # true (reference) n_sec = [1, 1]; predicted (rollout) n_sec = [1, 0]
-    labels, mat = _integer_confusion(np.array([1, 1]), np.array([1, 0]))
-    assert labels == ["0", "1+"]
-    assert mat.tolist() == [[0, 0], [1, 1]]  # row=true, col=pred
-
-
-def test_integer_confusion_caps_pathological_outliers():
-    labels, mat = _integer_confusion(np.array([0, 500]), np.array([0, 0]), max_bins=5)
-    assert labels[-1] == "4+"
-    assert mat.shape == (5, 5)
-    assert mat.sum() == 2
-
-
-def test_integer_confusion_explicit_cap_overrides_local_range():
-    # Even though this pair's own max is 1, an explicit shared cap forces a
-    # wider (and so cross-rollout-consistent) label set.
-    labels, mat = _integer_confusion(np.array([1, 1]), np.array([0, 1]), cap=3)
-    assert labels == ["0", "1", "2", "3+"]
-    assert mat.shape == (4, 4)
-
-
 def test_containment_depths_simple_ramp():
     # one event, edep concentrated in the first bin -> 90%/95% containment
     # depth is the first bin's right edge; a zero-energy event is dropped.
@@ -259,17 +239,25 @@ def test_containment_depths_simple_ramp():
     assert depths.tolist() == [1.0]
 
 
-def test_n_sec_confusion_spec(bundle):
-    spec = get_spec("n_sec_confusion")
+def test_sec_count_per_step_counts_empty_steps(bundle):
+    spec = get_spec("sec_count_per_step")
     r = spec.finalize([spec.compute_partial(bundle)], bundle.ctx)
-    assert r.payload["row_labels"] == r.payload["col_labels"] == ["0", "1+"]
-    assert r.payload["series"]["rollout"] == [[0, 0], [1, 1]]
+    # reference: 3 steps, two of which emit exactly one secondary
+    assert r.payload["reference"][:2] == [1, 2]
+    # rollout: 4 physical steps, one of which emits a single secondary
+    assert r.payload["series"]["rollout"][:2] == [3, 1]
+    assert sum(r.payload["reference"]) == 3
 
 
-def test_n_sec_confusion_shares_one_cap_across_rollouts(two_bundle):
-    spec = get_spec("n_sec_confusion")
-    r = spec.finalize([spec.compute_partial(two_bundle)], two_bundle.ctx)
-    assert list(r.payload["series"]) == ["flow", "wgan"]
-    # both rollouts share the same fixture data here, so their matrices (and
-    # the shared label set) must be identical.
-    assert r.payload["series"]["flow"] == r.payload["series"]["wgan"]
+def test_sec_count_per_step_by_species_zero_row_is_per_species(bundle):
+    spec = get_spec("sec_count_per_step_by_species")
+    r = spec.finalize([spec.compute_partial(bundle)], bundle.ctx)
+    cols = r.payload["col_labels"]
+    ref = r.payload["reference"]
+    g = cols.index(pdg_label(22))
+    # two reference steps emit one photon each; the third emits none
+    assert [row[g] for row in ref][:2] == [1, 2]
+    # every other species column is "no such secondary" on all 3 steps
+    for j, _ in enumerate(cols):
+        if j != g:
+            assert ref[0][j] == 3 and sum(row[j] for row in ref[1:]) == 0

tests/test_render.py

@@ -304,13 +304,15 @@ def test_render_one_of_each_kind(tmp_path: Path):
             "hm1",
             "secondaries",
             "heatmap",
-            "Confusion (single rollout)",
+            "Heatmap (single rollout)",
             "predicted",
             {
                 "series": {"flow": [[1, 0], [0, 1]]},
+                "reference": [[2, 0], [0, 1]],
                 "row_labels": ["0", "1+"],
                 "col_labels": ["0", "1+"],
                 "cbar_label": "count",
+                "log_color": True,
             },
         ),
     ]