Interface edep-sim h5 output by drinkingkazu · Pull Request #7 · DeepLearnPhysics/pysupera

drinkingkazu · 2026-03-10T16:08:12Z

This addresses #1

…. Voxelization of edep-sim output is added in the form of pre-processor.

…o handle very busy event topologies (which we plan to use for FM data generation). Algorithms (LEScatter merging ones) optimized for speed for similarly complex dataset (180s=>2s process time per event)

Copilot

Pull request overview

This PR implements interface support for reading EDepSim HDF5 simulation output files directly into pysupera, as requested in Issue #1. It adds an EDepSimHDF5Reader class and supporting infrastructure, renames ancestor_id to root_id throughout the codebase, introduces a VoxelizeProcessor preprocessing stage, and adds a Pipeline class for end-to-end event processing.

Changes:

New pysupera/readers/ package with EventReaderBase, EDepSimHDF5Reader, and Geant4 process-type classification logic for EDepSim HDF5 files
Rename of ancestor_id → root_id everywhere (data model, IO, partitioner, conditions, tests)
New VoxelizeProcessor preprocessor, Pipeline class in config.py, and performance improvements (bbox cache, parallel CC, KDTree-based LE scatter candidate building)

Reviewed changes

Copilot reviewed 27 out of 27 changed files in this pull request and generated 4 comments.

Show a summary per file

File	Description
`pysupera/readers/format_edepsim_h5.py`	New EDepSim HDF5 reader with G4 process classification
`pysupera/readers/example_edepsim_h5.py`	Draft/prototype script with critical executable module-level code
`pysupera/readers/base.py`	New abstract base class for all file-format readers
`pysupera/readers/__init__.py`	Package init exposing `EDepSimHDF5Reader`
`pysupera/data.py`	Renames `ancestor_id`→`root_id`; adds optional attributes and factory methods
`pysupera/utils.py`	Removes `+1` from `InteractionType` conversion, breaking native file reads
`pysupera/io.py`	Renames `ancestor_id`→`root_id` in HDF5 schema
`pysupera/preproc.py`	Adds `VoxelizeProcessor`; parallel CC; `last_stats`
`pysupera/config.py`	Adds `Pipeline`, `build_pipeline`, `build_voxelizer`, `build_reader`
`pysupera/app.py`	Adds EDepSim format selector, draw-mode toggle, cache; duplicates `_SEM_COLORS`
`pysupera/partitioner.py`	Bbox cache, per-condition timing, `batch_check_multi_cloud_proximity`
`pysupera/partitioner2.py`	Renames `ancestor_id`→`root_id` references
`pysupera/proxck.py`	Adds `batch_check_multi_cloud_proximity` to multiple backends
`pysupera/conditions/touching_em_shower.py`	Renames `ancestor_id`→`root_id`
`pysupera/conditions/combine_le_scatters.py`	KDTree-based candidate building
`pysupera/conditions/absorb_le_scatter.py`	KDTree-based candidate building
`pysupera/diagnostics.py`	Renames `ancestor_id`→`root_id` in strings/comments
`pysupera/_run.py`	Adds voxelizer support, switches to `partition_combined`
`pysupera/conf/config.yaml`	Adds voxelize config block and EDepSim reader default
`pysupera/conf/reader/edepsim_h5.yaml`	New reader config group for EDepSim HDF5
`tests/conftest.py`	Renames `ancestor_id`→`root_id` in `make_particle` helper
`tests/test_utils.py`, `test_preproc.py`, `tests/test_partitioner.py`, `tests/test_conditions.py`	Update tests to use `root_id`; add new tests for `VoxelizeProcessor` and `last_stats`
`README.md`	Documents new optional attributes, `root_id`, EDepSim format, draw modes

Comments suppressed due to low confidence (1)

pysupera/io.py:143

The comment at line 141-142 says "process_type is stored as raw int (0-based, as passed to SetSemanticType before the +1 inside that function)" but the +1 has been removed from SetSemanticType. This stale comment is now incorrect and misleading: it should be updated to reflect the new convention (whether the stored value is 0-based or 1-based). Moreover, if native HDF5 files still store 0-based values and SetSemanticType now expects 1-based, reading those files will fail.

            # process_type is stored as raw int (0-based, as passed to
            # SetSemanticType before the +1 inside that function).
            p_proc_type[j]   = int(p._process_type)

💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.

Copilot · 2026-03-10T16:15:14Z

+import time
+t0=time.time()
+parts = read()
+print(time.time()-t0,Particle.time)


The example_edepsim_h5.py file contains top-level executable code at module scope (lines 248-251) that runs read() on a hardcoded filename out_0100.h5 and prints timing. This will execute immediately upon import, which would fail with a FileNotFoundError (or NameError since h5 is not imported — only h5py is), crashing any process that imports from the pysupera.readers package. This file appears to be a development notebook/scratch script that was accidentally included in the package.

Additionally, the read function on line 215 uses h5.File(...) but h5 is never imported — only h5py is available. The np module used in search_parents (line 39) is also not imported at the top of the file (it is imported much later on line 207, after the functions are defined, but search_parents is called by get_parent_pdg and get_interaction_type which are defined before the import).

This file should either be removed, moved to an examples/notebooks directory outside the package, or at minimum guarded with if __name__ == "__main__": for the executable lines.

Copilot · 2026-03-10T16:15:14Z

    """

-    process_type = InteractionType(process_type+1)
+    process_type = InteractionType(process_type)


The change from InteractionType(process_type + 1) to InteractionType(process_type) in SetSemanticType breaks reading from native pysupera HDF5 files.

The native HDF5 format stores 0-based process_type integers (the comment at io.py line 141 confirms this: "process_type is stored as raw int (0-based, as passed to SetSemanticType before the +1 inside that function)"). The test conftest also confirms this: PT_TRACK = 0, PT_NEUTRON = 1, etc. (0-based mapping).

Since InteractionType uses auto() (1-based, starting at 1), calling InteractionType(0) — which occurs when reading any native particle with process_type=0 (kTrack) — will raise a ValueError.

If the intent is that _get_interaction_type in the EDepSim reader now returns 1-based .value integers (matching InteractionType's auto-assigned values), then the SetSemanticType function should keep process_type + 1 for the native format path, or the native format should be updated to store 1-based values, or the entire codebase needs to agree on a consistent convention.

The docstrings in both InteractionType (lines 9 and 12) and SetSemanticType (line 280) still say InteractionType(process_type + 1) but the code now does InteractionType(process_type) — this inconsistency confirms the change may be unintentional or incomplete.

Suggested change

process_type = InteractionType(process_type)

process_type = InteractionType(process_type + 1)

Copilot · 2026-03-10T16:15:14Z

        ``InteractionType`` is encountered.
    """

-    process_type = InteractionType(process_type+1)
+    process_type = InteractionType(process_type)
    if process_type == InteractionType.kInvalidProcess:
        return SemanticType.kUnknown
        raise Exception("'kInvalidProcess' particle process encountered\n")


The docstring for SetSemanticType (at line 280 outside the diff, and the InteractionType class docstring) still states that the process_type parameter is "0-based" and is converted via InteractionType(process_type + 1), but the implementation on line 309 has been changed to InteractionType(process_type). These docstrings need to be updated to reflect the new convention to avoid confusion.

Copilot · 2026-03-10T16:15:15Z

+
+_SEM_COLORS: dict[str, str] = {
+    "kShower":    "#e74c3c",
+    "kTrack":     "#3498db",
+    "kDelta":     "#2ecc71",
+    "kMichel":    "#f39c12",
+    "kLEScatter": "#9b59b6",
+    "kUnknown":   "#95a5a6",
+}
+


_SEM_COLORS is defined twice in app.py — first at lines 82-89 and then identically at lines 126-133. The second definition silently overwrites the first. One of these duplicate definitions should be removed.

Suggested change

_SEM_COLORS: dict[str, str] = {

"kShower": "#e74c3c",

"kTrack": "#3498db",

"kDelta": "#2ecc71",

"kMichel": "#f39c12",

"kLEScatter": "#9b59b6",

"kUnknown": "#95a5a6",

}

# (_SEM_COLORS is defined earlier in this module.)

… the unit tests, a remanent of C++ version of Supera. Fixed. Also process_type use was mixed with bare integer values, which is updated.

The visibility-filtering paths degraded silently when handed the wrong file in either slot: every particle came back with an empty point cloud, so a misconfigured run looked like a physics result rather than an error. The worst offender was n_actual defaulting to 0 via attrs.get(), which turned "this is not a seg volume" into "this volume has no segments". Six such branches now raise KeyError naming the offending path, what was found there, and which file belongs in that slot. Genuinely empty data (n_actual == 0, or readout planes present but nothing above threshold) is still allowed through -- only schema and pairing problems raise. Also accept deposit_to_group as an alias for segment_to_group, since newer JAXTPC output renamed that table, and import hdf5plugin so the Blosc-compressed seg/inst files can be decoded at all. hdf5plugin becomes a declared dependency; it was already required in practice for io.compression=lz4.

_CHUNK_POINTS was a fixed row count (262144) picked to make a 5 MiB chunk for the 5-column point array -- 5 MiB / (5 cols * 4 B) = 262144 exactly. The event clouds reused that same constant with 9 columns, so their chunks came out at 262144 * 9 * 4 = 9.4 MB: nearly double the intended budget, and large enough that reading one event's rows decompressed the entire dataset. Chunk shapes are now derived from a byte budget by _chunk_rows(), so the row count follows from the row width instead of being asserted independently: _CHUNK_PARTICLES 256 -> 16384 (64 KiB budget, 1-D metadata columns) _CHUNK_POINTS 262144 -> 13107 (256 KiB, 5-column point arrays) _CHUNK_CLOUD new -> 7281 (256 KiB, 9-column event clouds) _CHUNK_PARTICLES at 256 rows was 1 KB per chunk, small enough to hurt both compression ratio and per-chunk bookkeeping. The rdcc_nslots comment describing 5 MB chunks was stale and is updated; the formula itself already tracked _CHUNK_POINTS. Measured on a 5-event JAXTPC run: default LZ4 output 9.40 -> 8.05 MiB with random row-range reads 2.7x faster, and the gzip copy the WebGL viewer reads 6.79 -> 3.96 MiB and 5.2x faster. Also adds repack(), which rewrites a finished file with every chunk clamped to its dataset's real length. The streaming writer creates each dataset empty and grows it with resize(), so the final length is unknown when the chunk shape is fixed; the byte budget bounds the waste to at most one chunk per dataset, and repack removes even that. It walks groups and datasets generically rather than reconstructing events, so it preserves whatever layout the writer produced; copies in chunk-aligned blocks to keep memory bounded; and swaps the result in via os.replace so the original survives any failure. Third-party filters are reproduced by filter ID, because h5py reports LZ4 as compression == 'unknown', which is not a value it accepts on write.

Two config-level additions to the event loop. allow_empty_image (default false) halts on the first event whose particles all have empty point clouds. Such an event carries no reconstructable charge, and in practice that means a misconfigured reader rather than real physics -- the case that motivated it was a JAXTPC run that processed 5 events of nothing and reported success. Setting allow_empty_image=true downgrades the same condition to a single warning for the whole run and continues, for setups where empty events are genuinely expected (very low readout thresholds, heavily masked volumes). repack.enabled (default false) rewrites the output after the event loop finishes and every writer has closed, so chunk shapes can be clamped to the final dataset lengths. repack.compression optionally switches the filter in the same pass, so run_pysupera ... repack.enabled=true repack.compression=gzip yields a browser/h5wasm-readable file directly instead of requiring a separate pysupera-recompress step. The voxmap companion file is repacked too when one was written. The benefit is concentrated in short runs, where several datasets end up shorter than a single chunk: on a 5-event file 18 of 47 datasets are re-chunked (8.05 -> 6.86 MiB), while on long production runs almost none qualify. It costs one full read and write of the output, which is why it is opt-in.

