fix(series4): Thor/Micromate decoder is now per-sample exact

Verified against Thor's own CSV exports, which carry a per-sample four-column block beside every binary (CSV/<name>.IDFW.csv). Those 1,012 paired files were in the corpus all along; the decoder had been pinned to a superseded walker on the stated grounds that "Thor has no ASCII ground truth in the corpus and its geo scaling is separately suspect". Both premises were false. IDFW per-sample exact 39.1% -> 100.000% (1,057,536/1,057,536) IDFW files fully exact 0/153 -> 153/153 IDFW PPV median error -3.32% -> -0.002% IDFH within 2% of Thor PPV 51.1% -> 100.0% (858/858) prod IDFW, 8 units -3.3% -> -0.001% Four independent root causes: - Geo LSB was 0.0003, the 4-dp *display rounding* of the real 0.000310308 mistaken for the LSB, so every series-4 geophone sample read 3.3% low. Pinned to +-6e-11 by intersecting 991,415 rounding constraints; corroborated by the +-full-scale seed (+-32226) left in unwritten IDFH slots. IDFH had a separate, also wrong, 10.0/32768. - IDFH histograms were capped at 250 intervals: the segment validator required the interval counter's high byte to be zero, but the counter is a uint16 cumulative index, so every segment past interval 255 was rejected. Runs over ~4 hours lost their tail, often the peak. 540/858 corpus files affected. - Record mode 00 00 (raw int16, 10-byte header) was unhandled and fell through the dispatch, silently dropping each channel's first 512 samples -- the long-standing "loud events truncate" symptom. MODE_ABSOLUTE is now also accepted as a segment-0 preamble. - The body-offset search matched 00 02 00 *inside* record headers, selecting a candidate part-way down the chain and decoding a rotation-shifted body. It now anchors on record headers and takes the chain head (6 ms/file). Also fixes the separately tracked "UM-series decodes ~1000x low" bug. Series-3 re-verified unchanged at 14,338/14,338 exact after the shared waveform_codec change. Known open: 41/575 prod IDFW files (7%, mostly UM12947/UM20147) decode with unequal channel lengths and also fail metadata extraction -- a different header variant with no Thor export in the store. NOTE: this is a codec change; the Thor store owes a regeneration via scripts/backfill_thor_events.py. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01Ru8Lg9HkkYvX9VWWo65SmL
2026-09-10 18:06:14 +00:00
parent 91b9b4578c
commit 726c2ce1b5
8 changed files with 949 additions and 66 deletions
@@ -0,0 +1,237 @@
+"""Per-sample verification of the Thor / Micromate (series-4) IDF binary codec.
+
+Ground truth is Thor's own CSV export, written next to each binary by the
+Thor desktop application.  For waveforms the export carries a per-sample
+block of four columns (Tran, Vert, Long, Mic) in in/s and psi -- the
+series-4 equivalent of Blastware's ``_ASCII.TXT`` exports.
+
+The full-corpus harness is ``scratch/verify_thor_against_csv.py``; these
+tests pin the two constants that harness established so they cannot
+regress silently.
+"""
+from __future__ import annotations
+
+import csv
+import os
+import sys
+from pathlib import Path
+
+import pytest
+
+sys.path.insert(0, os.path.dirname(os.path.dirname(os.path.abspath(__file__))))
+
+from micromate.idf_file import (
+    _GEO_LSB_IPS,
+    geo_count_to_ips,
+    read_idf_file,
+)
+
+FIXTURES = Path(__file__).parent / "fixtures" / "thor-idf"
+IDFW = FIXTURES / "UM11719_20231219162723.IDFW"
+IDFH = FIXTURES / "UM11719_20231219162648.IDFH"
+
+GEO_CHANNELS = ("Tran", "Vert", "Long")
+
+# tests/fixtures/ is gitignored, so a fresh checkout has no sample data.
+# Skip rather than fail, matching test_idf_ascii_report.py.  To populate:
+#
+#   B="<thor-watcher>/example-data/THORDATA_example/THORDATA_example/UPMC Presby"
+#   mkdir -p tests/fixtures/thor-idf
+#   for f in UM11719/UM11719_20231219162723.IDFW \
+#            UM11719/UM11719_20231219162648.IDFH \
+#            UM13981/UM13981_20220207084555.IDFW \
+#            UM13981/UM13981_20220207183102.IDFH \
+#            UM13981/UM13981_20221202063059.IDFH; do
+#     cp "$B/$f" tests/fixtures/thor-idf/
+#     cp "$B/$(dirname $f)/CSV/$(basename $f).csv" tests/fixtures/thor-idf/
+#   done
+pytestmark = pytest.mark.skipif(
+    not FIXTURES.is_dir() or not any(FIXTURES.glob("*.IDFW")),
+    reason=f"Thor IDF fixtures not present under {FIXTURES}",
+)
+
+
+def _parse_export(path: Path):
+    """Split a Thor CSV export into (header dict, per-sample rows)."""
+    header, rows = {}, []
+    with path.open(newline="", encoding="utf-8", errors="replace") as fh:
+        for rec in csv.reader(fh):
+            if len(rec) == 2:
+                header[rec[0].strip()] = rec[1].strip()
+            elif len(rec) >= 3:
+                try:
+                    rows.append([float(x) for x in rec])
+                except ValueError:
+                    pass
+    return header, rows
+
+
+def _header_float(header, key):
+    return float(header[key].split()[0])
+
+
+@pytest.fixture(scope="module")
+def idfw_export():
+    return _parse_export(IDFW.with_suffix(".IDFW.csv"))
+
+
+# ─── The geo scale constant ────────────────────────────────────────────────
+
+
+def test_geo_lsb_matches_thor_quantisation():
+    """Thor's own export quantises geo samples to this LSB.
+
+    Derived by maximising exact-match count over 1,046,016 paired samples
+    (454 channel-events, 2 units); independently corroborated on 8
+    production units via their device-reported PPV.  The historical value
+    0.0003 read every series-4 geophone sample 3.3% low.
+    """
+    assert _GEO_LSB_IPS == pytest.approx(0.000310308, rel=1e-6)
+
+
+def test_geo_lsb_is_not_the_legacy_value():
+    # Guards against a revert to the truncated 0.0003 constant.
+    assert abs(_GEO_LSB_IPS - 0.0003) > 1e-6
+
+
+# ─── Per-sample fidelity ───────────────────────────────────────────────────
+
+
+def test_waveform_channel_lengths_match_export(idfw_export):
+    _header, rows = idfw_export
+    result = read_idf_file(IDFW)
+    for channel in GEO_CHANNELS:
+        assert len(result.samples[channel]) == len(rows), (
+            f"{channel} truncated: decoded {len(result.samples[channel])} "
+            f"samples, export has {len(rows)}"
+        )
+
+
+def test_waveform_samples_match_export_exactly(idfw_export):
+    """Every geo sample must reproduce Thor's exported value to 4 dp."""
+    _header, rows = idfw_export
+    result = read_idf_file(IDFW)
+    for index, channel in enumerate(GEO_CHANNELS):
+        decoded = result.samples[channel]
+        expected = [row[index] for row in rows]
+        mismatches = [
+            (i, geo_count_to_ips(c), v)
+            for i, (c, v) in enumerate(zip(decoded, expected))
+            if abs(geo_count_to_ips(c) - v) >= 5e-5
+        ]
+        assert not mismatches, (
+            f"{channel}: {len(mismatches)} of {len(expected)} samples differ; "
+            f"first three {mismatches[:3]}"
+        )
+
+
+def test_waveform_ppv_matches_export(idfw_export):
+    header, _rows = idfw_export
+    result = read_idf_file(IDFW)
+    for channel, attr in (
+        ("Tran", "transverse_ips"),
+        ("Vert", "vertical_ips"),
+        ("Long", "longitudinal_ips"),
+    ):
+        decoded = getattr(result.event.peaks, attr)
+        assert decoded == pytest.approx(
+            _header_float(header, f"{channel}PPV"), abs=5e-5
+        ), f"{channel} PPV disagrees with Thor's export"
+
+
+# ─── Histogram path shares the same scale ──────────────────────────────────
+
+
+def test_histogram_peaks_match_export():
+    header, _rows = _parse_export(IDFH.with_suffix(".IDFH.csv"))
+    result = read_idf_file(IDFH)
+    assert result.intervals, "IDFH decoded no intervals"
+    for channel, attr in (
+        ("Tran", "transverse_ips"),
+        ("Vert", "vertical_ips"),
+        ("Long", "longitudinal_ips"),
+    ):
+        decoded = getattr(result.event.peaks, attr)
+        expected = _header_float(header, f"{channel}PPV")
+        # Histogram peaks are stored per-interval, so the export's PPV is
+        # reproduced within one quantisation step rather than exactly.
+        assert decoded == pytest.approx(expected, abs=2 * _GEO_LSB_IPS), (
+            f"{channel} histogram peak {decoded} vs export {expected}"
+        )
+
+
+# ─── Regressions found 2026-09-10 ──────────────────────────────────────────
+
+IDFH_LONG = FIXTURES / "UM13981_20220207183102.IDFH"      # 719 intervals
+IDFH_SENTINEL = FIXTURES / "UM13981_20221202063059.IDFH"  # holds an unwritten slot
+IDFW_RAW16 = FIXTURES / "UM13981_20220207084555.IDFW"     # segment 0 is MODE_RAW16
+
+
+def test_histogram_decodes_past_250_intervals():
+    """The segment validator must not require a zero counter high byte.
+
+    The interval counter is a uint16 cumulative index.  Requiring its high
+    byte to be zero rejected every segment past interval 255, capping each
+    histogram at 250 intervals and truncating any run longer than ~4 hours —
+    frequently discarding the part that held the peak.
+    """
+    result = read_idf_file(IDFH_LONG)
+    header, _rows = _parse_export(IDFH_LONG.with_suffix(".IDFH.csv"))
+    expected = float(header["NumberOfIntervals"])
+    assert len(result.intervals) == 719
+    assert len(result.intervals) == pytest.approx(expected, abs=1.0)
+
+
+def test_histogram_ignores_unwritten_interval_slot():
+    """A never-written interval keeps its ±full-scale seed and must be dropped.
+
+    Counting it fabricates a 10.0 in/s peak on every channel, which then wins
+    the max-over-intervals and poisons the whole file's PPV.
+    """
+    header, _rows = _parse_export(IDFH_SENTINEL.with_suffix(".IDFH.csv"))
+    result = read_idf_file(IDFH_SENTINEL)
+    for channel, attr in (
+        ("Tran", "transverse_ips"),
+        ("Vert", "vertical_ips"),
+        ("Long", "longitudinal_ips"),
+    ):
+        decoded = getattr(result.event.peaks, attr)
+        assert decoded < 1.0, f"{channel} peak {decoded} looks like the ±FS seed"
+        assert decoded == pytest.approx(
+            _header_float(header, f"{channel}PPV"), abs=2 * _GEO_LSB_IPS
+        )
+
+
+def test_waveform_raw16_segment_zero_is_decoded():
+    """Segment-0 records can be raw int16 (MODE_RAW16, 10-byte header).
+
+    That mode was absent from the dispatch, so the record fell through
+    unhandled and the channel silently lost its first 512 samples.
+    """
+    rows = _parse_export(IDFW_RAW16.with_suffix(".IDFW.csv"))[1]
+    result = read_idf_file(IDFW_RAW16)
+    for index, channel in enumerate(GEO_CHANNELS):
+        decoded = result.samples[channel]
+        assert len(decoded) == len(rows), f"{channel} lost segment 0"
+        expected = [row[index] for row in rows]
+        bad = sum(
+            1 for c, v in zip(decoded, expected)
+            if abs(geo_count_to_ips(c) - v) >= 5e-5
+        )
+        assert bad == 0, f"{channel}: {bad} samples differ from Thor's export"
+
+
+def test_body_offset_search_is_not_quadratic():
+    """The body scan must stay cheap enough for bulk ingest.
+
+    MODE_RAW16 is (0x00, 0x00), so scanning for candidate *preambles* treats
+    every run of three zero bytes as a body start and trial-decodes each one
+    (~0.5 s/file measured).  The search anchors on record headers instead.
+    """
+    import time
+
+    start = time.perf_counter()
+    for _ in range(3):
+        read_idf_file(IDFW_RAW16)
+    elapsed = (time.perf_counter() - start) / 3
+    assert elapsed < 0.15, f"body-offset search took {elapsed*1000:.0f} ms/file"