ILE: label the exported time grid with the spacing the likelihood steps by (minimal, export path only)

[oshaughnessy-junior]

DRAFT. Minimal, export-path-only fix for a time-labelling defect that is present in
every GWTC-6 production run. Opened as a draft against oshaughn:rift_O4c because this is
the branch Reviewed-RIFT-20260807.sif (tag 0.0.17.13) is built from and it needs a
minimal path forward. No review requested.

The defect

resample_samples() builds the time-marginalization grid as

tvals = linspace(-t_ref_wind, t_ref_wind, int(t_ref_wind*2/P.deltaT))

whose spacing is 2W/(N-1). The likelihood does not step by that.
DiscreteFactoredLogLikelihoodViaArrayVectorNoLoop uses only tvals[0] — it sets
tfirst = t_det + tvals[0] — and then reads lnLt[k] from the precomputed Q_lm buffer at
index ifirst + k, i.e. at detector time tfirst + k*deltaT. Labelling those values with
a closed-interval linspace stretches the reported time about the start of the window:

t_reported = t_true + (t_true - (event_time - W)) * delta
delta      = (1 + f)/(N - 1),   where  2*W*srate = N + f

At the production settings (--srate 4096, --data-integration-window-half 0.075):
N = 614, f = 0.4, delta = 0.0022838, so

• the exported geocentre time is +171 µs late at the window centre (+80 µs at −40 ms
  rising to +263 µs at +40 ms), and
• any exported time posterior is 0.23 % too wide.

delta is not 1/srate — it depends on the fractional part of 2*W*srate:

--srate	N	delta	bias at window centre
4096	614	0.0022838	+171.3 µs
8192	1228	0.0014670	+110.0 µs
16384	2457	0.00065147	+48.9 µs

Neither existing option works around it. --interpolate-time does not (the cubic stencil
still steps exactly one sample), and --srate-resample-time-marginalization does not (the
dense export grid is built as tvals[0] + arange(n)/srate_resample and inherits the
mislabelling).

The likely origin is someone making the integer grid-size rounding come out right:
int(2W/deltaT) and tvals[0] are both already correct, only the spacing of the labels
is wrong.

The fix

One line: tie the label spacing to deltaT, the step the likelihood actually takes.

tvals = -t_ref_wind + float(P.deltaT)*xpy_default.arange(int((t_ref_wind)*2/P.deltaT))

tvals[0] and len(tvals) are unchanged, so the likelihood reads exactly the same
Q_lm samples. This relabels the times; it does not move the likelihood.