In JAXTPC mode the inst file acts as a visibility mask over the energy depositions, but nothing reported how much it removed -- a run could discard most of the event's charge without saying so. The reader now records per-event counts on last_mask_stats: deposits present in the seg file, how many cleared the readout threshold, how many were masked out, and how many were attached to an EDepSim particle. A separate "visible but unmatched" count covers segments the readout saw whose group_to_track ID is absent from the particle list (including the -1 sentinel for out-of-range groups); those reach no point cloud, so silently losing them would be easy. _run.py aggregates the counts into a segment-visibility table in the end-of-run report, with the masked fraction over the whole run, and prints them per event when verbose=true. Non-JAXTPC readers supply no stats and the section is omitted. None of this reads anything extra. n_total comes from the per-volume n_actual attributes already loaded to build the point cloud, n_visible sums the index arrays the mask already produced, and n_attached is the length of the finished cloud. Deriving the total from the EDepSim step array instead would have been the expensive route: pstep/lar_vol is an object dtype dataset, so indexing one event materialises the whole ~11 MB record array for that event. Measured on a 5-event run: 0.129 ms/event, against ~389 ms for the event read and ~2795 ms for full per-event processing -- 0.005% of runtime. Collection is therefore unconditional and only the printing is gated. Sanity check: the reported total for event 0 (171,094) matches the EDepSim pstep count exactly.

Three related defects, all rooted in trusting fields that turned out not to be trustworthy in real EDepSim output. 1. ancestor_track_id is not the primary ancestor. Some secondaries are written with ancestor_track_id == their own track_id despite having a real parent, and every descendant inherits that false root. In one 13,460-particle event 15 such entries put 4,423 particles (33%) on the wrong root. _get_root_id now derives the primary by walking parent_track_id, memoising along each chain to stay O(n), and is used for root_ids in both readers. 2. Interaction IDs were stranded by the same field. _get_interaction_id seeded from particles it judged to be roots and propagated down the ancestor field, so a false root failed the vertex match and left itself and its subtree at -1: 118 of 13,460 particles in that event, every one of them recoverable from its parent chain. A third pass now walks parent_track_id for anything still unassigned, the root seeds come from _get_root_id rather than the stored field, and a root reference missing from the event list no longer raises IndexError. 3. Geant4 track IDs are not a safe index. They are not guaranteed contiguous, since an upstream stage may drop particles before pysupera sees them, so they cannot serve as direct array indices. _to_index_space assigns pysupera's own per-event particle index -- id == arange(n) -- and remaps parent_id and root_id into that space; a particle whose parent is absent becomes its own parent, the convention every genealogy walk already terminates on. The original track IDs are preserved as Particle.geant4_id, plumbed through from_arrays and from_flat_arrays. Verified on a JAXTPC run: the chain pi+ 28 -> mu+ 3065 -> e+ 3066 now carries root_id 20 and interaction_id 7 throughout, where the positron previously had no interaction ID at all; zero unassigned interaction IDs across three events; and every id/parent_id/root_id is a valid index into its event's particle table.

Format 2.2.0 adds a geant4_id column beside id at all three levels. id is now pysupera's per-event particle index, so it can be used as a direct row index; geant4_id carries the Geant4 track ID for provenance and for joining back to the EDepSim input. Adding the name to _SCALAR_FIELDS propagates it to particles, particle_fragments and particle_instances at once. Read paths fall back to id when the column is absent, so 2.1.0 files still load: in data where track IDs happened to be contiguous the two values coincided, and EventStore returns identical particles from either version. repack() also becomes a console script, since fixing up an existing file otherwise meant re-running the event loop or reaching for pysupera-recompress, which cannot resize chunks. Three arguments were added to cover standalone use: dst write to a new file rather than rewriting in place rechunk let h5py pick chunk shapes instead of clamping the source's verify compare every dataset against the source afterwards; for an in-place repack this runs BEFORE the original is replaced Compression is deliberately three-way, because two of the cases are easy to conflate: omitting --compression preserves each dataset's existing filter, --compression gzip switches it, and --compression none stores the data uncompressed. An earlier draft mapped "none" onto the preserve path, silently doing the opposite of what it said. Measured on a 5-event file (8.05 MiB, lz4): pysupera-repack out.h5 -> 6.86 MiB (lz4 kept) pysupera-repack out.h5 vis.h5 -c gzip -> 3.96 MiB pysupera-repack out.h5 vis.h5 -c gzip --rechunk -> 3.19 MiB pysupera-repack out.h5 raw.h5 -c none -> 17.78 MiB The gzip form replaces the pysupera-recompress step for producing a browser-readable file, and unlike recompress it repacks the voxmap companion too, since it walks the file generically instead of replaying events. recompress is left in place as public API.

Four event-loop changes, all about what ends up in the output and whether it can be trusted. Validation (utils.py, on by default) validate_interaction_ids every particle must carry an interaction ID. The diagnostic prints each offender's parent chain with PDG and interaction ID at every link, which is normally enough to see where propagation stopped. validate_root_ids root_id must equal the primary the parent chain leads to, and must name a particle present in the event. Both run after resolve_orphans, since that legitimately re-roots particles whose references dangle, and both are gated by check_interaction_id / check_root_id -- true raises, false warns once per run. Producer links no longer dead-end _find_parent_idx returned the first ancestor found in the member lookup. When a representative's parent had been merged into the representative's own group -- routine for an EM shower whose initiator is a late photon -- that lookup answers with the group itself, so parent_inst_id pointed at its own row and the production history stopped there. It now skips its own group and keeps walking. On a sample event this repaired 11 of 535 instances (and the same 11 fragments); self-parent rows fall from 48 to 37, exactly the event's primary count, with zero non-primary dead ends and every chain terminating on a genuine primary. instance_output (default 'traceable') Writes instances with a non-empty point cloud, plus the zero-point ancestors needed to complete their genealogy -- a pi0 is the only link between its two visible photons -- plus the primaries of every interaction that has at least one visible instance. An interaction with no visible instance is dropped whole: it left no trace, so nothing represents it. Retains 207 of 535 instances on the sample event, and loses no point data, because an instance cloud is the exact union of its members' clouds and so a zero-point instance has only zero-point members. 'all' restores the previous behaviour. write_fragments (default false) The fragment level is an intermediate product whose point payload duplicates what the instance level already stores; skipping it saves ~12% of the file. An empty list is still appended per event so the fragment event offsets stay consistent with zero rows, the same way particles/ behaves under output_compact. Instance parent_frag_id is set to -1 rather than pointing into a table that was not written.

… edges Three changes to vis_force.html for format 2.2.0 and the new output defaults. The genealogy graph was empty for every default-output file. Edges were built by walking the level-0 particles parent map, and parent_inst_id was never read at all -- but particles/ is empty under output_compact, which is the default, so findParentInst had nothing to walk. The stored parent_inst_id is now preferred and the walk kept only as a fallback for older files. On a sample event that yields 176 edges across 207 instances where there were none. geant4_id is read separately from id, falling back to id on pre-2.2.0 files where the two coincided. Both the hover tooltip and the locked bar show the two distinctly, so a particle index is no longer mistaken for a Geant4 track ID. Fragment fields are gated on the level actually being present. With write_fragments=false the table exists but is empty and every parent_frag_id is -1, which would otherwise render as an em-dash and read as "this particle has no fragment" rather than "this file has no fragment level". probeFragLevel() checks the row count at open time, and the tooltip omits the field entirely when there is nothing to show. The open-time log now reports format version and which optional levels the file carries. Note for future edits: the worker source lives inside a template literal, so backticks in comments there terminate it. Syntax-check with node --check <extracted script block> which is how that mistake was caught in this change.

_STEP_FIELD_MAP declared the dx candidates as ("dx") -- a string, not a 1-tuple -- so the lookup loop iterated the characters 'd' and 'x', found no field 'd', then matched the very real field 'x'. Every dx value in the EDepSim reader was the point's x coordinate: col dx (5) first 5 : [-1467.1515 -1467.2798 ...] step["dx"] first 5 : [ 0.2999999 0.2999999 ...] The JAXTPC reader was unaffected; it assigns the column directly. dX is now also written. _PC_NDIM goes 5 -> 6, so the per-particle clouds are (x, y, z, time, dE, dX) and format becomes 2.3.0. Both merge stages operate within a single particle, so dE and dX both sum over a particle-voxel and dE/dX is column 4 / column 5 exactly. Sanity check on real data: median 0.195 MeV/mm from both readers independently, right for a MIP in liquid argon. Guard against dX == 0 -- EDepSim writes zero-length steps, about 7% of voxels on the non-JAXTPC path. This also settles the two tests that have been failing since before this work. They asserted that column 5 takes the maximum, under the name test_dedx_takes_maximum. Both premises were wrong: the column holds dX, not dE/dX, and dE/dX is not an aggregation rule at all -- it is derived. Renamed to test_dx_sums_within_a_particle and corrected to expect 1.5 + 3.2 = 4.7. The suite is green for the first time in this branch.

Three defects in defragmentation, all found by the validators added earlier in this branch once the non-JAXTPC path was finally run with them enabled. Dropped parents orphaned their children. When every fragment of a particle is small the particle itself is removed, but its children -- both the pieces split off from it and its real children elsewhere in the event -- still pointed at it. resolve_orphans then rewrote each dangling parent_id to the particle's own id while leaving root_id on the original ancestor, producing particles that are their own parent yet claim a foreign root: 292 of 13,589 in one event. DefragmentBase.process now redirects references to a dropped particle's own parent, following chains of consecutive drops. Split pieces invented a Geant4 track ID. _split_fragments never passed geant4_id, so the constructor defaulted it to the fresh pysupera id and 154 particles carried track IDs that exist in no input file. A split piece is part of the same Geant4 particle, so it now inherits the provenance: `id` distinguishes the pieces, `geant4_id` records their shared origin. The voxel mapping described the wrong particle list. It was accumulated in voxelizer.last_diagnostics before defragmentation ran, so it referenced 25 particles that no longer existed and omitted the 154 that were created. The mapping is now attached to the Particle, so it follows the particle: a dropped one takes its mapping with it, and a split piece gets its share via subset_voxmap(). That works because voxel_offsets is a CSR over the particle's voxels in point-cloud row order, so the same boolean mask that splits the cloud selects the matching voxels -- no re-derivation needed.

…ress Format 2.3.0. _PC_NDIM 5 -> 6 so the flat point arrays are (x, y, z, time, dE, dX); see the reader commit for why dX is worth keeping. Older 5-column files still load -- Particle accepts any width >= 3, and _normalize_cloud zero-pads on write -- verified by reading a 2.1.0 file. flat/input_ids in the voxmap narrows from int64 to int32. Those index points within a single event (max ~3e5 observed), unlike the two offset arrays, which are cumulative over the whole file and stay int64. Note this is not a disk win: LZ4 was already compressing away the zero high bytes, so stored size is unchanged at 670 KiB. It halves the in-memory array and the decompressed bandwidth, and matters when compression is off. pysupera-recompress is removed, along with recompress() and the _read_rep_level() helper that only it used. pysupera-repack supersedes it completely: recompress reconstructs the file event by event, so it cannot resize chunks and silently ignores the voxmap companion, while repack walks the file generically. Keeping two commands where the older one quietly does less is worse than having one. Neither symbol was exported from __init__.

