Fixes The Beam Accounting

[jzennamo]

From Claude with love:

Bugs fixed in BNB beam-quality retrieval / FOM

Found while auditing BeamSpillInfoRetriever against ICARUS flat CAFs (2 files, 99 triggered spills, 3,021 subrun spills; Run 2 + Run 4). Numbers below are from those files.

High impact

1. Missing devices are silently used as real values (getFOM.cpp)
Every device is pushed into a one-element std::vector, so all the .empty() checks (TOR875 fallback, "return 2/3 when BPM data missing") are dead code. A device that fails to read arrives as -999 and goes straight into the position/angle extrapolation, so the beam is placed ~1 m off target, FOM ≈ 0, and the spill ends up with FOM = -999 in the CAF.

• Fix: treat -999/non-finite as missing. Return 2 (horizontal) or 3 (vertical) when the needed BPM or BPM offset is missing, as the code intended.

2. Failed BPM-offset queries throw away good spills (POTTools.cpp)
The offsets (E_*S, BNB_BPM_settings) are slowly changing setpoints, but the query fails fairly often. Combined with (1), 624 / 3,021 subrun spills (21%) with full intensity (median TOR860 2.7e12) got FOM = -999 only because an offset was missing.

• Fix: keep a per-job cache of the last valid offsets and use it when a query fails (ReuseLastBPMOffsets, default true). In endSubRun (ICARUS and SBND BNB retrievers), back-fill offsets for spills recorded before the first valid reading and recompute their FOM. Emulated on the data, this recovers 624/624; the remaining 414 -999s all have a BPM reading missing.

3. Multiwire profiles matched to the wrong spill (ICARUS/SBND retrievers)
The MWR↔spill matching has no maximum time difference and falls back silently to index 0 when there is no best match. In the data, *_spill_time_diff reaches 22–27 s, so profiles from unrelated spills were used for the beam width.

• Fix: the three duplicated loops are replaced by one tested function, sbn::pot::matchMWRToSpill() (new MWRMatching.h/.cpp, no art deps). It returns -1 for no match and rejects matches beyond MWRMaxTimeDiff (default 0.0333 s, half the 15 Hz period; <= 0 disables).

4. Pre-fit FOM uses an unrelated / uninitialised chi2 (getFOM.cpp)
The "pre-fit" FOM (database widths M876HS/VS, M875HS/VS) applies the chi2 cut from whichever multiwire profile fit ran last, which has nothing to do with those widths. With no multiwire data (all of Run 2), chi2x/chi2y are uninitialised (undefined behaviour), so whether the pre-fit FOM is accepted depends on stack contents.

• Fix: the pre-fit FOM uses only the width window plus validity of the database widths.

Medium

5. No intensity requirement. Beam-off spills (TOR860 ≤ 0) get FOM ≈ 1; 13 such spills passed FOM > 0.95. Also there was no TOR860→TOR875 fallback, and a -999e12 intensity gives a negative emittance → sqrt(<0) → NaN.

• Fix: use TOR860, else TOR875; return -1 if neither is valid and positive.

6. Out-of-bounds read on short multiwire arrays. The profile is used if size() > 0, but 96 elements are read.

• Fix: require a complete 96-channel profile.

7. Fit result never checked; histogram leaked. TFitResultPtr is dereferenced without checking status/null; the TH1D is leaked on the empty-profile path and registered in gDirectory under a fixed name.

• Fix: local histogram with SetDirectory(nullptr); processBNBprofile() now returns bool (status 0, ndf > 0).

8. assert contradicting the code below it (POTTools.cpp). makeBNBSpillInfo asserted every MWR device was non-empty, then handled the empty case; debug builds would abort.

• Fix: assert removed; empty / unmatched devices give an empty profile and spill_time_diff = -999.

9. SBNDBNBZEROBIASRetriever spill selection. BrokenClock(times_temps[i], …) tested the previous best (initially entry 0) instead of the candidate times_temps[k]. When no spill qualified, entry 0 was used anyway (out of range for an empty list).

• Fix: test k; if no spill is found, skip (the event still gets an empty product).

Low / cosmetic

10. Angles are slopes in mm/m (= mrad) but were passed through atan() as if they were tangents. This is <1% for typical slopes but wrong. Now uses the slope directly.
11. FillExposure.cxx treated a physical FOM = 0 (beam fully off target) as a failure (> 0.0) and used bitwise &. Now [0, 1] with &&.
12. π = 3.14159 in the target integral (normalisation +8e-7, contributes to FOM saturating at exactly 1). Now M_PI.
13. Comment fixes: MWRtoroidDelay is in s (not ms); multiwire coverage is 24 mm (not 24 cm).

Behaviour changes to be aware of

• FOM values change slightly (max |ΔFOM| ≈ 7e-4 vs. existing CAFs) from removing atan (10) and the pre-fit chi2 fix (4).
• New flag values from getBNBqualityFOM: -1 (no beam), 2/3 (missing H/V BPM or offset). FillExposure maps these to -999 exactly as before. The 100 + FOM convention for the nominal-width FOM is unchanged.
• New fhicl parameters MWRMaxTimeDiff and ReuseLastBPMOffsets have defaults, so existing fcls run unchanged. Both are documented in icarusbnbspillinfo.fcl / sbndbnbdefaults.fcl.
• Affects ICARUS BNB, SBND BNB and SBND BNB zero-bias retrievers, plus CAFMaker for both experiments. EXT and NuMI retrievers are untouched.

Testing

• getFOM.cpp compiled (-Wall -Wextra) against a minimal ROOT fit stand-in and run on all 3,120 spills. It agrees with an independent Python implementation to 3e-8 on every spill. The same Python code in "legacy" mode reproduces the FOM stored in the CAFs for 3,120/3,120 spills, so the reference was validated against production first.
• matchMWRToSpill() was checked against the original loop on 388,644 random device–spill cases. It is identical except for the intended changes (no index-0 fallback, time cut).
• Not yet done: build of the art modules / POTTools in a full mrb environment, a test job with run_icarusbnbinfo_sbn.fcl, and any SBND data.

Open items (not fixed here)

• Target multiwire (MMBTBB) never gives a usable width: median χ²/ndf ≈ 260, with a dead channel in the middle of the vertical peak and alternating hot/cold wires. The FOM always falls back to M875/M876. This could be a channel-mapping problem in MWRData::unpackMWR or the hardware; it needs a beam-instrumentation look.
• Check whether BNBSpillInfo::POT() guards against a failed TOR860 (-999e12) being added to the subrun POT.
• The zero-bias retriever writes per event, so it cannot do the end-of-subrun offset back-fill. Spills before the first valid offset in a job still lose their FOM there.

[nathanielerowe]

Some of these things are probably issues in develop as well, can you create a PR for that as well?