Deliberately not the rift_O4d convention. rift_O4d (junior PR #162, issue #146) now
uses a centred grid, (arange(npts) - npts//2)*deltaT, so that t = 0 sits on the grid
and the JAX driver agrees sample-for-sample. That shifts tvals[0] by 0.2 samples, which
can round ifirst to a different integer and therefore change which data samples are
integrated. Keeping tvals[0] = -t_ref_wind here is what makes this patch provably
numerically inert, which is the right trade for the release branch.

Scope: why only this one site

There are seven grid sites in this file. This is the only one whose tvals values leave
the process as data (they become t_out, hence the exported t_ref). Checked on the
O4c file, not assumed — the branch lines have diverged:

• 4 sites (1583, 1632, 1738, 1782) feed DiscreteFactoredLogLikelihoodViaArrayVector
  / ...NoLoop, which consume only tvals[0] and len(tvals) and integrate with
  dx=deltaT. Relabelling them changes nothing at all, so a minimal patch leaves them.
• 2 sites (1553, 1835) feed FactoredLogLikelihoodTimeMarginalized, which integrates
  with dx=tvals[1]-tvals[0] (factored_likelihood.py:694). Fixing those would shift
  lnL by a constant -log(1+delta) = -0.00228 nats on the non-vectorized and ROM paths.
  Correct, but not numerically inert, so out of scope here.

Verification

Against the shipped likelihood on this branch (synthetic Q_lm carrying a feature at a
known detector arrival time, <h|h> cross terms zeroed so lnL(t) tracks the feature):

Q1  lnLt bitwise identical under the two grids : True   (max |dlnLt| = 0.000e+00)
Q2  time-marginalized lnL                      : -1.9389090783052292 under both, identical
Q5  npts at srate 1024/2048/4096/8192/16384    : unchanged (153/307/614/1228/2457)
    tvals[0] at every rate                     : unchanged

Recovered peak position of the injected feature:

t_true[ms]   legacy err[us]    fixed err[us]   predicted[us]
     -40.0            97.12            17.15           79.93
     -20.0           142.72            17.07          125.61
      -5.0           177.17            17.26          159.87
       0.0           188.37            17.04          171.29
       5.0           200.25            17.50          182.71
      20.0           234.05            17.05          216.97
      40.0           279.76            17.08          262.64

legacy − fixed matches the closed form above to better than 0.05 µs at every offset.
The residual +17 µs that survives is constant — it is the separate nearest-sample
detector-time snap for that sky position, not this defect.

--srate-resample-time-marginalization still lands on exactly 1/srate_resample:

srate_resample   n_dense (legacy -> fixed)   step exactly 1/srate_resample   stays <= tvals[-1]
       8192            1229 -> 1227                     yes                        yes
      16384            2458 -> 2453                     yes                        yes
      32768            4916 -> 4905                     yes                        yes

The dense grid loses a few points at the far edge because tvals[-1] moves from +W (a
label the likelihood never actually evaluated) to the last real sample. No coverage is
lost: the likelihood always evaluated npts samples stepping deltaT from -W, and the
old code merely mislabelled the last of them as +W.

On real data (recorded separately, not rerun here): re-running the extrinsic stage of
GWTC-6 S250331o over identical blocks of the final intrinsic grid with this correction
moves the exported posterior by −146.4 ± 22.6 µs against a −149.4 µs prediction fixed
before the run, with a null replicate (same arguments, different RNG) at +14.4 ± 23.8 µs.

How much does it matter in practice

JS divergence attributable to this defect alone, per event, over the 21 GWTC-6 events with
a bilby SEOBNRv5PHM counterpart (the bias is a deterministic map, so it is evaluated by
applying the inverse map to the production samples — no rerun needed):

median JS from the bug alone : 6.73e-04 bits
max    JS from the bug alone : 2.39e-02 bits  (S250331o)
median JS RIFT vs bilby      : 5.73e-03 bits
bug / total, median ratio    : 17.4 %

So it is negligible for poorly-localised events — median 6.7e-4 bits, below the
~2e-3 bit scale at which these comparisons are normally called indistinguishable, on the
broad (≈40 ms wide) events — and dominant for sharply-localised ones: 2.4e-2 bits on
S250331o (90 % width 1.65 ms), which is larger than the entire current RIFT-vs-bilby JS
for that event (1.2e-2 bits), and 1.1e-2 bits on S250213dg.

Tests

MonteCarloMarginalizeCode/Code/test/test_time_marginalization_grid_labels.py (13 tests,
no GPU, no data): the corrected grid's spacing, that tvals[0] and len(tvals) are
preserved, the closed-form bias at five offsets, the srate ladder, the width stretch, and
a drift guard that reads the shipped source line, evaluates it at five sample rates, and
fails if the export grid becomes a closed-interval linspace again.

Verified to fail for the right reason. Against the unpatched tree the two source-reading
guards fail and the other eleven pass. Mutation battery, with the unmutated control run
through the same scoring path:

control (unmutated)          : 13 passed
npts halved                  :  1 failed
npts + 1                     :  1 failed
origin sign flipped          :  2 failed
spacing doubled              :  1 failed
the original bug restored    :  2 failed

(An earlier revision of this patch introduced _npts_t on its own line; review caught that
the guard then never evaluated it and the "npts halved" mutant survived. npts is now
inlined into the guarded expression.)

Run with:

PYTHONPATH=MonteCarloMarginalizeCode/Code pytest -q \
  MonteCarloMarginalizeCode/Code/test/test_time_marginalization_grid_labels.py

Not addressed here

A second, independent defect in the same machinery — the nearest-sample detector-time snap
broadens the time posterior where the inter-detector delay is sharp (the constant +17 µs
above, and a +12.5 % width excess on S250331o). rift_O4d has --interpolate-time for it;
rift_O4c does not. Out of scope for a minimal release-branch patch.