id is documented as a per-event particle index usable as a direct row index, and that was true of what the reader produced but not of what got written: defragmentation spawns particles with fresh ids and drops others, leaving 25 gaps in one event and a maximum id past the row count. The property held on the JAXTPC path only because nothing happened to split in the event that was tested -- defragmentation runs there too. Particles are therefore renumbered to 0..n-1 in list order once preprocessing has finished, with parent_id and root_id remapped. This runs after resolve_orphans, so every reference already points at a particle that is present and the remap is total. id == arange(n) now holds on the splitting path as well. The voxel map moves with it. It was written immediately after voxelization, which is before defragmentation, so it described a particle list that no longer matched the output. It is now written after renumbering and rebuilt from the final particles, each carrying its own share of the mapping. The result covers exactly the stored particles -- 13,589 of 13,589, no orphans in either direction, where there were previously 154 missing and 25 stale -- and total input_energies equals total voxel dE, which is the check that the per-split subsetting partitions rather than duplicates. particle.voxelize.store_mapping now defaults to false. The map carries one row per *input* deposit rather than per voxel, making it the largest artifact pysupera produces: 9.7 MiB against 6.0 MiB of actual output on a 5-event run. Nothing is computed when it is off, not merely not written, and the file is a cache of a deterministic computation -- re-running reproduces it exactly -- so the default costs nothing and saves the common case from paying for a debugging aid it does not use.

…r output The output format had no documentation outside io.py docstrings, which made the two integer spaces in the file easy to conflate and the per-event id base (events/offsets[ev] + id) invisible. README gains an Output File Format section covering the top level, the three particle levels with every field and which space it lives in, membership arrays, point-cloud layouts, and CSR event indexing. ID conventions get four rules: id is pysupera's, not Geant4's; ids restart per event; a primary is its own parent; and geant4_id is one-to-many, because a particle split by defragmentation yields several rows sharing one track ID -- 129 of 13,589 on a sample event. That last one is the sharpest edge in the scheme and was recorded nowhere. Expert output is split into its own section, stating up front that all of it is off by default, exists to debug pysupera rather than to be analysed, and is regenerable so not worth archiving: the fragment level, the voxel map (with a worked read and its size), instance_output=all, and enable_diagnostics. The fragment entry also spells out the naming trap that caused real confusion -- a particle_fragment is a group of *merged* particles, while *defragmentation* splits one particle into several, moving the count the other way. Also adds tools/h5_read_bench.py, which compares HDF5 read speed across files with matching layout. It interleaves repeats so no file gets a warmer page cache, reports the first pass separately as cold, and reads identical ranges from every file via a fixed seed. It is what established that chunk shape, not compression filter, dominated read cost here.

The DBG button's tooltip reads "Log debug info [D]", but the global keydown handler never mapped 'd' -- it covers c, a, s, space, r, o, i, g, ? and Escape. The key did nothing, so the only way to dump the camera position and loaded point counts was to find the button.

The IO layer had no tests. It is the part of the codebase most likely to be changed by a format revision and the only part whose breakage is silent -- a wrong column reads back as plausible numbers -- so it was the wrong place to have zero coverage, particularly with a 3.0.0 schema change coming. 29 round-trip tests pinning format 2.3.0: scalar fields, stored sem_type (which must be trusted rather than re-derived on read), point-cloud values and width including narrow-pad and wide-truncate, zero-point particles, empty events, multi-event boundaries with per-event id reset, bulk and streaming writers agreeing, representative levels with their member lists and group parent pointers, an empty level keeping its event fencepost consistent, event clouds, fencepost monotonicity and sizing, all four compression filters producing identical data, repack in-place and to-dst with filter preservation, compression switching, rechunking, the dst==src guard, verification actually detecting corruption, and the pre-2.2.0 fallback when geant4_id is absent. Reading is checked both through EventStore, which is what consumers use, and through raw h5py, which pins the dataset layout other tools depend on -- the WebGL viewer reads these names directly. The first run failed 17 of 29, all on the same defect: write_events() never wrote particles/interaction_id, in any commit in this repository's history. read_events() requires that column, so the bulk writer has always produced files its own reader cannot open -- despite being exported and shown in the README quick start. It survived because everything real goes through the streaming EventWriter. geant4_id was missing from the same path. Both are added here, which is why io.py is part of this commit.

"root" said nothing about which root. The field is the primary the particle descends from, so name it that. This commit carries the files that contain only the rename; _run.py, io.py, preproc.py, utils.py, the readers and the viewer carry it too, alongside their own changes, and say so. No behaviour change.

Three unrelated defects in one layer. dX was never read. _STEP_FIELD_MAP had (PointFeature.dx, ("dx")) where ("dx") is the string "dx", not a tuple, so the loop iterated its characters and matched the field "x". Every dX in every output file since the column was added held an x coordinate. Now ("dx",). ancestor_track_id from EDepSim is wrong for decay products, so it is no longer read. _get_ancestor_id walks parent_track_id instead -- memoised and cycle-safe -- which is the field that is correct. On the old sample, 15 seeds stranded 4,423 of 13,460 particles in one event; across five events, 34,570 of 123,231 (28.1%). Reported upstream; see edepsim-ancestor-bug.md. The JAXTPC path had six branches that silently produced an empty point cloud when an expected dataset was missing -- a wrong --jaxtpc-seg-path gave zero points and no complaint. Each now raises KeyError naming what it looked for and what the file holds. It also accepts segment_to_group or deposit_to_group, imports hdf5plugin so Blosc-compressed inputs open at all, and records per-event mask statistics. get_visible_segment_indices_by_track grouped segments with a scan per track, O(tracks x segments): ~5k tracks each compared against 250k segments, 0.41 s of self time and the single largest cost in the pipeline -- in a function that is pure bookkeeping. Sorting once and splitting at the boundaries is O(N log N). End to end on a 38k-particle event: 3.45 s -> 3.04 s, 12%. Both files also carry their share of the root_id -> ancestor_id rename.

Defragmentation dropped a particle when all of its fragments fell under min_pc_size. That is the right call for a leaf and badly wrong for anything with children: the node vanishes and its whole subtree is stranded. The case that exposed it: a 534 MeV primary photon whose only deposits were two points 113 mm apart -- 61 eV of Compton and 3178 eV of argon K-shell relaxation at photoabsorption. Both clusters were size 1, so the photon was deleted and the 384-particle shower hanging off it re-rooted on a secondary electron. Over 100 events this removed 2,314 photons, every one of them the parent of an e+/e- pair, 77% of them heading a subtree larger than min_pc_size. Only one was a primary, so a primaries-only guard would have missed almost all of it. A fully split particle is now kept as a genealogy node with an empty cloud; its points live on in the spawned pieces. A neutral particle with no cloud of its own is the ordinary case the merging conditions already expect -- PhotonDecay is genealogical and needs no points -- so the node costs nothing. Nothing is dropped any more, which leaves the dropped-particle remap below it as defence rather than routine. Two provenance fixes in the same function. A spawned piece inherits the original's geant4_id rather than minting one from the fresh pysupera id, which previously fabricated track IDs present in no input file. And its parent_pdg is now the original's own PDG: parent_id is set to p.id, so copying p.parent_pdg left the two fields describing different particles. The dropped-particle remap likewise updates parent_pdg when it redirects parent_id, which it never did. After both, parent_pdg matches the parent's actual PDG in 0 of 433,593 particles over 20 events; it was wrong for every spawned fragment before. Also: dX is summed on merge (merging only ever happens within one particle, so path length adds), and subset_voxmap keeps the voxel mapping aligned when a cloud is split. Carries its share of the root_id -> ancestor_id rename.

resolve_orphans repaired parent_id and ancestor_id independently: if p.ancestor_id not in known_ids: p.ancestor_id = p.id Self-rooting is only correct when nothing survives. When a particle is removed but its descendants' parent links still terminate at a live particle, this scattered one shower across hundreds of false roots -- and contradicted validate_ancestor_ids, which walks the parent chain. In one event 373 particles were left claiming to be their own ancestor while their parents led to a live particle two links away. It now runs in two passes: repair every parent_id first, then set a dangling ancestor_id to the root of the repaired chain. _chain_root memoises along the way, since one shower shares one chain, and terminates on cycles. Over 100 events: 0 inconsistent ancestors, 0 missing interaction IDs. Three checks added, because both classes of defect were silent: validate_interaction_ids every particle must carry an interaction validate_ancestor_ids ancestor_id must equal what the parent chain says select_traceable_instances the instances worth writing, plus the zero-point ancestors needed to reach them PointFeature.id is documented rather than changed: it holds the input step's row index, in-memory clouds are (N, 7) while the file stores only columns 0-5, and the merge rules define no aggregation for it. Carries its share of the root_id -> ancestor_id rename.

geant4_id records where a particle came from, separately from the pysupera id, which is a dense per-event row index. It defaults to id so existing callers are unaffected, and it is what makes a join back to the input file possible -- the pysupera ids are assigned by the labeller and mean nothing upstream. voxmap carries the optional voxelization mapping so it can survive the splits and merges that follow. Carries its share of the root_id -> ancestor_id rename.

Chunk shapes were fixed row counts, so a chunk's size in bytes tracked the column count -- 262144 rows of a 9-wide float32 cloud is 9 MiB, far past anything useful, and HDF5 must read a whole chunk to return one row. Shapes now come from a byte budget: 256 KiB for bulk arrays, 64 KiB for metadata. recompress() only changed the filter and left chunk shapes alone, which is the half that matters least. repack() replaces it: it clamps every chunk to the dataset's real length -- the streaming writer cannot know that length when it creates the dataset empty -- and can change the filter in the same pass. It walks groups generically rather than reconstructing events, copies in chunk-aligned blocks so memory stays bounded, and verifies dataset equality before an in-place swap. Two options ported from the standalone h5_convert_compression.py, which this replaces: --compression auto flip the file's dominant filter between lz4 and gzip, so the round trip to a viewable copy needs no memory of which way round it currently is. Errors rather than guess on an uncompressed file. --scope all compress every chunked dataset, not only those already compressed (the default, "source"). On a 100-event file gzip is 65.1 MiB against LZ4's 88.8 MiB, and only gzip is readable by the browser viewer, which uses h5wasm. Also fixes write_events() never writing particles/interaction_id -- absent from every 2.x file ever produced -- and reads ancestor_id from a legacy root_id when present. Carries its share of the root_id -> ancestor_id rename.

New module, separating what to write from how to write it. io_v3 does no derivation; everything below happens here. Point ordering is interaction > instance > is_LE > fragment > particle which makes three things true at once: an instance's points are one contiguous block split into a non-LE and an LE run; a fragment's points are exactly two runs, one per side, never more; and LE-ness is positional, so no per-point label is stored. The query training issues most often -- the non-LE voxels of instance X -- is then a slice rather than a mask. snapshot_groups has to be called at the right moment and the reason is worth stating: merge_em_showers mutates the surviving representatives in place and returns those same objects as instances, so after step 2 a representative's member_ids and sem_type describe its instance and the fragment values are gone. Reading fragment membership afterwards silently misattributes most particles -- 10,299 in one event before this was understood. Fragments are snapshotted before the merge, instances after. parent_id is redirected to the nearest *stored* ancestor, so walking the genealogy in the output always terminates; ancestor_id is then derived from the redirected chain so the two agree. parent_pdg describes whatever parent_id points at. Because only ~17% of particles are stored, the true Geant4 parent is unrecoverable from the file alone, so it is kept explicitly as geant4_parent_trackid / geant4_parent_pdg. The interactions table keeps live interactions only, renumbered densely, with interaction_id as the index particles refer to and vertex_id as the surviving link back to the EDepSim vertex row.

2.x wrote three parallel particle tables -- raw, fragments, instances -- and a point cloud per level, so a point reachable at two levels was stored twice. 3.x writes one particle table and one point array, with every level expressed as index ranges into it. On a two-event file: 6404.6 KiB -> 1660.0 KiB, 26%. 3.1.0 settles the naming. The level qualifier is a prefix everywhere (frag_sem_type, frag_parent_id, inst_sem_type, inst_parent_id), point ranges carry pc at every level (frag_pc_start, inst_pc_start), and both levels describe their two runs the same way -- a non-LE range and an LE range. At instance level the runs are adjacent, so inst_pc_end == inst_pc_le_start; the pair is spelled out anyway rather than encoded as a single split index, so instance and fragment level read identically. geant4_id becomes geant4_trackid alongside the new geant4_parent_trackid. The interactions table gains its own dense interaction_id, and carries energy_sum, ke_sum and reaction through from EDepSim. 29 particle columns, 13 interaction columns, 6 point columns. 3.0.0 files still read: _read_legacy_particles maps the old names, synthesises the two instance ranges from the old split index, and fills the Geant4 parent columns with -1, which is honest -- the value cannot be recovered from the file. EventView rebases absolute point ranges to event-local on read, so every index a caller sees shares one origin. points_of(row, level, le) answers every combination in a single slice except fragment-level "all", where the two runs are genuinely not adjacent.

The event loop writes through layout + io_v3: fragments are snapshotted before merge_em_showers, ids renumbered densely after resolve_orphans, and the voxel map written against the renumbered ids. Three guards, all for failures that previously passed unnoticed: allow_empty_image an event where every particle has an empty cloud is almost always a misconfigured reader, not physics. Halts by default; set true for a one-time warning. check_interaction_id one broken ancestor link strands a particle and all its descendants at -1. check_ancestor_id ancestor_id must be what the parent chain says. An optional post-run repack block clamps chunks once the final lengths are known, and can emit gzip directly for the viewer. output_compact, write_fragments and write_event_clouds are retired: the single table makes the first two meaningless, and instance_output=traceable expresses the intent better -- keep every instance reachable from a visible one, which on a sample event is 207 of 535.

PhotonDecay merges a photon's PDG +/-11 children into its partition genealogically, with no proximity test, which is right: a photon's shower can start hundreds of mm from where the photon itself deposited. Defragmentation can also cut the photon's *own* cloud into pieces. Those pieces are PDG 22, so they fell outside the condition and could only rejoin by drifting within distance_threshold of a sibling electron. That is a coin flip. In the measured case the two pieces sat 0.01 mm and 5.30 mm from their nearest sibling against a 5.2 mm threshold: one rejoined, the other became a standalone one-voxel instance, losing by a tenth of a millimetre. A PDG 22 child is now claimed when it shares the photon's geant4_id, which means it *is* that photon rather than a separate particle. Matching on geant4_id rather than PDG alone keeps genuine photon-from-photon secondaries -- fluorescence X-rays, which are distinct tracks -- out of this path. In one event that is 94 split-off pieces against 6 fluorescence photons. With this and the preproc keep-alive, the 534 MeV photon of event 81 comes out as one kShower instance of 386 merged particles: 640 non-LE points from the electron shower it initiated, 357 LE points from its own stray deposits and its low-energy Compton products.

Supporting module for the training pipeline. memory_snapshot() reports host RSS -- including DataLoader workers, which is where reading and pre-processing actually allocate -- alongside the CUDA allocator's view and the driver's free/total. Enumerating child processes scans /proc and costs ~4.6 ms, against ~0.2 ms for the rest of the snapshot, so the handles are cached. An *empty* cache expires far sooner than a populated one: nothing in it could signal that a worker pool has since appeared, and without that a freshly spawned pool reads as zero worker memory for the whole refresh window. The psutil handle is also rebuilt when the pid changes, since Process() captures a pid at construction and a forked worker would otherwise report its parent's memory. WandbLogger collapses every reason not to log into the same silent no-op -- wandb absent, no active run, disabled, or a rank that should not log -- so instrumented code needs no conditionals. In a distributed job only rank 0 logs by default, since every rank sees a different shard and letting them all write to one run interleaves unrelated series.

A batch is a *list* of events rather than a padded tensor -- point clouds differ in length, sparse and transformer backbones both prefer the ragged form, and batch_size then means exactly what it says. Model input is one row per voxel, but a voxel can receive energy from several instances. Energy is summed so no charge is lost; labels come from the winning contribution, by semantic priority (track > shower > delta > michel > LE) and then earliest time within a class. On a sample event ~10% of raw points share a cell with another instance's. LE-ness is a property of a point, not of its instance: a point in the inst_pc_le range is labelled kLEScatter whatever its host is, which makes it the lowest-priority contributor to a contested voxel, while its instance label stays the host's so an absorbed deposit remains attributed to the shower that absorbed it. LE points are excluded by default -- around a third of the total and usually noise to a segmentation model -- and the flag governs input and labels together, so the model never sees a point the labels do not describe. StreamMonitor separates two kinds of timing that are easy to confuse. read_s/preprocess_s/collate_s are worker-side *costs*, running concurrently with the training step, so their sum routinely exceeds wall time. wait_s, stage_s and consume_s are main-process wall time and partition one iteration; stall_fraction = wait_s / step_s is the number that says whether input is the bottleneck. summary() measures wall_s with one clock spanning the loop and reports unaccounted_s against the per-step breakdown, which is how the 27% overhead of the uncached psutil walk was caught. It also warms the CUDA transfer path on construction: the first host-to-device copy in a process pays ~35 ms of lazy initialisation against ~1 ms warm, and allocating on the device is not enough to trigger it.

The viewer reconstructs the old per-level view from the single table: instance blocks give ownership, LE-ness is read positionally from the inst_pc_le range, and fragment labels take two runs per side. h5wasm reads the raw file, so it gets no help from the Python reader's compatibility layer and has to know the 3.0.0 column names itself. V3_ALIAS plus v3Name() resolve either spelling, including the old single split index, so files written before the 3.1.0 rename still open. Also wires the [D] shortcut the DBG button has always advertised.

hdf5plugin was a hard requirement and an undeclared one: it registers the Blosc and LZ4 HDF5 filters, so without it JAXTPC inputs cannot be opened at all and io.compression=lz4 cannot be written. torch and wandb are genuinely optional -- pysupera.torchdata is not imported from __init__, and WandbLogger is a no-op without wandb -- so they are extras. psutil rides with torch, since the worker-memory figures are the part of the snapshot that needs it. Drops the pysupera-recompress entry point, whose function repack replaced.

Documents what the file actually contains: the point ordering and the three consequences that follow from it, which particles are stored and why, the full column list at all three levels, and the ID conventions -- ids are dense per-event row indices, Geant4 track IDs are provenance, and you join back to the input on geant4_trackid, never on id. Adds sections for the training pipeline and its profiling, including the distinction between worker-side cost and main-process wall time that the metric names alone do not convey, and the measured numbers behind the advice (persistent_workers is worth more than the worker count; discard the first epoch when benchmarking). Also documents repack's auto and scope options, and the hdf5plugin, torch and wandb dependencies.

A correctly drawn muon track was reported as PDG 2112, a neutron. The geometry was right because it comes from point ranges, read per row; the PDG was wrong because it came from whatever instance happened to sit at that array slot. The assumption dates from 2.x, where instances had their own table and a row index equalled an instance id. In 3.x instances are a sparse subset of one shared table, ordered by the point layout rather than by id, so the two are unrelated: a 2,902-row event carries ids up to 21,254, unsorted. Five places used one where the other was meant -- * the per-point PDG buffer, keyed by node index and read by id * the hover panel's node lookup * the point-count and interaction patch onto graph nodes * the highlight uniform, which the shader compares against a point's id while the graph passed a node index * the graph's parent edges, in both copies of the builder, which used inst_parent_id -- an id -- directly as a node index -- and of 850 instances in one event only 89 ids even landed inside the node array. 85 of those returned another particle's PDG; the other 761 ran off the end and displayed nothing. _instNodeIdxById makes the conversion explicit at every boundary between the cloud (which speaks ids) and the graph (which speaks node indices). The README said id "is 0..n-1 over the event's stored particles, so it indexes rows directly", which is what the viewer believed. It is false -- ids are dense over the *full* particle list, assigned before the output subset is chosen -- and it is now corrected, with the map to build instead, and the warning extended to parent_id, frag_id, inst_id and the two group parents.

Answers "which readout hits belong to this reconstructed object" without going through voxels. A hit is the projection of a JAXTPC group, and group_to_track gives each group one Geant4 track -- but defragmentation splits a track across several pysupera particles, so track -> particle is one-to-many and cannot regroup hits. What rescues it is that a group never straddles a split: a group is itself a tight spatial cluster and defragmentation clusters at distance_threshold. Measured over eight events, none of 354,959 groups spanned two particles, and group_to_track agreed with the owning particle's geant4_trackid every time. So group -> particle is a function, and hit -> group -> particle -> fragment is a chain of functions. Two arrays per group are stored, and no more. groups/fragment_id The particle would be the natural key and does not work: a group's owner is usually an ordinary member, and only about a sixth of particles get a row -- 6,605 of 6,777 owners were unresolvable when this was first keyed on particle id. Every fragment with points is stored by construction (0 missing of 29,201), so the fragment is the finest key that always resolves. groups/is_le Irreducible. LE-ness is a pysupera-particle attribute and a fragment mixes LE and non-LE members by design -- that is what its two point ranges are for -- so 26-30% of fragments contain both. The bit costs 0.17 bytes/group; making it derivable instead, by storing the owning particles as rows, costs about 20x more. The instance is *not* stored: a fragment's points lie inside exactly one instance's block, so it follows by containment, verified exact on 29,201 groups. instance_of_groups() does it. The interaction follows from the instance. Note this cannot be a column join -- inst_id and frag_id are self-markers, equal to id on a representative and -1 elsewhere -- which is why a first attempt at a join disagreed on 12,798 of 29,201 groups. check_group_ownership (default true) verifies the one measured assumption per event rather than trusting it, naming the group and both particles and pointing at store_mapping=true -- the exact deposit-level record the table compresses -- as the alternative. groups/ costs about 21% of the output under LZ4. Going via voxels would be worse as well as larger: a group's deposits land in 2-4 different voxels 38% of the time, so hit <-> voxel is many-to-many and would need charge apportionment. Supporting changes: * the JAXTPC reader keeps the real deposit index instead of overwriting it with a running counter, so the mapping reaches the input at all. It rides as an int64 array beside the cloud, not in it: the cloud is float32 and exact only to 2**24, which a hundredfold rise in deposits per event would reach. * computing the CSR is separated from writing it (track_provenance), so the group table can be built in memory while the deposit-level map stays opt-in and off by default. * input_energies is dropped from the voxmap. Nothing consumed it but one test, and once input_ids are real deposit indices the energy is recoverable from the seg file -- storing it duplicated upstream data. The test now looks the energies up through input_ids, which is what a caller has to do, and a second test pins that the CSR partitions the inputs. The voxmap file shrinks 18%. * subset_voxmap gains tests. It is called on every defragmentation split and had none, which is how it shipped briefly returning a stale 3-tuple. groups_of_particles drops NO_OWNER from the request rather than matching it: a caller passing an unset id would otherwise get back every unclaimed group, confidently and wrongly.

vis_hits.html puts the 3D voxels next to the 2D readout planes and colours both by the same owner id, so a colour means the same reconstructed object in either view. Hover a 3D voxel and its hits light up across every plane image; hover a 2D pixel and its voxels light up in 3D. That is the whole demonstration: the link is a single integer, and the two views agree. A separate file rather than an addition to vis_force.html, which is already 122k characters and answers a different question. Colour by instance, fragment or interaction. Only the fragment is read from the file -- the instance is derived by containment and the interaction from the instance -- so the viewer demonstrates the same derivation the format relies on rather than depending on stored duplicates. Unselected points are genuinely transparent, not merely darkened: the cloud uses RGBA vertex colours with depthWrite off, so the selection is never hidden behind the points being de-emphasised. A slider sets how transparent. "hide LE" drops low-energy points from both views, filtering at build time so the picker cannot select something hidden. A per-point readout shows x/y/z, time, dE, dX, dE/dX, LE and owner in 3D, and wire, time, peak charge, group, LE and owner in 2D. It refuses to draw rather than show something false. h5wasm decodes gzip and nothing else, and an LZ4 dataset yields garbage instead of raising -- which cost an afternoon: a viewer showing a muon track labelled PDG 2112, from a file that looked repacked but had kept its original filter. The viewer now checks the filter on open, and also banners a missing groups/ table or an event where every voxel resolves to one owner. The 3D canvas is absolutely positioned, sized from getBoundingClientRect with a window fallback, tracked by a ResizeObserver, and drawn by a rAF loop. A zero-sized canvas renders nothing *and* swallows every mouse event, which is indistinguishable from a broken viewer. pysupera-hits-subset prepares the second input: the JAXTPC hits file is 1.3 GB, almost all of it per-tick arrays the viewer never opens. Keeping the four arrays it does read gives 1279.7 MiB -> 2.73 MiB for three events. Verified against a three-event sample: 22,067 voxels and 87,268 hits, every one owned, all 196 instances present in both views. Rendering itself is unverified -- there is no browser in the build environment. The README gains a Browser viewers section documenting vis_force.html and vis_hits.html side by side -- what each is for, what it needs, its controls and shortcuts. vis_force.html had never been documented at all, only pysupera-app had. Writing it up turned up two labels still saying "Root Track ID" and "Root tracks" from before the root_id -> ancestor_id rename; both now say ancestor.

…ource One change rather than seven because the pieces landed on the same files and cannot be separated after the fact. What follows is what is in it. Hit provenance ============== Answers "which readout hits belong to this reconstructed object" without going through voxels. A hit is the projection of a JAXTPC group, and group_to_track gives each group one Geant4 track -- but defragmentation splits a track across several pysupera particles, so track -> particle is one-to-many. What rescues it is that a group never straddles a split: measured over eight events, none of 354,959 groups spanned two particles. So hit -> group -> particle -> fragment is a chain of functions. Two arrays per group, and no more. groups/fragment_id, because the particle is the natural key and does not work -- only about a sixth of particles get a row, and 6,605 of 6,777 owners were unresolvable when this was first keyed on particle id, where every fragment with points is stored by construction. And groups/is_le, which is irreducible: a fragment mixes LE and non-LE members by design and 26-30% contain both. The instance is derived by containment and the interaction from the instance; neither is stored. That cannot be a column join -- inst_id and frag_id are self-markers, equal to id on a representative and -1 elsewhere -- which is why a first attempt disagreed on 12,798 of 29,201 groups. groups/ costs about 21% under LZ4. Supporting: the reader keeps the real deposit index instead of overwriting it with a running counter, riding beside the cloud as int64 rather than in it (float32 is exact only to 2**24). Computing the CSR is separated from writing it, so the group table can be built while the deposit-level map stays opt-in. input_energies is dropped from the voxmap -- nothing consumed it, and once input_ids are real deposit indices the energy is recoverable upstream; the file shrinks 18%. subset_voxmap gains the tests it never had, which is how it shipped briefly returning a stale 3-tuple. Wire and pixel readout ====================== JAXTPC simulates both LArTPC geometries: wire records three 2-D projections (U/V/Y planes), pixel records one natively 3-D image (a single Pixel plane). Only wire had ever been exercised. The reader needed nothing. Visibility is a question about a *group* -- a group is above threshold or it is not, whatever geometry saw it -- and the filter already discovered plane subgroups structurally and read only group_ids. Rather than add a branch that would have to be kept honest, this documents why there is none and tests it: the same groups fed through both geometries must give the same visible segments. What does differ is anything that draws a hit. pysupera-hits-subset hardcoded U/V/Y and center_wires, so handed a pixel file it wrote an empty one and said nothing; it now finds planes structurally and keeps whichever centres they store. vis_hits.html gains a second 3-D pane for pixel readout -- x = drift tick, y = py, z = pz -- in place of the 2-D canvases. Two bugs found there, both of which made the feature dead rather than merely wrong. The load gate tested groups/instance_id, a column groups/ has never had by design, so every output ever produced reported "no groups/ table" and drew nothing. And group numbers were shifted twice, gids + base - ga, where the owner array is already sliced to the event; event 0 survived because ga is 0 there, and event 1 resolved 1,686 of 12,450 hits instead of all of them. The detected image as the point source ====================================== reader.point_source=hits replaces each particle's point cloud with the pixel hits themselves. A point stops being a Geant4 deposit that survived the readout and becomes a fired pixel, so the cloud carries pixelation, diffusion, threshold and the drift-time ambiguity, and every algorithm downstream sees what the detector saw. The assignment stays exact: a CSR entry belongs to one group and group_to_track gives each group one track. y and z invert exactly. x is drift time, and a tick measures t_drift + t0 -- t0 being precisely what a real detector must determine separately. So hit_x_from chooses: 'nominal' assumes the beam time, leaving x offset from truth by drift_direction * v * t0; 'true_t0' puts the real t0 back. On this batch t0 is constant within an interaction (max spread 0.53 us) and differs across the interactions of one event by up to 1094 us -- 175 cm. So nominal is a rigid translation per interaction, and the two conventions partition almost identically (4,723 against 4,727 over two events), which is itself the answer to a question worth asking. points/true_x_shift keeps the truth anyway, as a displacement rather than a second coordinate: true_x = x + shift, zero under true_t0 and in deposit mode so the relation holds for every file. It costs 0.14% of the output, because it is drift_direction * t0 * v and takes one value per (interaction, volume) pair -- 28 in a test event against 4.18M hits. Per point rather than per particle because the two exceptions are cathode-crossing tracks carrying 15% of the points. Geometry is configured and then checked, never measured. An earlier round fitted the pitch, velocity, tick length and drift direction from the truth deposits and checked the result against the same deposits, which is circular in the direction that matters: a geometry wrong but linear cannot fail such a check. Now the anode comes from config/volume_ranges, the velocity and sampling period from the sensor file, and the pitch and drift direction from config. A wrong pitch of 4.0 instead of 4.32 shows up as a 248 mm residual; a mirrored drift direction as 277 mm. Before, all of them fitted cleanly. What the readout window hides ============================= num_time_steps is 2701 and a tick is 0.8 mm, so the recordable ticks span exactly one full drift and a nominal x cannot leave the volume. The t0 displacement instead pushes charge off the end of the window, and that turns out to be nearly the whole story of what "visible" means here: per event and volume, the fraction of masked deposits that the window never reached is 97.0%, 98.1%, 99.7%, 99.4%. A deposit is invisible because the readout had closed, not because it was below threshold. Reported per run. A group can also straddle a fragment split now -- 1.84%, against 0 of 354,959 for truth deposits -- because the readout leaves gaps the deposit cloud does not have. build_group_owners grows resolve= and report=, hit mode settles by majority, and check_group_ownership must be off. Classification by extent, not by row count ========================================== min_pc_size decides whether a particle is a low-energy scatter, and it counted rows. Row count is a property of the input sampling, not of the particle: deposits arrive every 0.3 mm and pixel hits on a 4.32 mm grid, so the same object is 1 row in one and 40 in the other. With point_source=hits the size rule stopped contributing almost entirely and 3,818 particles the truth input calls low-energy scatters came out as showers. It now counts occupied voxels, which is what defragmentation already splits on, so the two stages finally speak one unit. The value follows from a physical rule: a chunk is worth calling a trajectory only when its direction can be read off it, so the threshold is the cell count of a blob one track-width across. Truth deposits are 4.9 mm wide -> 5, which is the long-standing default and reassuring about where it came from. reader.hit_charge_threshold drops hits below a given induced charge in electrons, on the magnitude, matching JAXTPC's own encoder -- the response is bipolar and 16% of hits are negative signal. JAXTPC thresholds at 1 electron, which is none, so this is where a readout threshold enters. It interacts with the extent threshold, because cutting the halo thins the tracks: hit_charge_threshold track width min_pc_size agrees with truth@5 ---------------------------------------------------------------------- 0 e- 13.1 mm 85 99.0% 500 e- 6.5 mm 13 95.5% 500 e- 6.5 mm 85 (mixed) 85.2% Measured over ten events; the per-event best cut ranged 81-90 at threshold zero. Agreement means the same particles, not merely the same fraction. Mixing the pair is worse than either, and no single config group could express that they belong together -- min_pc_size is under particle, the threshold under reader, distance_threshold and check_group_ownership at the top level. Hence preset/, an optional group listed after _self_ so it can override the defaults while a command-line flag still overrides it. preset=cubic_pixel_hits carries the matched set and a sensor-hit run is now just file paths. Note this also moves truth-deposit labelling, 85.8% -> 95.3% LE, because a particle spread over a handful of rows inside one or two voxels is now measured as what it is. A correction, but a change. Configuration ============= conf/reader/ held one 197-line file covering three unrelated inputs, almost all of it nulls that mattered to one of them. It is now edepsim_h5 (plain EDepSim, and the base), jaxtpc_wire and jaxtpc_pixel, the latter two inheriting the first through a within-group defaults list. jaxtpc_wire carries no point_source or pixel-geometry keys at all: a wire plane has no 3-D image, and options that cannot work should not look as though they might. Selecting a JAXTPC config without its input paths used to fall through to plain EDepSim in silence -- OmegaConf's .get() returns None for an unfilled mandatory value and `key in cfg` is False for one, so both obvious guards fail open. You would ask for reader=jaxtpc_pixel and get truth steps. It now refuses, naming the flags. The output also records point_source, hit_energy, hit_x_from and readout_type as root attributes, because nothing in the point columns says whether the energy column holds true dE or measured charge, and inferring it from negative values works on a busy event and fails on a quiet one. Viewers ======= vis_drift.html puts true x beside inferred x for the same points, framed on the union so a displacement reads as one. vis_truth.html puts the Geant4 deposits beside the pysupera cloud, painting each deposit with the label of whatever object claimed it -- the label travels backwards down deposit -> group -> fragment, and a fragment's truth deposits sit a median 2.2 mm from that same fragment's hits. Deposits whose group was never detected are shown as unowned rather than dropped: 34% of event 0. Both link by hover, colour by energy / semantic type / fragment / instance / interaction / PDG, and carry a per-pane display threshold showing the absolute cut and what survives. On the sensor side half the points cost 0.4% of the charge; on the truth side dropping the faintest 60% of deposits costs 61% of the dE, which is the halo stated in one number. Colour was chosen by rule and the data pushed back twice. The energy column is charge in hit mode and 13% of voxels are negative, so signed data about a meaningful zero wants a diverging ramp, not a sequential one. And linear is useless across twelve decades, so the default is sign-preserving log with the range clipped to the 2nd-98th percentile. Six categorical hues cannot clear the all-pairs colour-vision floor -- no six in one lightness band can -- so identity is never colour alone: a legend is always up, the hover names the class, and the filter can isolate any pair. Two rendering problems only a browser could have shown. A third of a million points with depthWrite off is a pale wash, every overlap blending instead of occluding; and once the bulk writes depth, the selection hides behind whatever faint points sit in front of it, so it became a second pass with depthTest off. Also: vis_truth was composed from vis_drift and syncFrom came across still naming the other file's panes, so the camera-sync checkbox did nothing at all there. It now copies to every view rather than a named sibling. Panning existed on right-drag and nobody finds it, so the arrow keys pan the pane under the cursor. pysupera-truth-subset writes the four arrays vis_truth needs -- positions, de, charge, deposit_to_group -- into one gzip file, because both JAXTPC files are Blosc and h5wasm decodes gzip only; 65 MB becomes about 1 MB per event. Verification ============ 417 tests pass, 91 of them new. Everything visual was checked in a real Chromium against real output rather than by node --check: all six colour modes, both thresholds, the hover on both panes, the linked highlight, the semantic filter, camera sync in both directions and with the box off. The container needed DISPLAY unset before ANGLE would stop trying to reach an X server and fall back to swiftshader, which is why earlier rounds of these viewers went unrendered.

…s/flat For JAXTPC runs the model's input is the JAXTPC hits file itself, so its targets now have to be one label per hit, in the hits file's own order. Before this, pysupera wrote the hits it had partitioned into points/flat, together with truth/, points/hit_index and points/true_x_shift, so that a loader could map them back. Wire and pixel did this differently. Now both readouts write the same thing: hit_labels/volume{V}/plane{P}/ (attrs: source U|V|Y|Pixel, volume, plane) fragment_id, instance_id int32, -1 = untraced is_le bool on_sensor bool, pixel only offsets per-event CSR Entry k labels CSR entry k of that plane in the JAXTPC hits file, so no reordering or join is needed. For wire, a hit belongs to its group's owner, and that is exact because a group's deposits are never split. For pixel, the hits are what was partitioned, so every hit is read directly. -1 marks exactly the hits off the sensor image, and check_pixel_labels raises if that does not hold. Over 10 pixel events, 24.9M of 28.9M hits (86.1%) carry a label, and the rest are all off the sensor. points/flat now holds the true energy deposits behind the hits: every EDepSim segment whose JAXTPC group left hits, read from pstep (matched 1:1 with the JAXTPC segments, per volume and in order). By default they are voxelized on 3 mm (io.store=voxels; io.store=points keeps the raw segments). There are nine columns, named in points/columns: x, y, z [mm], t [us], dE [MeV], dX [cm], theta, phi [rad], p [MeV/c]. Each voxel is its cell centre, with the earliest t, dE and dX summed, and theta/phi/p taken from the earliest step, i.e. the particle entering the cell (voxelize_truth). The generic voxelizer only merges columns 3-5, and left theta, phi and p at zero in every voxel. Also: - particles gains frag_inst_id, the instance a fragment row belongs to. groups/ is kept. - truth/, points/hit_index and points/true_x_shift are retired, and so is the pixel "deposits" reader option. reader=jaxtpc_pixel means the sensor-masked hits, and its cross-group settings live in conf/readout/jaxtpc_pixel.yaml (min_pc_size 63, distance_threshold 6.2, check_group_ownership false), loaded after _self_ so that a flag still wins. - The run records the pixel geometry as the root attribute pixel_geometry (JSON, one entry per volume, with x_anode_mm and mm_per_tick), so a viewer can place a hit from its tick and pixel index. - points/flat is stored in column chunks (16384 x 1) with bitshuffle + LZ4 (io.points). On a pixel-hit file this is 45.6 MB against 78.7 MB as rows. energy_mantissa_bits is an optional rounding, off by default. - pysupera-hits-subset now decodes every CSR entry into hit_tick, hit_charge and hit_wire or hit_py/hit_pz, in CSR order, so row k lines up with hit_labels entry k (--centres-only keeps the old contents). - pysupera-truth-subset is removed together with vis_truth.html. points/flat is now the truth, so the companion file had nothing left to add. Two speed-ups in preproc.py, which also carries the voxelizer change above. Both leave the output bit-for-bit identical, checked on 10 events of each mode (65/79/55 datasets compared): - _unique_rows: the voxelizer's np.unique(keys, axis=0) packs (particle, vx, vy, vz) into one int64 and sorts in 1-D, in the same order. On 3.3M rows this takes 0.14 s instead of 2.46 s. Keys wider than 63 bits fall back. - ScipyDefragmenter._get_labels_many: one KD-tree and one connected- components pass for all clouds, not one per particle (14k per pixel event). Clouds are kept apart by a fourth coordinate k*(eps+1). Within a cloud the distances are unchanged, and labels are shifted to start at 0, so they equal the per-cloud labels exactly. Over 10 pixel events (3 alternating repeats), voxelize went from 21.6 to 6.0 s, defragment from 10.1 to 4.0 s, and wall time from 68.3 to 49.6 s (-27%). Wire is -18% and EDepSim about -13%. README: the output format section is rewritten for 3.3.0, and so is the browser-viewer section, now one viewer per run mode (the viewers themselves are the next commit). Tests: 460 passed, 1 skipped.

The viewers were written for 3.2, where pixel hits were in points/flat and truth/ sat beside them. Under 3.3.0 points/flat is the truth, and the hits stay in the JAXTPC hits file, labelled by hit_labels/. So each viewer reads the output and, for JAXTPC, a hits subset (pysupera-hits-subset), which must be the same file the run read. EDepSim only vis_force.html ?file=out_gz.h5 JAXTPC wire vis_hits.html ?out=out_gz.h5&hits=hits_small.h5 JAXTPC pixel vis_drift.html ?out=out_gz.h5&hits=hits_small.h5 vis_force.html: only needed fixes to its format probes. The fragment level and geant4_trackid are now detected from the 3.x particles table (with them wrong, the log said "fragment level: no"). vis_hits.html (wire), rewritten: 3-D points/flat on the left and one 2-D image per hit_labels plane (wire vs tick) stacked on the right. Hovering a voxel lights its hits on every plane. Hovering a hit lights its voxels and its hits on the other planes. The link is the label itself (instance, fragment or interaction), so the views cannot disagree. Checked on event 0 of the sample: 27,762 voxels and 1.66M hits. An instance in volume 1 lights only volume 1's planes. vis_drift.html (pixel), rewritten: truth at true x on the left, and every hit on the right at the x inferred from its tick with t0 at the beam time, placed using the output's pixel_geometry attribute. The cameras are synced, both views are framed on the union, and hovering either side lights the object in both. Colour by energy (dE and ADC on separate ramps), semantic type, fragment/instance/interaction or PDG. There is a percentile threshold per pane, and an off-sensor toggle for the unlabelled hits (drawn grey). Event 0: 24,572 voxels and 3.61M on-sensor hits. Instance 29 lit 4,015 voxels and 475k hits. Each viewer puts up a banner when given the other readout's file, a file without hit_labels, or labels and hits of different lengths. vis_truth.html is removed (with pysupera-truth-subset in the previous commit), since vis_drift's left pane now shows the truth it was built for.

…ling Summary: the README had fallen behind the code in its defaults, point columns, pixel hit-mode numbers and output sizes, and said nothing about performance. This fixes those, and adds a performance section (per-mode timings, the checker comparison, and the before/after for the preproc speed-ups: pixel -27% wall, output identical), plus algorithm notes for voxelization and batched defragmentation. Documentation only; no code changes. A pass over the README against the current algorithms, defaults and data structures, beyond the 3.3.0 sections of the previous commits. Corrected: - Config defaults: distance_threshold 5.2 (6.2 for pixel), verbose false, checker cpu_single (the table said gpu), min_pc_size 5 (it said -1), and the voxelize rows added. - Preprocessing is three stages in order (merge duplicates, voxelize, defragment), not two. Merge duplicates sums dX rather than taking a max, and is skipped when the voxelizer runs. Defragmentation only examines kShower/kDelta/kMichel/kLEScatter by default. - Point columns: six for EDepSim only, nine for JAXTPC (theta, phi, p added), named in points/columns. - What the energy column holds: in points/flat it is always true dE. The point_source / hit_energy / hit_x_from attributes describe what the partitioning ran on, not what is stored. The hit_x_from section says the same. - Pixel hit mode after the sensor mask and the uint16 group_sizes fix: min_pc_size 63 from an 11.9 mm track width (the old 12.8 mm / 85 table was measured on unmasked, mis-decoded hits); groups straddling a split are 70 of about 574k over 10 events (0.01%), not 1.84%. - Output sizes per 10 events, remeasured with real theta/phi/p (the old numbers were taken when those columns were all zero and compressed away): pixel 16.3 MB (points/flat 3.8, hit_labels/ 11.1), wire 9.7 MB (5.7, 2.2). - Hit provenance points to hit_labels/ first, with groups/ as the group-level view behind it. The viewers named are the current ones. - Training data: torchdata reads points/flat, so for JAXTPC it loads truth voxels. A hits-file + hit_labels/ loader is not written yet. - Stale mentions of a pixel "deposit mode" are gone. New: - Compute Backends: the checker comparison. partition_combined on 3 pixel events took 3.4 s with cpu_single, 7.8 s with cpu_multi, 16.9 s with cell_hash_cpu_multi and 24.9 s with numba, so the default is the fastest. - Pipeline performance: wall and per-stage times per run mode (10 events: EDepSim 19.7 s, wire 29.0 s, pixel 49.6 s), the before/after table for the two preproc speed-ups (pixel voxelize 21.6 -> 6.0 s, defragment 10.1 -> 4.0 s, wall -27%), the profile that found them, and a note that identical runs on this machine varied by up to 35%. - Algorithms: Voxelization (the packed-key unique) and Defragmentation (the batched connected components, and why the labels come out identical). All internal links resolve.

verify_volume_geometry took its reference-tick anchor from the intercept of a least-squares fit of truth x on (1, tick, t0). When every deposit shares nearly the same non-zero t0, the t0 column is the constant column scaled, the fit cannot separate intercept from t0 slope, and it splits them arbitrarily. The check then subtracted t0 with the *stated* velocity from an anchor that had absorbed a different one, and rejected a correct geometry by hundreds of mm. Seen on test_doraemon_2026_v0.1 (run 93, event 0, volume 1): x off by 634 mm median with y, z within 1.1 mm and mm/tick exactly the stated 0.8000, while the fit claimed v = 0.23 mm/us, dt = 3.5 us. t0 spans 18-1300 us across the event, but the 932 groups that reach volume 1's readout are all one interaction: t0 = 346.8 us, spread 0.047 us. On those same pairs the new check passes with x within 0.20 mm. Synthetic groups on the exact stated geometry reproduce it: t0 = 460 +- 0.01 us fails by 97 mm; t0 ~ 0, or 460 +- 5 us, passes. The local test batches all have t0 ~ 0, which is why it never showed. The anchor is now the median over groups of the tick each would have had at the anode at t0 = 0, by the stated constants -- a one-parameter offset, so a wrong pitch, tick length or drift direction still fails as before. The fit is also no longer trusted for a velocity it cannot see. A new _fit_drift computes the standard error of the t0 slope; above 5 % relative, _measure_volume reports v and dt as NaN, with velocity_measured False and the t0 spread. The error message, the run summary in _run.py and the "config is missing" hint in format_jaxtpc.py say "not measurable" instead of printing a number that looks like a measurement. calibrate_volume (reader.pixel_geometry_from_fit=true) had the same degeneracy, and on the same event it fails too: the arbitrary split gave a t0 slope of the opposite sign to the drift, which the code assumes cannot happen, so x came out 2 * 0.229 mm/us * 346.8 us = 159 mm off. Where the signs happen to agree it would instead pass and store whatever drift velocity, time step and reference tick the split chose. It now uses the same _fit_drift test. When the t0 slope is not measured it takes the stated velocity (or, failing that, the stated time step) from config or the sensor file and refits only the intercept and mm per tick. With neither stated it refuses, naming what to set. A stated value is only a fallback: when t0 varies, the fit still wins. The report records which via velocity_source. Tests: a correct geometry with a shared 460 us t0 passes verification for both drift directions (fails on the old anchor). A wrong tick length is still rejected. Calibration refuses without a stated drift (did not raise before), recovers v, dt and the reference tick from either stated value in both directions, and prefers its own measurement when t0 varies.

The empty-image guard stops a run on any event with no points, because that is nearly always a wrong seg/inst path. test_doraemon_2026_v0.1 has a second cause. The readout window is exactly one full drift (2701 x 0.5 us = 1350.5 us) and interaction times run to ~1300 us, so an event whose interactions all happen late reaches the anode only after the last tick. JAXTPC zeroes that charge and the event has no hits at all. Run 725 event 11 is one: 42,805 deposits with t0 860-1137 us, earliest arrival 1369 us in volume 1 and 1996 us in volume 0, 0 hit groups. It is the only such event among run 725's 100 (others: 2,199 to ~47,000 groups); run 93 has none, though its sparsest event keeps only 18-25 % of its charge. Supera job 40076964_16 halted on it. The reader already counts, per event, the deposits that land before and after the window. When an event is empty and every one of its deposits is outside the window, the run now prints that, writes the event out empty -- so event numbering stays aligned with the input -- and counts it in the hit summary ("Events wholly outside it"). Every other empty event still halts unless allow_empty_image is set. The decision is _wholly_outside_window, kept apart so it can be tested without a run. Checked on run 725, events 0-12: the run completes, event 11 is written with 0 points and 0 particle rows, and the other twelve are written as usual (6,200 to 31,296 points each). The production-side fix -- a longer window or a narrower t0 range, so that a typical event stops losing half its charge (50.1 % of deposits past the window over these 13 events) -- is separate.

The stated pixel geometry was verified once, on the first event, for every volume at the same time. A volume with fewer than 8 groups there stopped the run. With interaction times up to ~1300 us against a 1350.5 us window, a volume can be empty in an event: every deposit in it reaches the anode after the last tick. Run 1045 event 0 is one -- volume 0 has 22,241 deposits at t0 1104-1111 us, earliest arrival 1601 us, and no hits, while volume 1 has 2,591 groups. Supera job 40110241_11 stopped on it after 4 s. It is the only such event in the file. Had event 0 had no hits anywhere, the job would have died here before reaching the empty-event handling of the previous commit. Stated constants do not depend on the data, so they are now applied from the first event, and each volume is checked on the first event that gives it MIN_GROUPS (8) groups. Nothing about the events before the check changes once it runs. A wrong geometry still stops the run -- on the event where the check runs, rather than on event 0. Until checked, a volume's entry in pixel_geometry_report is None. The run summary now reads the report at the end, says which event each volume was checked on, and warns (it does not stop) about a volume that no event in the file gave enough groups. The fit path (pixel_geometry_from_fit) is unchanged: it needs data to produce the geometry at all. Checked on run 1045, events 0-2: the run completes; volume 1 is checked on event 0 (2,591 groups, residual x 0.22 y 1.07 z 1.06 mm) and volume 0 on event 1 (11,047 groups, x 0.22 y 1.07 z 1.08 mm), both measuring v 1.6003 mm/us against the stated 1.6. Event 0 is written with 1,955 points, all from volume 1. Tests drive _check_pending_volumes on a reader with synthetic pairs: a volume empty in the first event is checked on the next, and checked only once; a wrong geometry is still rejected when its check runs; a volume with MIN_GROUPS - 1 groups waits.

…tics LE now marks blobs: a point is LE when its fragment is a kLEScatter fragment, a cluster with no trajectory. A low-energy scatter that AbsorbLEScatter merged into a track or shower fragment has become part of that trajectory, so it is no longer LE, although the particle's own sem_type stays kLEScatter. An LE fragment grouped into an instance through parentage stays LE. points/flat, hit_labels/*/is_le and groups/is_le all follow this one rule (layout.fragment_le). That rule only holds if proximity is the one way an LE particle enters a non-LE fragment, and it was not: PhotonDecay merged every e+- child of a photon into the photon's fragment, Compton electrons and photoelectrons included -- 23-31% of an event's points, many of them isolated blobs metres from anything (3.3.0: 20% of LE points sat over 100 mm from their fragment's non-LE points). A shower starts at the pair production, so: - PhotonDecay merges only the kConversion pair (and the photon's own split-off pieces), and promotes the photon to kShower only then. - merge_em_showers lets LE fragments join a shower instance but never head one, and lets a Compton electron or photoelectron climb to its photon only when that photon is itself inside a shower. A top-level photon's Compton electrons, from before it converted or instead of converting, are their own instances; the Compton blobs of a brem photon inside a shower join it, as LE. LE is 22-29% of points and 6-8% of the energy on three sample events per mode. Photon fragments that were several distant blobs go away: fragments that are not one connected blob at distance_threshold drop from 182 of 2366 (3.3.0, wire, 10 events) to 15 of 1005 (wire) and 9 of 1484 (EDepSim). Two more rules put low-energy pieces where they came from: - AbsorbLEScatter, when an LE partition touches several neighbours, takes the one holding its nearest ancestor if that is among them, instead of whichever the proximity search listed first. A delta ray goes to its track rather than to a shower crossing it: in pixel mode 505 of 523 LE deltas now join their parent's fragment, against 454. The partitioner publishes its live rep_lookup so post_filter can ask. - A stand-alone piece that defragmentation cut off a kTrack, kMichel or kDelta particle joins that particle's instance, and its fragment takes the particle's fragment type, so it is not LE (merge.find_split_pieces, relabel_split_pieces, attach_split_pieces). A piece proximity absorbed elsewhere stays there. Tracks are not defragmented by default, so this acts on delta and Michel pieces unless preprocessor.sem_types says otherwise. Since a fragment is wholly LE or not, it is one run: frag_pc_le_start/end are gone. An instance is three fenceposts, inst_pc_start <= inst_pc_le_start <= inst_pc_end, so inst_pc_le_end is gone, and inst_pc_end now ends the whole instance where in 3.x it ended the non-LE run. The column keeps its name but changes its meaning, so the reader takes 4.x only and refuses a 3.x file with a note to regenerate it. The 3.0 compatibility code goes with it. _check_ranges also refuses a fragment that mixes the two sides. Kinematics, all float64 (M, 4): *_pos [x, y, z, t] (mm, ns), *_mom [px, py, pz, E] (MeV, E total), plus mass. float64 because E - mass would lose a slow nucleus's KE in float32: about 2 keV of resolution for 40Ar. - mass, start_* and end_* are the row's own Geant4 particle's, joined on geant4_trackid, so a split piece shares its track's. Geant4 records a converting photon's end momentum as zero; end_mom is its momentum just before converting instead -- start less its Compton electrons', exact up to binding energy (Rayleigh scattering is not seen). - first_* and last_* describe an instance, the thing a model reconstructs, and sit on instance rows only. pysupera.kinematics makes them once the instances exist: each Geant4 step goes to the particle holding its deposit (the voxel map with EDepSim, the segment's group with JAXTPC, as the truth cloud does), and of an instance's counted steps -- visible with JAXTPC, all with EDepSim; root attribute kinematics_steps -- those of Geant4 particles with a non-LE member in it count, so a blob absorbed from an unrelated particle cannot be the last point. A conversion pair counts however small; an instance whose non-LE particles left no counted step uses all of its steps. first is the start of the earliest step, last the end of the latest: a step record holds its midpoint (median residual 0.0005 mm, against 0.097 mm for the start-of-step reading), so the ends move dx/2 along the step, and t by dx/2 / (beta c). - A kShower instance has no last point (NaN). When its photon's pair made the first counted step, first_mom is the pair's summed creation momentum: the photon's energy exactly, its direction up to the nuclear recoil (median 0.56 MeV against end_mom). When the readout missed the pair, the first point is a later product's, with that step's momentum. The readers hand over their step table (last_steps) instead of first/last; the JAXTPC reader's step-to-segment matching moves into _matched_steps, shared by the truth cloud, the visibility mask and the step groups. torchdata labels the [inst_pc_le_start, inst_pc_end) run as LE, and its instances table carries the kinematics. The three viewers read 4.x, take LE from frag_sem_type, and show start/end and, on instance rows, the first and last point with |p| and E. vis_hits' "hide LE" box becomes an LE selector (show all / hide LE / hide non-LE), and vis_drift gets the same one. They parse under node --check but were not run in a browser for this change. README updated throughout. Checked on three events in each mode (EDepSim, JAXTPC wire, JAXTPC pixel): no overlapping ranges; every instance row with points has a first point and no other row has one; no shower has a last point; at the conversion, pair-sum E matches the photon's end_mom E to 1e-7. Tests: 511 passed, 1 skipped.

Every supera_pixel job of test_doraemon_2026_v0.1 failed with "non-LE particle(s) have no instance after step-2 merging". In edepsim_j001995-001995.h5, event 50, a stopped mu+ (36 visible points, so kLEScatter) was absorbed by proximity into the shower of a photon that its own Michel positron radiated. The walk finds a fragment's parent through whichever fragment holds the parent particle, so the Michel's parent fragment was that shower, and the shower's was the Michel. find_initiator stops on a fragment it has seen, so each of the two named the other as the initiator. The fragment that headed no group of its own disappeared with its members (the Michel and an LE electron), and 35 member ids came out twice. Initiators are now followed to a fixed point, and a loop of them shares one instance, headed by its lowest-id member that may head a shower. A chain without a loop groups as before. Checked: all 100 events of that file run, with no member lost or duplicated by merge_em_showers or attach_split_pieces, and a test replays the case. Tests: 512 passed, 1 skipped.

kPrimary was meant to keep a generator primary out of kLEScatter: it maps to kTrack, or kShower for an electron, positron or photon, whatever the size. But _get_interaction_type marked a particle kPrimary only when track_id == parent_track_id, and EDepSim gives a primary a parent of -1, so nothing was ever kPrimary. Instead, primaries were classified by proc_start/subproc_start, which mean nothing for a primary. On the production file edepsim_j001995-001995.h5, 1,879 of its 5,480 primaries came out kLEScatter (816 photons, 449 pi0, 364 protons, 82 muons), and 61 primary electrons and positrons came out kUnknown (kInvalidProcess). The 1.7 MeV mu+ behind the step-2 failure of event 50 was one of them. The readers now read the input's primary table (reader.primary_key, primary/geant4), and a particle in it is kPrimary. Neutrons and nuclei keep their own rules, as before. The table is checked against parent_track_id == -1 every event (primary_mask), and the readers expose the result as last_primary_check. A file without the table falls back on parent_track_id == -1. Of 5,480 primaries in that file, 1,940 change type (all 1,879 LE ones, and the 61 kUnknown ones become kShower); on data/kazu.h5, 206 of 526 with truth deposits and 334 with pixel hits. No other particle changes. The table and parent -1 agree on every event of both files.

A warning printed among a production job's log lines is not seen. So each check that lets the run go on now records what it found under a fixed integer code (pysupera.issues.Issue), with the event, the Geant4 track and a message. The output stores the record as issues/: event, code, geant4_trackid, message, and a codes/ table naming every code. The root attribute n_issues holds the total and issues.attrs the count per code. The group is always written, so an empty one says the checks ran and found nothing. The end of the run log lists what was found, and EventStoreV3.issues() reads it back. Codes are never renumbered. 1 PRIMARY_TABLE_MISMATCH the primary table and parent_track_id == -1 disagree about a track; one row per track 2 NO_PRIMARY_TABLE the input has none; parent -1 taken (once) 3 EMPTY_IMAGE an all-empty event that allow_empty_image let through, every one now, not only the first 4 OUTSIDE_READOUT_WINDOW an event written empty because nothing reached the readout window dprod reads the group into each task's summary. The README documents the group, and the primary rule of the previous commit. Checked: the production file and the three local samples run with an empty group; a run on a synthetic file records each of codes 1-3. Tests: 531 passed, 1 skipped.

…on's pieces stay apart A primary is the interaction's generator particle, whatever it is. But in _get_interaction_type the neutron (pdg 2112) and nucleus (pdg > 1e9) rules came before the primary one, so a primary neutron was kNeutron, hence kLEScatter, and a primary nucleus kNucleus. The primary rule now comes first: such a particle is kPrimary and so kTrack, like a primary pi0. A neutron or nucleus that is not a primary keeps its rule. A primary neutron usually deposits nothing itself, leaving an instance with no voxels; those are common already (primary pi0, and protons and pions that leave no deposits). But some physics lists put sub-cut recoil energy on the neutron's own track, and find_split_pieces sends the split-off pieces of a kTrack particle back to its instance. That would gather blobs from metres apart into one instance. A neutron's pieces are now never sent back, whatever its semantic type (_NOT_A_TRAJECTORY_PDG). The production file edepsim_j001995-001995.h5 has no primary neutron or nucleus. On its 100 events (pixel), every particle's interaction type, semantic type, fragment and instance are the same as before this change.

…s a kDelta row With pixel hits a low-energy delta ray often cannot be told apart from the muon it came off: diffusion leaves its hits on the track or hugging it. It still headed its own fragment and instance, so the labels asked for something the image does not show. delta_min_standoff (mm, 0 = off) merges such a delta ray into its track. It applies to a non-primary kDelta fragment whose parent sits in a track-like fragment and which holds nothing else but LE and the delta ray's own pieces. The delta ray is merged when none of its voxels is that far or more from the nearest voxel of the track's fragment. It is measured on the clouds step 1 groups, so for pixel on the voxelized hits. Only the points move, into the track's fragment and so into its instance. The cut is distance and not a count of voxels off the track. A delta ray running along the track can have many voxels beside it and still be invisible, while a short one that leaves the track is plain to see. In a hand scan of event 0 (sdf_data_03, pixel, 3 mm voxels), how far each delta ray of muon 74 got from it matched the call for all 11: hidden 482 511 447 525 519 468 481 0 - 7.3 mm visible 507 479 461 480 9.0 - 20.3 mm A count of voxels off the track did not separate them: 507 (visible) has 19, the same as 481 (hidden), and 519 (hidden) has 21. readout/jaxtpc_pixel.yaml sets 8 mm, about two 4.32 mm pixels. Over 3 events that merges 25-30 delta rays an event. Wire and edepsim leave it off. On the true deposits the same delta rays are resolved (519 gets 12.4 mm from the muon), and their delta rays are decided as before, by min_pc_size and AbsorbLEScatter. A merged delta ray keeps a row, whichever rule merged it. A delta ray smaller than min_pc_size is typed kLEScatter, and AbsorbLEScatter puts it in the track it touches. Before this it got no row unless it lay on someone's ancestry, and then as kLEScatter. Now any delta ray in its parent track's fragment, merged by either rule, is written as sem_type kDelta, with its own pc_start/pc_end inside the track's fragment range and frag_id/inst_id -1. So sem_type == kDelta selects every delta ray, kept or merged. Its steps do not count toward the track instance's first/last point (buried_in_track), as its LE type kept them out before. One merged into anything else, or left as an LE blob, stays kLEScatter. This only adds rows: +5% (edepsim) and +9% (pixel) particle rows over 3 events, about 490 and 340 an event. Points, groups, hit labels, fragment/instance heads and instance first/last are identical with and without it, checked on 3 events of both. The run report prints both counts ("Delta rays into track", "LE delta rays in track"). The README has both rules and the scan table.

… colouring and filters Following one object through a busy event by hover alone is fragile: the highlight moves as soon as the cursor does. All three viewers now let you pin objects by clicking: click pin the object under the cursor; click more to pin several together, click a pinned one to unpin. While anything is pinned, hovering changes nothing. x / Esc / button unpin everything h / checkbox hover highlighting on or off p / checkbox left-drag pans instead of rotating; right-drag and shift-drag always pan, and so do the arrow keys (shift+arrows in vis_force) A drag never pins, so rotating does not change what is lit. vis_force.html pins instances in the shader (uPins, up to 8 at once). Its toolbar gains HOVER, UNPIN and PAN, and its help lists the keys. vis_hits.html (wire) was behind vis_drift.html, so it is rebuilt on vis_drift's point-set model, keeping its six plane images. New in it: - colour by energy/charge, semantic type, fragment, instance, interaction or PDG, separate from what a hover highlights (highlight: instance, fragment or interaction; it was "link by") - the semantic-type strip, turning each type on or off in every view - log scale, and the show >= percentile cut - a legend The energy ramp and the cut are one for the voxels (dE, MeV) and one shared by all planes (hit charge), so the planes compare. LE filtering now hides points as it draws instead of rebuilding the event, and the pickers skip whatever is hidden. Colours are worked out once per colouring or filter change, not on every hover, since a plane holds ~280k hits. Plane headers give how many hits are shown and how many are highlighted. Checked in the browser on 3 wire events (sdf_data_03_wire, 29k voxels and 1.70M hits in event 0): - colour by semantic type, with shower and LE hidden - hover on a hit (instance 69, a track) - pinning 69 and 76 together, lighting both in the 3-D view and on all six planes - x, h and p - energy colouring with the 50% cut There were no console errors. vis_drift and vis_force were checked the same way. The README documents vis_hits' controls, with one Highlighting table for both JAXTPC viewers.

…d id A production writes one EDepSim file per job, and JAXTPC can merge the events of several jobs into one file (--merge-output). pysupera took one EDepSim file, matched JAXTPC events by position (event_NNN), and kept no record of which input event an output event was. io.input_path may now be a list. Each reader.jaxtpc_*_path is then a list of the same length (one JAXTPC file per input), or one path all inputs share (a merged JAXTPC file). config.reader_inputs works out the pairs. build_reader builds one reader per input and chains them (readers/chain.py, ChainReader). Events are written in input order. The chained readers share one pixel geometry, so each volume is checked once per run. A single path works as before. JAXTPC format 2 numbers its groups entry_NNNNNN by position in the file only, so the reader pairs each EDepSim event with its entry by identity: - by job_id/job_event against the EDepSim file's /dprod/events, when both have them; - otherwise by source_event_idx, and by the file name in config/source_files (an error if two inputs share a name). An EDepSim event with no entry stops the run, saying so. The step and sensor files must hold the same event under that entry, or the run stops. Format-1 files (event_NNN) are read as before. Format 4.1.0 adds, and changes nothing else (4.0.0 files read as before): inputs/{edepsim,jaxtpc_step,jaxtpc_hits,jaxtpc_sensor} one row per input, "" when not used events/input which input each event came from events/input_event its index in that EDepSim file events/jaxtpc_entry its JAXTPC entry, -1 without one and, when every input carries /dprod/events, its identity copied unchanged: events/event_id, job_id, job_event event_id is the same in the EDepSim and JAXTPC files, so it pairs an event across all of them. If only some inputs carry the table, the run stops. EventStoreV3.event_ids() and .event_sources() read these columns (pysupera/identity.py). pysupera-hits-subset also takes entry_ groups and keeps each one's identity attributes. tests/test_inputs.py has 12 tests: input lists, the chain, entry pairing (by id, by name, a missing entry, mismatched files, format 1) and the columns written. Full suite: 559 passed, 1 skipped. Checked end to end on inputs built from sdf_data_03_wire. A.h5 holds EDepSim events 0-2 as job 1995, and B.h5 holds events 3-5 as job 1996. A merged step/hits file in format 2 stores B's entries first. Running [A, B] against it gave the right input, input_event, jaxtpc_entry 3-5 then 0-2, and the event_ids. Every other dataset is identical to a single-input run over events 0-5 of the original files. pysupera-repack keeps the new groups.

…ouching search Three changes take 16-26% off a run's wall time. The output does not change by a bit. provenance: the true deposits behind the hits were voxelized per particle, about 7,000 np.unique(axis=0) calls an event. voxelize_truth_many does all particles in one pass, with the packed-key sort of the voxelizer on (particle, cell). Within a particle the cells and the steps in each cell come in the same order, and dE/dX are summed in the same order, so each particle's rows are what voxelize_truth gives it alone. That function stays as the definition, and a test holds the two equal, bit for bit. proxck: CPUSingleThreadChecker built a KD-tree, a float32 cloud, a bounding box and a centroid for every particle at start (about 15,000 an event). Partition-level proximity, the default, builds its own trees and read none of them. They are now made on first access, through a dict that fills a missing key. The particle-level path gets the same objects when it asks. AbsorbLEScatter: the touching (LE rep, other rep) pairs came from query_ball_tree and a Python loop over every point pair, about 1M an event in pixel. One cKDTree.sparse_distance_matrix and NumPy dedupe replace them. The order of the candidates decides merges. When an LE partition touches several reps and none is its home, the first pair wins, and "first" is the iteration order of the candidate set, which depends on insertion order. query_ball_tree returns each point's list sorted. So the pairs are inserted by LE point, then other point, each rep pair at its first occurrence, as before. The set and the candidate list are therefore identical. A test checks this against the old code on a grid where many pairs are exactly D apart. On 10 real events per mode, all 117 calls gave the same list, element for element. Checked on sdf_data_03 (pixel, EDepSim only) and sdf_data_03_wire, all 100 events of each. Every dataset and attribute is identical to the code before (85, 99 and 75 datasets). Timings are for 20 events, mean of 2 alternating serial runs, which agreed within 2%: pixel 131.7 s -> 109.7 s (-17%) wire 80.5 s -> 59.4 s (-26%) EDepSim 69.8 s -> 58.4 s (-16%) partitioner_init drops from 5-9 s to under 1 s in every mode. The truth voxels halve wire's write_event (29.9 -> 14.0 s) and cut pixel's from 23.9 to 13.6 s. AbsorbLEScatter's search is 33% faster in pixel on the same calls, and 12-18% faster in the others. Wire's partition_combined rose from 17.4 to 18.4 s, consistently but unexplained. Full suite: 562 passed, 1 skipped (3 new tests). The README describes each change where the algorithm is documented and adds the timing table under Pipeline performance.

drinkingkazu added 3 commits March 7, 2026 10:17

edep-sim interface added. Auxiliary variables for Particle also added…

bd27f18

…. Voxelization of edep-sim output is added in the form of pre-processor.

ready to interface edep-sim H5 files. visualization options tweaked t…

f2b55b4

…o handle very busy event topologies (which we plan to use for FM data generation). Algorithms (LEScatter merging ones) optimized for speed for similarly complex dataset (180s=>2s process time per event)

minor README update

afe1c57

drinkingkazu requested a review from Copilot March 10, 2026 16:08

Copilot started reviewing on behalf of drinkingkazu March 10, 2026 16:09

View session

Copilot AI reviewed Mar 10, 2026

View reviewed changes

drinkingkazu added 24 commits March 10, 2026 09:25

a bug fix: process_type enum assumed to start from 0 in some parts of…

7234f95

… the unit tests, a remanent of C++ version of Supera. Fixed. Also process_type use was mixed with bare integer values, which is updated.

add input point <=> voxelized points mapping

839d941

supporting jaxtpc stuff

d174bef

drinkingkazu added 28 commits September 20, 2026 22:50

Navigation Menu

Repository navigation

Interface edep-sim h5 output - #7

Interface edep-sim h5 output#7
drinkingkazu wants to merge 55 commits into
DeepLearnPhysics:developfrom
drinkingkazu:develop

drinkingkazu commented Mar 10, 2026

Uh oh!

Copilot AI left a comment

Uh oh!

Copilot AI Mar 10, 2026

Uh oh!

Copilot AI Mar 10, 2026

Uh oh!

Copilot AI Mar 10, 2026

Uh oh!

Copilot AI Mar 10, 2026

Uh oh!

Reviewers

Uh oh!

Assignees

Labels

Projects

Milestone

Development

Uh oh!

2 participants

	process_type = InteractionType(process_type)
	process_type = InteractionType(process_type + 1)

Conversation

drinkingkazu commented Mar 10, 2026

Uh oh!

Copilot AI left a comment

Choose a reason for hiding this comment

Pull request overview

Reviewed changes

Uh oh!

Copilot AI Mar 10, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 10, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 10, 2026

Choose a reason for hiding this comment

Uh oh!

Copilot AI Mar 10, 2026

Choose a reason for hiding this comment

Uh oh!

Reviewers

Uh oh!

Assignees

Labels

Projects

Milestone

Development

Uh oh!

2 participants