diff --git a/CHANGELOG.md b/CHANGELOG.md index 5284dc6..591a5de 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,6 +8,84 @@ All notable changes to seismo-relay are documented here. --- +## v0.27.0 — 2026-08-28 + +**Per-sample decoder verification at scale, plus the offset investigation.** +The series-3 codec is now verified sample-by-sample against **14,338** preserved +Blastware ASCII exports — 1,249 waveform and 13,089 histogram, spanning 45 units +and files back to 2018. That is 11x the ground truth the production store +carried, and it found one real codec bug (below). + +### Fixed +- **Sub-minute histograms with a partial final block decoded to nothing** + (`histogram_codec.detect_multi_interval_stride`). The stride search confirmed + itself on a third block header whenever the body was long enough to hold one — + but a body can exceed two strides and still contain only two real blocks, because + a *partial* final block leaves trailing padding. BE18193 `T193L0XM.CI0H` (51 + intervals at 2 s = one full 30-interval block plus a 21-interval remainder, in a + 2787-byte body) therefore had its correct stride of 612 discarded and produced an + empty decode. A missing third header now means end-of-stream rather than + disqualification; the block-counter check, which is what actually prevents the + false positives that once mis-dispatched 9,082 files, is unchanged. + + Found by decoding the full DL2 archive against its preserved Blastware ASCII + exports. Across **63,535 unique** histogram binaries the fix recovers **4 files** — + `K440HJCN.3C0H` and `K557IF1U.8K0H` (stride 252), `T191HVNP.0S0H` (92) and + `T193L0XM.CI0H` (612) — with **zero** files regressed. Verification over all + 14,340 archive pairs goes 14,337 → 14,338 exact, the only remainder being two + series-4 IDF files that belong to a different codec. + + (The DL2 export keeps a byte-identical `Sent/` mirror of its root, so a naive + walk double-counts every binary — 127,035 paths are 63,535 distinct files. The + ASCII exports are *not* mirrored, so the 14,340 pair count is already distinct.) + + ⚠ Prod stores hold `.h5` files generated before this fix. Those 4 events stay + empty until `backfill_sidecars.py` is re-run — not worth a two-hour prod backfill + on its own; fold it into the next one. + +- **Histogram/waveform twin matching is now interval-based** (`find_twins`). A real + trigger is recorded twice — as a triggered waveform (stamped at the trigger instant) + and inside the scheduled histogram whose interval contains it (stamped at the 7am/7pm + interval start) — so the two twins can be **hours apart**. The old ±5-minute window + silently missed them, which broke review propagation (flagging one twin didn't flag its + twin). Twins are now matched by same serial + identical `peak_vector_sum` + opposite + record type + the waveform falling within the histogram's interval (bounded by the next + same-serial histogram). `window_seconds` is retained but ignored. Fixes terra-view #102 + sub-task 2. + +- **`/health` reported a hard-coded `0.1.0`** instead of the real service version. + `sfm/server.py` now derives its version from `minimateplus.event_file_io.TOOL_VERSION`, + making that constant the single source of truth for the service version and the + sidecar stamp alike — one place to bump at release. + +- **`CLAUDE.md` had 793 NUL bytes appended** after its last line, which made `grep` + treat the file as binary and silently skip it. Present since at least v0.21.0. + Stripped. + +### Added +- **`docs/offset_investigation.md`** — a dated journal of the "offset" hardware + fault: base rate, detector design, per-unit case files, ruled-out hypotheses + (each kept with the evidence that killed it), and Instantel's own autozero + procedure with its 2027–2069 acceptance window. +- **`scratch/verify_against_ascii.py`** — decodes a corpus of BW binaries and + diffs every sample against the paired `_ASCII.TXT`. Includes a saturation + carve-out: BW clamps clipped events to the range maximum and writes `OORANGE`, + while the decoder faithfully reports counts past nominal full scale. +- **`scratch/offset_scan3.py`** — offset detector. Measures the resting floor in + the *pre-trigger* window (definitionally quiet) and requires it to hold across + pre / middle / end. Result: **5 of 45 units (11%)**, stable across a 2x + threshold range. Supersedes `offset_scan.py` and `offset_scan2.py`, both kept + as the reasoning trail. + +### Verified +- **19,244 healthy channel-events sit at a pre-trigger floor of exactly 0.000 + (62.7%), 94.5% within ±1 quantisation unit, median +0.0000.** No systematic + zero-point bias in the decoder — an independent confirmation of the + 32000-count geo full scale, arrived at from a different direction than the + ASCII sample comparisons. + +--- + ## v0.26.0 — 2026-08-27 **Series-3 decode correctness.** Two body-model rewrites, a systematic diff --git a/CLAUDE.md b/CLAUDE.md index f2c7e75..2e3ae32 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -2,21 +2,24 @@ Ground-up Python replacement for **Blastware**, Instantel's Windows-only software for managing MiniMate Plus seismographs. Connects over direct RS-232 or cellular modem -(Sierra Wireless RV50 / RV55). Current version: **v0.26.0**. +(Sierra Wireless RV50 / RV55). Current version: **v0.27.0**. --- -## Where things stand (updated 2026-08-27) +## Where things stand (updated 2026-08-28) Read this first when picking the project back up. -- **Series-3 decode is correct and verified.** All 11,603 series-3 binaries in - the prod snapshot pass every check (channel lengths, peaks vs the device's - own reported PPV, nothing above full scale, length vs declared record time). - Ground truth: 1,211/1,211 histograms exact per-interval and 75/75 waveform - sample counts exact against preserved Blastware ASCII exports. - ⚠ That is per-sample proof on 11% of files and peak-only consistency on the - other 89% — see `docs/instantel_protocol_reference.md` §7.6.1. +- **Series-3 decode is verified per-sample at scale (v0.27.0).** The full DL2 + archive decodes **14,338 / 14,338** paired files exactly against their + preserved Blastware ASCII exports — 1,249 waveform + 13,089 histogram, 45 + units, files back to 2018. That is 11x the ground truth the prod store + carried, and it supersedes the old "per-sample on 11%, peak-only on 89%" + caveat. Harness: `scratch/verify_against_ascii.py` (note its saturation + carve-out — BW clamps clipped events, the decoder reports true counts). + Independent corroboration of the 32000-count scale: 19,244 healthy + channel-events sit at a pre-trigger floor of exactly 0.000 (62.7%), 94.5% + within ±1 quantisation unit, median +0.0000 — no zero-point bias. - **Series-4 (Thor / Micromate) is NOT verified.** UM-series sits at ~48% against device peaks with a ~1.7% systematic bias and a near-zero tail. Thor IDFW is pinned to `decode_waveform_legacy` deliberately. @@ -24,11 +27,23 @@ Read this first when picking the project back up. (= 10.0/1.25) units discrepancy; `scripts/backfill_sidecars.py --force` also inserts DB rows for store files that have none (one-time per store) and the dry-run does not report that count. -- **After any codec change, regenerate the store** — `backfill_sidecars.py - --force` then `backfill_event_shape.py`, DB backup first. Stored `.h5` files - do not update themselves. -- **Parked:** the "offset" hardware-fault investigation (Appendix E of the - protocol reference) pending the multi-year BW archive. +- **After any codec change, regenerate the store** — `backfill_sidecars.py` + then `backfill_event_shape.py`, DB backup first. Stored `.h5` files do not + update themselves. No `--force` needed as long as `TOOL_VERSION` was bumped + (it gates regeneration). ⚠ On the office NAS this takes **~2 hours** + (~1.5 files/sec vs 85/sec on the dev box — gzip-4 in `sfm/event_hdf5.py` + against a Synology CPU). Budget it up front. + **v0.27.0 owes prod a backfill:** the partial-final-block fix recovers 4 + histograms that are still empty in the store. +- **The "offset" hardware fault has its own journal** -- + `docs/offset_investigation.md`. **5 of 45 units (11%)**, and the fault is + **persistent** — it stays until the geophone is serviced. Detect it with + `scratch/offset_scan3.py`: the resting floor in the **pre-trigger** window, + required to hold across pre/middle/end. Never score only the dominant-peak + axis and never use the mean — both produce false recoveries (see the + retraction banner in the journal). Instantel's autozero procedure and its + 2027-2069 acceptance window are recorded there too. Best open lead is + `SUB 0x0E` (unimplemented), which may carry those very numbers. When new information about the protocol is discovered, please update the instantel_protocol_reference.md with the findings in addition to this document @@ -1803,4 +1818,4 @@ body) because writing a dial string may require DLE escaping for embedded contro To parse BW TX captures: use `bridges/captures/` scripts or adapt the `find_write_frames()` pattern in `/tmp/analyze_write_payload.py` — it correctly handles `0x10 0x03` DLE-escaped ETX bytes -inside write frame data (the naive parser terminates early at the escaped `0x03`). \ No newline at end of file +inside write frame data (the naive parser terminates early at the escaped `0x03`). diff --git a/README.md b/README.md index 02f7644..4cb9676 100644 --- a/README.md +++ b/README.md @@ -1,4 +1,4 @@ -# seismo-relay `v0.26.0` +# seismo-relay `v0.27.0` A ground-up replacement for **Blastware** — Instantel's aging Windows-only software for managing seismographs. Supports both the **MiniMate Plus diff --git a/docs/instantel_protocol_reference.md b/docs/instantel_protocol_reference.md index 4fecbad..a75393a 100644 --- a/docs/instantel_protocol_reference.md +++ b/docs/instantel_protocol_reference.md @@ -3524,6 +3524,11 @@ of them was, for a while. ### E.1 The "offset" fault +> **The full investigation now lives in `docs/offset_investigation.md`** -- +> base rate, detector definition, per-unit case files, ruled-out hypotheses, +> and Instantel's own autozero procedure with its 2027-2069 acceptance +> window. This appendix is kept as the protocol-side summary. + **Symptom.** One geophone channel's baseline steps away from zero and stays there. The trace still carries the real AC signal, but it rides on a DC pedestal of a few tenths of an in/s. Operators call this an diff --git a/docs/offset_investigation.md b/docs/offset_investigation.md new file mode 100644 index 0000000..20da227 --- /dev/null +++ b/docs/offset_investigation.md @@ -0,0 +1,514 @@ +# The "offset" fault — investigation journal + +> ## ⚠ CORRECTED 2026-08-28 (same day) — the v1 detector was wrong +> +> Brian pushed back on the finding that offsets "come and go": in the field, +> once a unit develops one it stays broken until the geophone is replaced. +> He was right, and the challenge exposed **two real flaws** in the v1 detector: +> +> 1. **It scored only the axis with the largest peak.** A real event on one axis +> hid a persistent pedestal on another. BE12599 on 2026-08-21 read "clean" +> solely because Long had a 1.065 in/s event — Tran was sitting at +> **+0.4732 in/s** at that moment and was never examined. +> 2. **It used the MEAN**, which a real transient perturbs. The **median** is the +> resting baseline — most samples sit at it, so a blast does not move it. +> Same event, Long channel: mean **+0.0783** vs median **-0.0050**. +> +> Both flaws manufactured false recoveries. The corrected detector +> (`scratch/offset_scan2.py`, per-channel median) shows the pedestal is +> **persistent**, exactly as the field experience says. See §2b and §3b. +> +> **Then Brian proposed a better detector still** — measure the floor during +> the *pre-trigger* window, and require it to hold across pre/middle/end. +> That is now the detector of record (§2c). Final answer: **5 of 45 units +> (11%)**, stable across a 2x threshold range. +> +> Sections below that were written against v1 are marked; v1 numbers are kept +> for the reasoning trail, not as current fact. + +A running record of the **offset** hardware fault on Instantel Series III +seismographs: a geophone channel whose trace sits displaced from zero rather +than centred on it. + +This is a *journal*, not a spec. Findings are dated, dead ends are kept with +the reason they died, and every number says where it came from. When something +here is superseded, strike it and say why rather than deleting it — the point +is that a future session can tell what was actually established from what was +merely believed at the time. + +Companion material: +- `scratch/offset_scan.py` — the detector +- `scratch/verify_against_ascii.py` — decoder verification harness +- `docs/instantel_protocol_reference.md` — wire protocol, incl. the + unimplemented `SUB 0x0E` this investigation now wants + +--- + +## TL;DR (current state, 2026-08-28) + +- **It is real device data, not a decode bug.** Settled early and confirmed + against Blastware's own ASCII exports. +- **Base rate: 5–6 of 45 units (11–13%)** across the full DL2 archive, + 2018–2026. This *confirms* the earlier 2-of-21 (9.5%) estimate from the much + smaller Terra-View DB — survivorship bias from deleted events had **not** + concealed a wave of cases. +- **The fault is bimodal, not a drift continuum.** A unit is either clean or + grossly off. Loosening the amplitude threshold 11× adds no new units. +- **The unit's own sensor check cannot see it.** 102 offset events, zero + sensor-check failures. Do not try to use it as a screen. +- **Cause is still unsettled.** Instantel's autozero fixes the minority of + cases; the rest are hardware. We cannot yet tell which is which remotely. +- **Best open lead:** `SUB 0x0E` (channel sensor data, 8 channels × 10 bytes, + unimplemented) may carry the very numbers Instantel says to check against + **2027–2069**. Untested. + +--- + +## 1. What the fault looks like + +A healthy geophone trace is centred on zero. An offset channel is parked away +from zero, so the channel **mean approaches its own peak**. In Blastware the +signature is "parallel lines above or below the zero line" (Instantel's own +wording). + +Consequences observed in the field: +- The unit can **self-trigger on its own offset** when the displacement exceeds + the geo trigger level, producing streams of junk events with no ground + motion. Instantel has a separate FAQ for this symptom (13-0-22, *"Unit + triggers continuously without activity"*). +- Recorded PPV for that channel is meaningless while the fault persists. + +--- + +## 2. The detector + +Implemented in `scratch/offset_scan.py`. Operates on raw BW binaries only — no +DB, no sidecars. + +``` +for each series-3 waveform binary: + decode -> per-channel ADC counts + dominant axis = channel with the largest |peak| + flag when |mean| / peak > 0.70 + and |mean| >= 0.90 x the unit's geo trigger level +episodes = per-serial runs of flagged events, split on a >12 h gap +``` + +Why each term: + +| term | purpose | +|---|---| +| `\|mean\|/peak > 0.7` | the discriminator. A DC-parked trace has mean ≈ peak. | +| `\|mean\| >= 0.9 × trigger` | amplitude floor — suppresses quiet traces where mean and peak are both tiny and the ratio is meaningless. | +| dominant axis only | the fault is per-channel; scoring all three dilutes it. | +| 12 h episode gap | separates deployments/visits rather than counting events. | + +Trigger level comes from a paired `_ASCII.TXT` when one exists, else the +per-serial median learned from that unit's ASCII files, else 0.2 in/s. + +**Known limitation.** Event traces contain real ground motion, so this can only +see offsets large enough to *dominate* the trace. A mild offset on a real blast +is invisible. Instantel's A/D-mode check (§5) is the only thing that sees the +mild end. Our base rate is therefore a **gross-offset** rate. + +### 2b. Detector v2 — per-channel median (CURRENT) + +`scratch/offset_scan2.py`. Supersedes the above. + +``` +for each series-3 waveform binary: + for each geo channel independently: + pedestal = median(samples) # resting baseline, robust to blasts + flag the CHANNEL when |pedestal| >= 0.025 in/s (5 A/D counts) +a unit has a real fault when a channel is flagged on >=3 CONSECUTIVE events +``` + +Why median: a DC pedestal shifts every sample, so it moves the median. A real +event moves only a minority of samples, so it does not. This removes the need +for the `m/p` ratio guard entirely — that guard existed only to compensate for +using the mean. + +Why per-channel: the fault is on one geophone axis. Scoring only the dominant +axis means any event with motion elsewhere hides it. + +Why "3 consecutive": the 0.025 in/s floor is only ~2x a healthy channel's +resting median (observed 0.010-0.015), so isolated flags are noise. Persistence +is the discriminator — and it is what the field experience predicts. + +--- + +## 3b. Archive results, corrected (v2) + +| | v1 (dominant axis, mean) | **v2 (per-channel median)** | +|---|---|---| +| units with any flagged event | 6 of 45 | 19 of 45 | +| **units with a sustained pedestal (>=3 consecutive)** | — | **8 of 45 (18%)** | +| runs of >=3 consecutive | — | 29 | +| runs of 1-2 events (noise) | — | 69 | + +Units with a sustained pedestal: **BE9558, BE10895, BE11007, BE11529, BE12599, +BE13117, BE18003, BE18438**. BE10895 and BE18003 were invisible to v1. + +**The affected channel is most often Vert**, which v1 got wrong — it named +whichever axis had the largest peak. BE13117 and BE18438 are both Vert faults. + +Longest / clearest runs: + +| unit | ch | span | events | median in/s | +|---|---|---|---|---| +| BE13117 | Vert | 2023-05-03 → 05-04 | 194 | 0.035 → **1.915** | +| BE18438 | Vert | 2026-02-25 → 02-26 | 75 | 0.180 → 0.370 | +| BE9558 | Vert | 2020-02-11 (6 h) | 33 | 0.065 → 0.090 | +| BE12599 | Tran | 2026-08-14 → 08-23 | 8 | 0.030 → **0.565** | +| BE18003 | Vert | 2021-03-17 → 06-11 | 3 | 0.040 → 0.060 | + +BE12599 began **2026-08-14**, not 08-17 as v1 reported, and was still faulting +at the last event in the archive. + +### 2c. Detector v3 — PRE-TRIGGER floor + constant-floor test (CURRENT) + +`scratch/offset_scan3.py`. Brian's method, and better than v2 for a reason +worth naming: **the pre-trigger window is definitionally quiet** — it is the +buffer captured before the trigger fired — whereas a whole-record median is +merely *robust* to the event. `pretrig_samples` comes from the STRT record. + +``` +per channel: + pre = median of the first pretrig_samples samples + mid = median of the middle third + end = median of the final third + spread = max(pre,mid,end) - min(pre,mid,end) + + offset when |pre| >= floor AND spread <= 0.02 in/s +real fault when a channel is flagged on >=3 CONSECUTIVE events +``` + +A DC offset is a **constant floor** — present before the trigger, during, and +after. The spread test rejects transients (settling, handling, a long event +tail) that move one segment relative to the others, which is what v2's +whole-record median could not do. + +**The empirical noise floor justifies the threshold.** Across 19,244 +non-flagged channel-events the pre-trigger floor distributes as: + +| floor | share | +|---|---| +| −1 unit (−0.005) | 18.4% | +| **0.000** | **62.7%** | +| +1 unit (+0.005) | 13.4% | + +**94.5% within ±1 quantisation unit; median exactly +0.0000, mean −0.0008.** +So there is **no systematic zero-point bias in the decoder** — an independent +confirmation of the 32000-count scale. A healthy channel really does read +0.000, and "any constant floor that is not 0.000" is the right signal, with +±1 unit of slack for quantisation. + +**The result is threshold-insensitive**, which is what distinguishes a real +signal from a tuned one: + +| floor | units flagged | sustained units | +|---|---|---| +| 2 units (0.010) | 34 | 15 ← into the noise | +| 3 units (0.015) | 26 | 8 | +| **4 units (0.020)** | 17 | **5** | +| **5 units (0.025)** — Instantel's | 12 | **5** | +| **8 units (0.040)** | 8 | **5** | + +### FINAL RESULT: 5 of 45 units (11%) + +**BE9558, BE11529, BE12599, BE13117, BE18438.** + +Unchanged across a 2x threshold range. BE11007 and BE10895 drop out — the +spread test identifies them as transients, not pedestals. + +The 11% headline happens to match v1's, but the reasoning and the unit list +differ: v1 included BE11007 and named the wrong *channel* on most units. + +--- + +## 3. Archive results (2026-08-28) + +Source: DL2 event export, 6,577 **unique** series-3 waveforms, 45 units. +See [`dl2-archive`](#8-data-and-tooling) for the `Sent/` mirror trap. + +**283 suspect events, 15 episodes, 6 of 45 units (13.3%).** +Excluding BE11007 (§4, likely not an offset at all): **5 of 45 = 11.1%**. + +### Threshold sensitivity — the bimodality result + +Re-scoring the same corpus at a range of amplitude floors, with two +ratio cut-offs (1 A/D count = 0.005 in/s, see §5): + +| \|offset\| floor | m/p > 0.7 | m/p > 0.9 | +|---|---|---| +| 5 cts (0.025 in/s) — *Instantel's own* | 333 ev / 6 units | 279 ev / **5 units** | +| 10 cts (0.050) | 294 / 6 | 274 / 5 | +| 20 cts (0.100) | 250 / 5 | 244 / 4 | +| 40 cts (0.200) | 209 / 5 | 203 / 4 | +| 80 cts (0.400) | 152 / 4 | 148 / 2 | +| 160 cts (0.800) | 144 / 2 | 141 / 1 | + +Relaxing the floor by 11× (0.27 → 0.025 in/s) adds ~14% more events and **no +new units**. There is no population of mild offsets hiding below our threshold +*in event data*. Either a unit is clean or it is grossly off. + +--- + +## 4. Per-unit case files + +Ordered by severity. `m/p` medians are on the offending channel. + +### BE13117 — one violent day, never again +`145 / 454 events (32%)`, **1 episode**, 2023-05-04, 6.8 h. +Offset climbed **0.393 → 1.875 in/s within the episode**. `m/p` median +**0.996** — the trace is almost pure DC. No recurrence in the rest of its 454 +events. No ASCII files in the archive, so no calibration history. + +### BE18438 — recurring, months apart +`87 / 293 (30%)`, **2 episodes**: 2025-11-15 (1.2 h, n=12, 0.279 → 0.369) and +2026-02-25 (**28.8 h**, n=75, 0.183 → 0.366). `m/p` median 0.967. +Clean across all 196 events preceding its 2025-08-12 calibration. + +### BE9558 — six years apart +`38 / 196 (19%)`, **4 episodes**: 2020-02-11 (6.3 h, n=33, but only +0.068 → 0.086 — very mild), then 2026-04-14, 2026-04-29, 2026-05-04 +(0.28–0.45). `m/p` median 0.919. Calibrated 2026-06-26; 0/7 events flagged +after, but n=7 is far too small to call it fixed. + +### BE12599 — the live case ⚠ +`6 / 77 (8%)`, **6 single-event episodes, one per day at exactly 05:00**, +2026-08-17 → 2026-08-23. Offset rose 0.383 → 0.565 then fell back to 0.345. +`m/p` ≈ 0.965, geo trigger 0.3 in/s — **the offset exceeds the trigger level, +so the unit is triggering on its own fault**. Last calibrated 2025-08-12. + +This is the most recent and the most useful: a currently-faulting unit is the +natural experiment for the re-zero-vs-repair question (§7). + +### BE11529 — marginal +`4 / 99 (4%)`, 1 episode 2025-07-08, 0.4 h, offsets only 0.051 → 0.058 in/s. +`m/p` median 0.959, so DC-dominated, but the magnitude is near the noise of +this method. Treat as unconfirmed. + +### BE11007 — probably NOT an offset +`3 / 70 (4%)`, 1 episode 2022-01-17, offsets 5.500 → 6.904 in/s — by far the +largest. But `m/p` is only **0.719–0.738** against ≥0.9 for every other unit, +and the peaks are 7.6–9.4 in/s on a 10 in/s range. That reads as a **large +real blast with asymmetric ground motion**, not a parked trace. Excluded from +the headline base rate. + +--- + +## 5. Instantel's own procedure and thresholds + +From two Instantel technical-support FAQs supplied 2026-08-28 +(answers **13-0-21** *"How to determine offsets"* and **12-0-10** *"Removing +offsets on an Instantel Series III monitor"*; created 2008/2007, last updated +2009-03-06). + +### Identifying (13-0-21) + +1. Create or use an event with the **manual minimum trigger** set for the + connected geophone and microphone — i.e. an event that recorded no real data. +2. Save it and open in Blastware. +3. An offset shows as **parallel lines above or below the zero line**. +4. Put the unit in **A/D mode** — on Series III, press and hold `OPTION`, then + press `START MONITOR`. +5. **Display counts higher than 5**, with no vibration or overpressure present, + indicate an offset. + +### Removing — the autozero (12-0-10) + +1. Be in a **quiet area with low vibration**. +2. Power on the Blastmate III / Minimate Plus. +3. Connect the geophone and microphone — **LINEAR mic only**. + ⚠ *Do not connect an "A" weight microphone, regardless of what the monitor + displays.* +4. Press `Test`. +5. Wait for the **Sensor Check** results to appear. +6. Press `OPTION` and `START MONITOR` **simultaneously**. +7. `Performing Autozero` appears; press `Enter`. +8. Confirm the sensors are properly connected; press `Enter`. +9. Wait for the autozero to complete. +10. Press `Enter` twice → Main Menu, *Ready To Monitor*, offset corrected. + +### The go/no-go number — 2027 to 2069 + +> When you perform an Autozero on any Series III unit, the lists of numbers in +> the **X1 and X8 gains should all be between 2027 and 2069**. If not, repeat +> the Autozero. **If the numbers are extremely out of the specified range, then +> the unit should be sent in for repair.** +> +> If this process does not remove the offset problem, return the unit **and +> sensors** to Instantel for repair. + +This is the documented explanation for the field experience (Brian's dad, +2026-08-28) that **a re-zero works maybe 10% of the time** — the autozero only +recovers units whose zero reference is still near-correct. + +### Scale derivation (inference, well-supported — not proven) + +2048 is 12-bit midscale. Our codec's geo full scale is 32000 internal counts = +10 in/s, with 1 decoder unit = 16 counts = exactly 0.005 in/s +(`geo-full-scale-is-32000-counts`). ±2000 A/D counts about 2048 therefore maps +to ±10 in/s at **0.005 in/s per A/D count**. That makes: + +- Instantel's ">5 counts" threshold ≈ **0.025 in/s** +- the 2027–2069 window = **±21 counts = ±0.105 in/s** of tolerated zero error + +Consistent and mutually corroborating, but we have not confirmed the A/D-count +scale directly from a device reading. + +--- + +## 6. Ruled out — keep these dead + +### Condensation / humidity — DEAD (2026-08-25) +Proposed, then killed by its own controls: BE18438 stayed flat across a 10-hour +overnight gap, and only 2 of 21 units showed the fault while 19 sat in the same +weather. The apparent "diurnal cycle" was an artifact of binning by hour-of-day +across two days. See `waveform-dc-offset-is-real-device-data`. + +### Clipping as a false-positive source — RULED OUT (2026-08-28) +A rail-hitting trace would fake an offset (mean → peak). It isn't happening: +median suspect peak is only **10% of full scale**, p90 is 18.6%. Only BE11007's +3 events exceed 50% FS, and none reach 98%. + +### The sensor check as a predictor — DOES NOT WORK (2026-08-28) +Tested on 102 offset events across 4 units: + +| unit | state | n | failed | median ratio | median freq | +|---|---|---|---|---|---| +| BE11529 | offset | 4 | **0** | 3.90 | 7.6 | +| BE11529 | clean | 14 | 0 | 3.80 | 7.5 | +| BE12599 | offset | 6 | **0** | 4.00 | 7.4 | +| BE12599 | clean | 13 | 0 | 4.00 | 7.6 | +| BE18438 | offset | 87 | **0** | 3.70 | 7.6 | +| BE18438 | clean | 25 | 0 | 3.80 | 7.5 | +| BE9558 | offset | 5 | **0** | 3.90 | 7.8 | +| BE9558 | clean | 44 | 0 | 3.80 | 7.5 | + +Zero failures on either side and indistinguishable ratios/frequencies. The +swing test measures geophone frequency response and damping — it never examines +DC zero. **A grossly offset unit passes its own self-check.** This is why the +fault goes unnoticed until somebody looks at waveforms. + +### "Offsets are transient / come and go on their own" — RETRACTED 2026-08-28 +v1 reported episodes lasting hours that ended spontaneously. **This was an +artifact of the v1 detector** (see the banner at the top). With the per-channel +median, the pedestal persists. Every clear case reads clean again only after a +multi-day-to-multi-month gap consistent with service: BE13117 6 days, BE18438 +24 days, BE9558 63 days **with a confirmed Instantel calibration inside the +gap**. BE12599 never reads clean — it is still faulting at the end of the +archive. This matches the operational experience: once a unit develops an +offset it stays broken until the geophone is replaced. + +### "Offsets develop N months after calibration" — CONFOUNDED, NOT A FINDING +Tempting, and it looked strong: + +| unit | suspect before latest cal | after | +|---|---|---| +| BE18438 | 0 / 196 | 87 / 97 | +| BE12599 | 0 / 62 | 6 / 15 | +| BE11529 | 0 / 82 | 4 / 17 | +| BE9558 | 38 / 189 | 0 / 7 | + +But bucketing suspects by months-since-calibration gives **one unit per bucket**: +`0–3mo={BE11529}`, `3–6 & 6–9mo={BE18438}`, `9–12mo={BE9558}`, +`12–15mo={BE12599}`. The apparent "51% failure rate at 6–9 months" is entirely +BE18438's single February 2026 episode. Five units with roughly one episode +each cannot support a population trend. **Do not re-derive this.** + +Also note: all affected units are calibrated on a **~12–13 month cadence**, so +"sent to Instantel" is the routine annual schedule, not evidence of a +fault-driven return. + +--- + +## 7. Open questions + +### Q1 — Is it a latched bad zero or analog degradation? +The question that decides everything. A latched zero is correctable (possibly +over the wire); degradation means a repair. Instantel's 2027–2069 rule implies +*both* populations exist, with the split roughly 10/90 in the field. + +**BE12599 is the natural experiment** — faulting as of 2026-08-23. Read its +values, run the autozero, read them again. + +### Q2 — Can we read the autozero numbers over the wire? (best lead) +Instantel says to check *"the lists of numbers in the **X1 and X8 gains**"* — +4 sensors × 2 gains = **8 channels**. The protocol reference already documents +an unimplemented command with exactly that shape: + +``` +SUB 0x0E -> RSP 0xF1 "channel sensor data" + 2-step read; channel selector in params[6:8] = 0x0000..0x0007 + data length 0x0A (10 bytes) per channel +``` + +Blastware's *Unit Channel Test* sequence: +`POLL×N → 0x15 → 0x01 → 0x08 → 0x01 → 0x0E×8 → 0x98×2 → 0x0E×8` +— note the **second `0x0E` pass carries live ADC readings**. + +**Hypothesis (untested):** `0x0E` returns the numbers Instantel wants compared +against 2027–2069. If true, SFM could diagnose an offset remotely *and* predict +whether a re-zero will succeed — converting a 10%/90% shipping gamble into a +decision made before packing a box. + +**How to test.** `bridges/ach_mitm.py` is a generic TCP proxy: + +```bash +python bridges/ach_mitm.py --bw-host --bw-port 9034 --listen-port 9999 +``` + +Point Blastware at the proxy and run **Unit Channel Test**. +⚠ In this topology the output filenames are reversed — the tool labels the +*connecting* side "unit", so `raw_s3_*.bin` holds Blastware's bytes and +`raw_bw_*.bin` the unit's. + +Capture priority: (1) BE12599 while faulting, (2) a known-good unit as control, +(3) before/after an autozero on the same unit. Eight 10-byte payloads with an +expected value near 2048 is a very constrained puzzle. + +### Q3 — What is the mild-offset rate? +Unmeasurable from event files (§2). Only the A/D-mode check sees it. Would +need a fleet sweep in A/D mode, or Q2 to succeed. + +### Q4 — Does an offset recur on the same unit after service? +BE9558 shows episodes in 2020 and 2026; BE18438 twice in four months. Suggestive +of recurrence, but service records aren't in the data — only calibration dates. + +--- + +## 8. Data and tooling + +| what | where | +|---|---| +| detector | `scratch/offset_scan.py` | +| current results | `/home/serversdown/dl2-archive/offset_archive.csv` | +| earlier candidate list (Terra-View DB, 274 events) | `scratch/offset_candidates.csv` | +| archive working copy | `/home/serversdown/dl2-archive/files/` | +| archive source | NAS `DeathStar` 10.0.0.2, `/volume1/Uploads/TMI/DL2-Event-backup-8-25-26/Event/autocall home/` | + +⚠ **The DL2 export keeps a byte-identical `Sent/` mirror of its root.** 13,077 +waveform paths are 6,577 distinct files. Always dedupe by basename — this +doubled two reported figures before it was caught. + +--- + +## 9. Chronology + +| date | event | +|---|---| +| 2026-08-25 | Reported as a *waveform decode bug* — traces with a DC offset. Investigation shows the offset is **real device data**; the decoder is correct. | +| 2026-08-25 | Brian relays his dad's description: a known hardware fault called an "offset"; usually sent to Instantel. | +| 2026-08-25 | First detection pass over the Terra-View DB: **2 of 21 units**, 274 events, 11 months. Flagged as vulnerable to survivorship bias — flooded events were routinely deleted. | +| 2026-08-25 | Condensation hypothesis proposed, then **killed by its own controls**. | +| 2026-08-25 | Parked pending the multi-year archive. | +| 2026-08-28 | DL2 archive pulled (33 GB, 546k files; 6.6 GB working set). | +| 2026-08-28 | Archive scan: **6 of 45 units**, 283 events, 15 episodes. Prior base rate **confirmed**, not overturned. | +| 2026-08-28 | Clipping ruled out; `m/p` established as the discriminator; BE11007 reclassified as probably a real blast. | +| 2026-08-28 | Calibration-timing correlation attempted and **rejected as confounded**. | +| 2026-08-28 | Instantel FAQs supplied: autozero procedure, the **2027–2069** window, the **>5 counts** threshold. Explains the ~10% re-zero success rate. | +| 2026-08-28 | Bimodality established; sensor check proven **blind** to offsets; `SUB 0x0E` identified as the best open lead. | +| 2026-08-28 | **v1 detector retracted.** Brian challenged the "come and go" finding against field experience. Two flaws found: dominant-axis-only scoring and mean-instead-of-median. Corrected detector shows persistent pedestals on **8 of 45 units**, and the gaps are service windows. | +| 2026-08-28 | **Detector v3 (Brian's method):** pre-trigger floor + pre/mid/end consistency. Healthy channels proven to sit at 0.000 +/-1 unit (94.5%), confirming no decoder zero-point bias. Final: **5 of 45 units (11%)**, threshold-insensitive. | diff --git a/minimateplus/event_file_io.py b/minimateplus/event_file_io.py index b0393d8..09a5bf1 100644 --- a/minimateplus/event_file_io.py +++ b/minimateplus/event_file_io.py @@ -50,7 +50,7 @@ SIDECAR_KIND = "sfm.event" # bumped without a `pip install` re-run — leading to confusing stale # version stamps in sidecars. Bump this constant and CHANGELOG.md # together at release time. -TOOL_VERSION = "0.26.0" +TOOL_VERSION = "0.27.0" try: # Best-effort: prefer the installed metadata when it's NEWER than the diff --git a/minimateplus/histogram_codec.py b/minimateplus/histogram_codec.py index d60f853..a24d281 100644 --- a/minimateplus/histogram_codec.py +++ b/minimateplus/histogram_codec.py @@ -413,10 +413,18 @@ def detect_multi_interval_stride(body: bytes) -> Optional[int]: if (_ctr(stride) - _ctr(0)) & 0xFFFF != 1: continue - # confirm on a third block when the body is long enough - if 2 * stride + _MULTI_HEADER_LEN <= len(body): - if not _is_multi_header(body, 2 * stride): - continue + # Confirm on a third block WHEN ONE IS ACTUALLY PRESENT. A body can + # be longer than two strides and still hold only two real blocks: a + # final *partial* block leaves trailing padding. E.g. 51 intervals at + # 2 s = one full 30-interval block + a 21-interval remainder, in a + # 2787-byte body — long enough to demand a third header at 1224 that + # does not exist. Requiring it unconditionally threw away the correct + # stride and the file decoded to nothing (BE18193 T193L0XM.CI0H). + # The block-counter check above is the decisive anti-false-positive + # test; this one is corroboration, so a missing third header means + # end-of-stream, not disqualification. + if (2 * stride + _MULTI_HEADER_LEN <= len(body) + and _is_multi_header(body, 2 * stride)): if (_ctr(2 * stride) - _ctr(stride)) & 0xFFFF != 1: continue return stride diff --git a/pyproject.toml b/pyproject.toml index 24a4e1b..6153c0f 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -4,7 +4,7 @@ build-backend = "setuptools.build_meta" [project] name = "seismo-relay" -version = "0.26.0" +version = "0.27.0" description = "Python client and REST server for MiniMate Plus seismographs" requires-python = ">=3.10" dependencies = [ diff --git a/scratch/offset_scan.py b/scratch/offset_scan.py new file mode 100644 index 0000000..1b4e377 --- /dev/null +++ b/scratch/offset_scan.py @@ -0,0 +1,171 @@ +#!/usr/bin/env python3 +"""Scan series-3 waveform binaries for the 'offset' hardware fault. + +A healthy geophone trace is centred on zero. An offset unit sits displaced, +so the channel mean approaches its own peak. Detector (unchanged from the +2026-08-25 run, see memory note `offset-archive-analysis-backlog`): + + dominant-axis |mean| / peak > 0.7 + AND |mean| >= 0.9 * the unit's geo trigger level + +Trigger level is read from a paired _ASCII.TXT where one exists, otherwise +from a per-serial median learned across that unit's ASCII files, otherwise +--default-trigger. + +Serial is decoded from the BW filename: prefix letter encodes thousands +(chr(ord('B') + n)), next 3 digits the remainder -- T193 -> BE18193. + +Usage: + python scratch/offset_scan.py --dir [--jobs N] --out offsets.csv +""" +from __future__ import annotations + +import argparse, csv, json, re, sys +from collections import defaultdict +from concurrent.futures import ProcessPoolExecutor, as_completed +from pathlib import Path + +sys.path.insert(0, str(Path(__file__).resolve().parent.parent)) + +from minimateplus.event_file_io import read_blastware_file +from minimateplus.bw_ascii_report import parse_report + +GEO = ("Tran", "Vert", "Long") +_GEO_FS_COUNTS = 32000.0 +_WAVE_RE = re.compile(r"\.[A-Za-z0-9]{2}0[Ww]$") +_STEM_RE = re.compile(r"^([B-Z])(\d{3})") + +MEAN_OVER_PEAK_MIN = 0.7 +TRIGGER_FRACTION = 0.9 + + +def serial_from_name(name: str): + m = _STEM_RE.match(name) + if not m: + return None + letter, digits = m.group(1), m.group(2) + return f"BE{(ord(letter) - ord('B')) * 1000 + int(digits)}" + + +def counts_to_ips(c, gr): + return c * (gr or 10.0) / _GEO_FS_COUNTS + + +def scan_one(path_str: str, default_trigger: float) -> dict | None: + p = Path(path_str) + try: + gr, trig = 10.0, None + ap = p.with_name(p.name.replace(".", "_", 1) + "_ASCII.TXT") \ + if False else p.parent / (p.stem + "_" + p.suffix.lstrip(".") + "_ASCII.TXT") + if ap.exists(): + rep = parse_report(ap.read_text(errors="replace")) + gr = rep.geo_range_ips or 10.0 + trig = rep.geo_trigger_level_ips + ev = read_blastware_file(p) + s = ev.raw_samples or {} + if not all(s.get(c) for c in GEO): + return None + best = None + for ch in GEO: + arr = s[ch] + n = len(arr) + if n == 0: + continue + mean = sum(arr) / n + peak = max(abs(v) for v in arr) + if peak == 0: + continue + ratio = abs(mean) / peak + if best is None or peak > best["peak_counts"]: + best = {"channel": ch, "mean_counts": mean, + "peak_counts": peak, "ratio": ratio} + if best is None: + return None + ts = ev.timestamp + return { + "serial": serial_from_name(p.name) or "?", + "timestamp": (f"{ts.year:04d}-{ts.month:02d}-{ts.day:02d}T" + f"{ts.hour:02d}:{ts.minute:02d}:{ts.second:02d}") if ts else "", + "filename": p.name, + "channel": best["channel"], + "offset_ips": round(counts_to_ips(best["mean_counts"], gr), 4), + "peak_ips": round(counts_to_ips(best["peak_counts"], gr), 4), + "mean_over_peak": round(best["ratio"], 3), + "trigger_level_ips": trig if trig is not None else "", + "geo_range_ips": gr, + } + except Exception: + return None + + +def main(): + ap = argparse.ArgumentParser() + ap.add_argument("--dir", required=True) + ap.add_argument("--jobs", type=int, default=4) + ap.add_argument("--limit", type=int, default=0) + ap.add_argument("--default-trigger", type=float, default=0.2) + ap.add_argument("--out", required=True) + a = ap.parse_args() + + # The DL2 export keeps a byte-identical `Sent/` mirror of the root, so + # enumerate paths but keep only the first occurrence of each basename — + # otherwise every event is counted twice. + seen = set() + files = [] + for q in sorted(Path(a.dir).rglob("*")): + if q.is_file() and _WAVE_RE.search(q.name) and q.name not in seen: + seen.add(q.name) + files.append(q) + if a.limit: + files = files[: a.limit] + print(f"waveform binaries to scan: {len(files)}", flush=True) + + rows = [] + with ProcessPoolExecutor(max_workers=a.jobs) as ex: + futs = [ex.submit(scan_one, str(p), a.default_trigger) for p in files] + for n, f in enumerate(as_completed(futs), 1): + r = f.result() + if r: + rows.append(r) + if n % 2000 == 0: + print(f" {n}/{len(files)}", flush=True) + + # learn per-serial trigger levels from the rows that had an ASCII + by_serial = defaultdict(list) + for r in rows: + if r["trigger_level_ips"] != "": + by_serial[r["serial"]].append(float(r["trigger_level_ips"])) + med = {} + for k, v in by_serial.items(): + v.sort() + med[k] = v[len(v) // 2] + + for r in rows: + if r["trigger_level_ips"] == "": + r["trigger_level_ips"] = med.get(r["serial"], a.default_trigger) + r["suspect"] = int( + r["mean_over_peak"] > MEAN_OVER_PEAK_MIN + and abs(r["offset_ips"]) >= TRIGGER_FRACTION * float(r["trigger_level_ips"]) + ) + + cols = ["serial", "timestamp", "filename", "channel", "offset_ips", "peak_ips", + "mean_over_peak", "trigger_level_ips", "geo_range_ips", "suspect"] + with open(a.out, "w", newline="") as fh: + w = csv.DictWriter(fh, fieldnames=cols) + w.writeheader() + w.writerows(rows) + + sus = [r for r in rows if r["suspect"]] + print(f"\nscanned {len(rows)} decodable waveforms") + print(f"suspect events: {len(sus)}") + per = defaultdict(int) + for r in sus: + per[r["serial"]] += 1 + print(f"units with >=1 suspect event: {len(per)} of {len({r['serial'] for r in rows})}") + for s, n in sorted(per.items(), key=lambda x: -x[1])[:20]: + print(f" {s:10} {n}") + print(f"\nwrote {a.out}") + + +if __name__ == "__main__": + main() diff --git a/scratch/offset_scan2.py b/scratch/offset_scan2.py new file mode 100644 index 0000000..5fe2e62 --- /dev/null +++ b/scratch/offset_scan2.py @@ -0,0 +1,108 @@ +#!/usr/bin/env python3 +"""Offset detector v2 — per-channel MEDIAN pedestal. + +Supersedes the dominant-axis / mean detector in offset_scan.py, which had two +flaws that manufactured false "recoveries": + + 1. It scored only the axis with the largest peak, so a real event on one axis + hid a persistent pedestal on another. BE12599 2026-08-21 read "clean" + because Long had a 1.065 in/s event, while Tran sat at +0.47 in/s. + 2. It used the MEAN, which a real transient perturbs. The median is the + resting baseline: most samples sit at it, so a blast does not move it. + Same event, Long: mean +0.0783 vs median -0.0050. + +Flags a CHANNEL when |median| >= --floor in/s (default 0.025 = 5 A/D counts, +Instantel's own criterion; 1 A/D count = 0.005 in/s). + +Emits one row per (event, channel) so persistence can be tracked per channel. +""" +from __future__ import annotations +import argparse, csv, re, statistics, sys +from concurrent.futures import ProcessPoolExecutor, as_completed +from pathlib import Path + +sys.path.insert(0, str(Path(__file__).resolve().parent.parent)) +from minimateplus.event_file_io import read_blastware_file + +GEO = ("Tran", "Vert", "Long") +K = 10.0 / 32000.0 # ADC counts -> in/s at the 10 in/s range +_WAVE_RE = re.compile(r"\.[A-Za-z0-9]{2}0[Ww]$") +_STEM_RE = re.compile(r"^([B-Z])(\d{3})") + + +def serial_from_name(n): + m = _STEM_RE.match(n) + return f"BE{(ord(m.group(1))-ord('B'))*1000+int(m.group(2))}" if m else "?" + + +def scan_one(ps): + p = Path(ps) + try: + ev = read_blastware_file(p) + s = ev.raw_samples or {} + if not all(s.get(c) for c in GEO): + return None + ts = ev.timestamp + stamp = (f"{ts.year:04d}-{ts.month:02d}-{ts.day:02d}T" + f"{ts.hour:02d}:{ts.minute:02d}:{ts.second:02d}") if ts else "" + out = [] + for ch in GEO: + a = s[ch] + out.append({ + "serial": serial_from_name(p.name), "timestamp": stamp, + "filename": p.name, "channel": ch, + "median_ips": round(statistics.median(a) * K, 4), + "mean_ips": round(statistics.fmean(a) * K, 4), + "peak_ips": round(max(abs(v) for v in a) * K, 4), + }) + return out + except Exception: + return None + + +def main(): + ap = argparse.ArgumentParser() + ap.add_argument("--dir", required=True) + ap.add_argument("--jobs", type=int, default=4) + ap.add_argument("--floor", type=float, default=0.025) + ap.add_argument("--out", required=True) + a = ap.parse_args() + + seen, files = set(), [] + for q in sorted(Path(a.dir).rglob("*")): + if q.is_file() and _WAVE_RE.search(q.name) and q.name not in seen: + seen.add(q.name); files.append(str(q)) + print(f"unique waveform binaries: {len(files)}", flush=True) + + rows = [] + with ProcessPoolExecutor(max_workers=a.jobs) as ex: + for n, f in enumerate(as_completed([ex.submit(scan_one, p) for p in files]), 1): + r = f.result() + if r: rows.extend(r) + if n % 2000 == 0: print(f" {n}/{len(files)}", flush=True) + + for r in rows: + r["offset"] = int(abs(r["median_ips"]) >= a.floor) + + cols = ["serial","timestamp","filename","channel","median_ips","mean_ips","peak_ips","offset"] + with open(a.out, "w", newline="") as fh: + w = csv.DictWriter(fh, fieldnames=cols); w.writeheader(); w.writerows(rows) + + from collections import defaultdict + ev_flagged = {(r["serial"], r["filename"]) for r in rows if r["offset"]} + ev_all = {(r["serial"], r["filename"]) for r in rows} + per = defaultdict(set) + for r in rows: + if r["offset"]: per[r["serial"]].add(r["filename"]) + tot = defaultdict(set) + for r in rows: tot[r["serial"]].add(r["filename"]) + print(f"\nfloor = {a.floor} in/s ({a.floor/0.005:.0f} A/D counts)") + print(f"events with >=1 offset channel: {len(ev_flagged)} of {len(ev_all)}") + print(f"units affected: {len(per)} of {len(tot)}") + for s in sorted(per, key=lambda s: -len(per[s])): + print(f" {s:9} {len(per[s]):4} / {len(tot[s]):4} events") + print(f"\nwrote {a.out}") + + +if __name__ == "__main__": + main() diff --git a/scratch/offset_scan3.py b/scratch/offset_scan3.py new file mode 100644 index 0000000..2c2b953 --- /dev/null +++ b/scratch/offset_scan3.py @@ -0,0 +1,95 @@ +#!/usr/bin/env python3 +"""Offset detector v3 — pre-trigger floor, with pre/mid/end consistency. + +Brian's method, and better than v2's whole-record median for one reason: the +pre-trigger window is *definitionally* quiet (it is the buffer captured before +the trigger fired), whereas a whole-record median is merely robust to the event. + +Per channel: + pre = median of the first `pretrig_samples` samples (STRT record) + mid = median of the middle third + end = median of the final third + spread = max(pre,mid,end) - min(pre,mid,end) + +A DC offset is a *constant floor*: |pre| at or above the floor AND a small +spread. A transient (settling, handling, a long-tailed event) moves one segment +relative to the others and is rejected by the spread test. + +Floor default 0.025 in/s = 5 A/D counts (Instantel's own criterion; 1 count = +0.005 in/s). Quantisation is 0.005 in/s, so `spread` is measured in units of it. +""" +from __future__ import annotations +import argparse, csv, re, statistics, sys +from concurrent.futures import ProcessPoolExecutor, as_completed +from pathlib import Path +sys.path.insert(0, str(Path(__file__).resolve().parent.parent)) +from minimateplus.event_file_io import read_blastware_file + +GEO=("Tran","Vert","Long"); K=10.0/32000.0 +_WAVE=re.compile(r"\.[A-Za-z0-9]{2}0[Ww]$"); _STEM=re.compile(r"^([B-Z])(\d{3})") + +def serial_of(n): + m=_STEM.match(n) + return f"BE{(ord(m.group(1))-ord('B'))*1000+int(m.group(2))}" if m else "?" + +def scan(ps): + p=Path(ps) + try: + ev=read_blastware_file(p); s=ev.raw_samples or {} + if not all(s.get(c) for c in GEO): return None + pre_n=ev.pretrig_samples + ts=ev.timestamp + stamp=(f"{ts.year:04d}-{ts.month:02d}-{ts.day:02d}T" + f"{ts.hour:02d}:{ts.minute:02d}:{ts.second:02d}") if ts else "" + out=[] + for ch in GEO: + a=s[ch]; n=len(a); t=n//3 + pre = a[:pre_n] if (pre_n and 0 < pre_n < n) else a[:t] + mid, end = a[t:2*t], a[2*t:] + if not pre or not mid or not end: continue + v=[statistics.median(x)*K for x in (pre,mid,end)] + out.append({"serial":serial_of(p.name),"timestamp":stamp, + "filename":p.name,"channel":ch, + "pretrig_n": pre_n or 0, + "pre":round(v[0],4),"mid":round(v[1],4),"end":round(v[2],4), + "spread":round(max(v)-min(v),4), + "peak":round(max(abs(x) for x in a)*K,4)}) + return out + except Exception: + return None + +def main(): + ap=argparse.ArgumentParser() + ap.add_argument("--dir",required=True); ap.add_argument("--jobs",type=int,default=4) + ap.add_argument("--floor",type=float,default=0.025) + ap.add_argument("--max-spread",type=float,default=0.02) + ap.add_argument("--out",required=True) + a=ap.parse_args() + seen=set(); files=[] + for q in sorted(Path(a.dir).rglob("*")): + if q.is_file() and _WAVE.search(q.name) and q.name not in seen: + seen.add(q.name); files.append(str(q)) + print(f"unique waveform binaries: {len(files)}",flush=True) + rows=[] + with ProcessPoolExecutor(max_workers=a.jobs) as ex: + for i,f in enumerate(as_completed([ex.submit(scan,p) for p in files]),1): + r=f.result() + if r: rows.extend(r) + if i%2000==0: print(f" {i}/{len(files)}",flush=True) + for r in rows: + r["offset"]=int(abs(r["pre"])>=a.floor and r["spread"]<=a.max_spread) + cols=["serial","timestamp","filename","channel","pretrig_n","pre","mid","end","spread","peak","offset"] + with open(a.out,"w",newline="") as fh: + w=csv.DictWriter(fh,fieldnames=cols); w.writeheader(); w.writerows(rows) + from collections import defaultdict + per=defaultdict(set); tot=defaultdict(set) + for r in rows: + tot[r["serial"]].add(r["filename"]) + if r["offset"]: per[r["serial"]].add(r["filename"]) + print(f"\nfloor={a.floor} in/s ({a.floor/0.005:.0f} counts) max spread={a.max_spread}") + print(f"units affected: {len(per)} of {len(tot)}") + for s in sorted(per,key=lambda s:-len(per[s])): + print(f" {s:9} {len(per[s]):4} / {len(tot[s]):4} events") + print(f"\nwrote {a.out}") + +if __name__=="__main__": main() diff --git a/scratch/verify_against_ascii.py b/scratch/verify_against_ascii.py new file mode 100644 index 0000000..da14211 --- /dev/null +++ b/scratch/verify_against_ascii.py @@ -0,0 +1,221 @@ +#!/usr/bin/env python3 +"""Verify the series-3 decoder against preserved Blastware ASCII exports. + +Pairs each `__ASCII.TXT` with its binary `.`, decodes the +binary with the production codec, and compares against BW's own export: + + waveform — per-channel sample counts, then every sample value + histogram — interval count, then every per-interval channel peak + +ADC counts convert as ips = counts * geo_range_ips / 32000 (1 decoder unit = +16 counts = 0.005 in/s at the 10 in/s range; see CLAUDE.md). + +Usage: + python scratch/verify_against_ascii.py --dir [--limit N] [--jobs N] + [--out results.json] [--kind w|h|all] +""" +from __future__ import annotations + +import argparse, json, re, sys, traceback +from concurrent.futures import ProcessPoolExecutor, as_completed +from pathlib import Path + +sys.path.insert(0, str(Path(__file__).resolve().parent.parent)) + +from minimateplus.event_file_io import read_blastware_file +from minimateplus.bw_ascii_report import parse_report_file, parse_report + +GEO = ("Tran", "Vert", "Long") +_GEO_FS_COUNTS = 32000.0 +_ASCII_SUFFIX_RE = re.compile(r"_ASCII\.TXT$", re.IGNORECASE) + + +def binary_for(ascii_path: Path) -> Path: + """H907KXOW_WC0H_ASCII.TXT -> H907KXOW.WC0H""" + stem = _ASCII_SUFFIX_RE.sub("", ascii_path.name) + if "_" not in stem: + return ascii_path.with_name(stem) + head, _, ext = stem.rpartition("_") + return ascii_path.with_name(f"{head}.{ext}") + + +def counts_to_ips(counts, geo_range_ips): + r = geo_range_ips if geo_range_ips else 10.0 + return counts * r / _GEO_FS_COUNTS + + +def agrees(got, exp, geo_range_ips, tol=0.0006): + """True when decoded `got` matches BW's exported `exp`. + + Saturation carve-out: when an event clips, BW clamps its export to the + channel's range maximum (and writes OORANGE for the summary PPV), while + the decoder faithfully reproduces raw counts that can sit a decoder unit + or two past nominal full scale (32016 counts observed = 10.005 in/s on + the 10 in/s range). Same sign and both at/above the ceiling is agreement, + not a decode error. + """ + if abs(got - exp) <= tol: + return True + r = geo_range_ips if geo_range_ips else 10.0 + if abs(exp) >= r - tol and abs(got) >= r - tol and (got >= 0) == (exp >= 0): + return True + return False + + +def parse_interval_table(text: str): + """Histogram interval rows: time, Tpk, Tfq, Vpk, Vfq, Lpk, Lfq, PVS, ..., micdB, micfq""" + rows = [] + seen_header = False + for line in text.splitlines(): + if "\t" not in line: + continue + cols = [c.strip().strip('"') for c in line.split("\t")] + cols = [c for c in cols if c != ""] + if not seen_header: + if any(c in ("Tran", "Vert", "Long") for c in cols): + seen_header = True + continue + if len(cols) < 7: + continue + if not re.match(r"^\d{1,2}:\d{2}:\d{2}$", cols[0]): + continue + def num(s): + try: + return float(s) + except ValueError: + return None + rows.append({"time": cols[0], "Tran": num(cols[1]), + "Vert": num(cols[3]), "Long": num(cols[5])}) + return rows + + +def check_one(ascii_path_str: str) -> dict: + ap = Path(ascii_path_str) + bp = binary_for(ap) + res = {"ascii": ap.name, "binary": bp.name, "status": "?", + "kind": None, "detail": ""} + try: + if not bp.exists(): + res["status"] = "no_binary" + return res + text = ap.read_text(errors="replace") + rep = parse_report(text, parse_samples=True) + ev = read_blastware_file(bp) + gr = rep.geo_range_ips + res["kind"] = kind = ("histogram" + if (rep.event_type or "").lower().startswith(("full histogram", "histogram")) + else "waveform") + samples = ev.raw_samples or {} + dec_n = {c: len(samples.get(c) or []) for c in GEO} + + if kind == "histogram": + rows = parse_interval_table(text) + res["n_ascii"] = len(rows) + res["n_decoded"] = dec_n["Tran"] + if not rows: + res["status"] = "no_ascii_table" + return res + if dec_n["Tran"] == 0: + res["status"] = "decode_empty" + return res + if dec_n["Tran"] != len(rows): + res["status"] = "count_mismatch" + res["detail"] = f"decoded {dec_n['Tran']} vs ascii {len(rows)}" + return res + bad = 0 + worst = 0.0 + for i, row in enumerate(rows): + for ch in GEO: + exp = row[ch] + if exp is None: + continue + got = counts_to_ips(samples[ch][i], gr) + if not agrees(got, exp, gr): + bad += 1 + worst = max(worst, abs(got - exp)) + res["worst_abs"] = round(worst, 6) + res["status"] = "exact" if bad == 0 else "value_mismatch" + if bad: + res["detail"] = f"{bad} interval-channel values off" + return res + + # waveform + asc = rep.samples or [] + res["n_ascii"] = len(asc) + res["n_decoded"] = dec_n["Tran"] + if not asc: + res["status"] = "no_ascii_table" + return res + if dec_n["Tran"] == 0: + res["status"] = "decode_empty" + return res + if len({dec_n[c] for c in GEO}) != 1: + res["status"] = "channel_len_mismatch" + res["detail"] = str(dec_n) + return res + if dec_n["Tran"] != len(asc): + res["status"] = "count_mismatch" + res["detail"] = f"decoded {dec_n['Tran']} vs ascii {len(asc)}" + return res + bad = 0 + worst = 0.0 + for i, quad in enumerate(asc): + for j, ch in enumerate(GEO): + exp = quad[j] + got = counts_to_ips(samples[ch][i], gr) + if not agrees(got, exp, gr): + bad += 1 + worst = max(worst, abs(got - exp)) + res["worst_abs"] = round(worst, 6) + res["status"] = "exact" if bad == 0 else "value_mismatch" + if bad: + res["detail"] = f"{bad} sample values off" + return res + except Exception as e: + res["status"] = "error" + res["detail"] = f"{type(e).__name__}: {e}" + return res + + +def main(): + ap = argparse.ArgumentParser() + ap.add_argument("--dir", required=True) + ap.add_argument("--limit", type=int, default=0) + ap.add_argument("--jobs", type=int, default=8) + ap.add_argument("--kind", choices=["w", "h", "all"], default="all") + ap.add_argument("--out", default=None) + a = ap.parse_args() + + root = Path(a.dir) + files = sorted(p for p in root.rglob("*") + if p.is_file() and p.name.upper().endswith("_ASCII.TXT")) + if a.kind != "all": + want = "0W" if a.kind == "w" else "0H" + files = [p for p in files + if _ASCII_SUFFIX_RE.sub("", p.name).upper().endswith(want)] + if a.limit: + files = files[: a.limit] + print(f"pairs to check: {len(files)}", flush=True) + + out = [] + from collections import Counter + tally = Counter() + with ProcessPoolExecutor(max_workers=a.jobs) as ex: + futs = {ex.submit(check_one, str(p)): p for p in files} + for n, f in enumerate(as_completed(futs), 1): + r = f.result() + out.append(r) + tally[(r["kind"], r["status"])] += 1 + if n % 500 == 0: + print(f" {n}/{len(files)}", flush=True) + + print("\n=== results ===") + for (kind, status), n in sorted(tally.items(), key=lambda x: -x[1]): + print(f" {str(kind):10} {status:22} {n}") + if a.out: + Path(a.out).write_text(json.dumps(out, indent=1)) + print(f"\nwrote {a.out}") + + +if __name__ == "__main__": + main() diff --git a/sfm/database.py b/sfm/database.py index 6986208..04c5d24 100644 --- a/sfm/database.py +++ b/sfm/database.py @@ -603,56 +603,107 @@ class SeismoDb: ).fetchall() return [dict(r) for r in rows] - def find_twins(self, event_id: str, *, window_seconds: int = 300) -> list[dict]: + def find_twins(self, event_id: str, *, window_seconds: int | None = None) -> list[dict]: """ - Find this event's histogram/waveform twin(s): rows sharing the same - serial and an identical peak_vector_sum, whose timestamp falls - within ``window_seconds`` of this event's timestamp. Excludes the - event itself. Returns [] if the event or any required field - (serial / peak_vector_sum / timestamp) is missing. + Find this event's histogram/waveform twin(s): the SAME physical event + recorded both as a scheduled histogram and as a triggered waveform. - Caveat: identical-PVS matching is a proxy for "same physical event - recorded twice," not a guarantee. In the rare case where the device - clamps/saturates PVS (clamped to sqrt(3) * geo_range), two distinct - saturated events on the same serial within the window can share the - same clamped PVS value and be matched as twins even though they are - different events. This is harmless in practice — false_trigger/ - reviewed_real are derived/index columns re-derivable from the - sidecar source of truth — but worth knowing if twin counts look - surprising on a saturated/clamped run. + A real trigger is captured twice — once as a triggered waveform (stamped + at the trigger instant) and once inside the scheduled histogram whose + interval contains it (stamped at the histogram's interval start, e.g. the + 7am/7pm call-in). The two can be HOURS apart in time yet report the same + serial and identical peak_vector_sum. Twins are therefore matched by: + + * same serial, + * identical peak_vector_sum, + * OPPOSITE record type (one histogram, one waveform), and + * the waveform's timestamp falls within the histogram's interval — + from a histogram's timestamp up to the next histogram (same serial). + + This replaces the old ±``window_seconds`` heuristic, which silently + missed twins more than a few minutes apart (a histogram's interval-start + stamp and the trigger instant routinely differ by hours). ``window_seconds`` + is still accepted for backward compatibility but is ignored. + + Returns [] if the event or a required field (serial / peak_vector_sum / + timestamp) is missing. + + Caveat: identical-PVS matching remains a proxy for "same physical event" + — if the device clamps/saturates PVS (to sqrt(3) * geo_range), two + distinct saturated events could share a PVS. The added opposite-type and + interval constraints make a false pairing far less likely than the old + time-window match, and false_trigger/reviewed_real stay re-derivable from + the sidecar source of truth. """ + def _parse(ts): + if not ts: + return None + try: + return datetime.datetime.fromisoformat(str(ts).replace(" ", "T")) + except ValueError: + return None + + def _is_hist(rt): + return str(rt or "").lower().startswith("hist") + row = self.get_event(event_id) if not row: return [] - serial = row.get("serial"); pvs = row.get("peak_vector_sum"); ts = row.get("timestamp") - if serial is None or pvs is None or not ts: + serial = row.get("serial"); pvs = row.get("peak_vector_sum") + t_target = _parse(row.get("timestamp")) + if serial is None or pvs is None or t_target is None: return [] - try: - t = datetime.datetime.fromisoformat(ts.replace(" ", "T")) - except ValueError: - return [] - lo = (t - datetime.timedelta(seconds=window_seconds)).isoformat() - hi = (t + datetime.timedelta(seconds=window_seconds)).isoformat() - with self._connect() as conn: - rows = conn.execute( - "SELECT * FROM events WHERE serial=? AND id!=? AND peak_vector_sum=? " - "AND timestamp BETWEEN ? AND ?", - (serial, event_id, pvs, lo, hi), - ).fetchall() - return [dict(r) for r in rows] + target_hist = _is_hist(row.get("record_type")) - def propagate_review_to_twins(self, event_id: str, *, window_seconds: int = 300) -> list[str]: + with self._connect() as conn: + cand_rows = [dict(r) for r in conn.execute( + "SELECT * FROM events WHERE serial=? AND id!=? AND peak_vector_sum=?", + (serial, event_id, pvs)).fetchall()] + hist_ts = [r["timestamp"] for r in conn.execute( + "SELECT timestamp FROM events WHERE serial=? AND lower(record_type) LIKE 'hist%'", + (serial,)).fetchall()] + + # Histogram interval-start times for this serial, sorted, to bound intervals. + starts = sorted(x for x in (_parse(t) for t in hist_ts) if x is not None) + + def _interval_end(h_start): + # The next histogram strictly after h_start bounds the interval; else open-ended. + for x in starts: + if x > h_start: + return x + return None + + def _covers(h_start, w_time): + end = _interval_end(h_start) + return h_start <= w_time and (end is None or w_time < end) + + twins = [] + for c in cand_rows: + if _is_hist(c.get("record_type")) == target_hist: + continue # twins are strictly cross-type (one histogram, one waveform) + c_time = _parse(c.get("timestamp")) + if c_time is None: + continue + h_start, w_time = (t_target, c_time) if target_hist else (c_time, t_target) + if _covers(h_start, w_time): + twins.append(c) + return twins + + def propagate_review_to_twins(self, event_id: str, *, window_seconds: int | None = None) -> list[str]: """ Copy this event's `false_trigger`/`reviewed_real` columns onto each of its histogram/waveform twins (see `find_twins`), so flagging one twin flags both. Returns the list of twin ids updated. + + ``window_seconds`` is accepted for backward compatibility but ignored; + twin matching is now interval-based (see `find_twins`). """ row = self.get_event(event_id) if not row: return [] ft = 1 if row.get("false_trigger") else 0 real = 1 if row.get("reviewed_real") else 0 - twins = self.find_twins(event_id, window_seconds=window_seconds) + twins = self.find_twins(event_id) moved = [] with self._connect() as conn: for tw in twins: diff --git a/sfm/server.py b/sfm/server.py index 0004702..eeeae73 100644 --- a/sfm/server.py +++ b/sfm/server.py @@ -67,6 +67,7 @@ from minimateplus.blastware_file import write_blastware_file, blastware_filename from minimateplus.client import _decode_a5_metadata_into, _decode_a5_waveform, _decode_event_count from minimateplus.framing import build_bw_write_frame, SESSION_RESET, POLL_PROBE, POLL_DATA from minimateplus.protocol import SUB_STOP_MONITORING +from minimateplus.event_file_io import TOOL_VERSION as SFM_VERSION # single source for the service version (release-bumped) from sfm import event_hdf5 from sfm.cache import SFMCache, get_cache from sfm.database import SeismoDb @@ -90,7 +91,7 @@ app = FastAPI( "Implements the minimateplus RS-232 protocol library.\n" "Proxied by terra-view at /api/sfm/*." ), - version="0.26.0", + version=SFM_VERSION, ) # Allow requests from the waveform viewer opened as a local file (file://) @@ -371,7 +372,7 @@ def _backfill_events(events: list, info: "DeviceInfo") -> None: @app.get("/health") def health() -> dict: """Service heartbeat. No device I/O.""" - return {"status": "ok", "service": "sfm", "version": "0.1.0"} + return {"status": "ok", "service": "sfm", "version": SFM_VERSION} @app.get("/", response_class=FileResponse) diff --git a/tests/test_find_twins.py b/tests/test_find_twins.py index d0b356c..db501c1 100644 --- a/tests/test_find_twins.py +++ b/tests/test_find_twins.py @@ -1,33 +1,71 @@ -import datetime +import sqlite3 from sfm.database import SeismoDb from minimateplus.models import Event, Timestamp -def _ins(db, key, serial, pvs, ts): +def _ins(db, key, serial, pvs, ts, record_type="Waveform"): ev = Event(index=0) ev._waveform_key = bytes.fromhex(key) ev.timestamp = ts - # peak_vector_sum comes from peak_values; simplest: insert then UPDATE pvs directly db.insert_events([ev], serial=serial) row = [r for r in db.query_events(serial=serial) if r["waveform_key"] == key][0] - import sqlite3 with sqlite3.connect(db.db_path) as c: - c.execute("UPDATE events SET peak_vector_sum=? WHERE id=?", (pvs, row["id"])) + c.execute("UPDATE events SET peak_vector_sum=?, record_type=? WHERE id=?", + (pvs, record_type, row["id"])) return row["id"] -def test_find_twins_matches_same_serial_pvs_near_time(tmp_path): +def _ts(hour, minute, second=0, day=25): + return Timestamp(raw=b"", flag=0x10, year=2026, unknown_byte=0, month=2, day=day, + hour=hour, minute=minute, second=second) + + +def test_histogram_and_waveform_twin_across_hours(tmp_path): + # The real UM12947 case: histogram stamped at its 7pm interval start, the + # triggered waveform 75 min later — same serial + identical PVS. The old + # ±5-min window missed this; interval matching catches it, both directions. db = SeismoDb(tmp_path / "s.db") - base = Timestamp(raw=b"", flag=0x10, year=2026, unknown_byte=0, month=2, day=25, hour=20, minute=19, second=5) - twin = Timestamp(raw=b"", flag=0x10, year=2026, unknown_byte=0, month=2, day=25, hour=20, minute=19, second=45) - far = Timestamp(raw=b"", flag=0x10, year=2026, unknown_byte=0, month=2, day=25, hour=21, minute=0, second=0) - # d needs a timestamp distinct from `twin` (UNIQUE(serial, timestamp) would - # otherwise collide with b and UPSERT onto its row instead of inserting a - # new one) while staying near `base` in time. - near = Timestamp(raw=b"", flag=0x10, year=2026, unknown_byte=0, month=2, day=25, hour=20, minute=19, second=44) - a = _ins(db, "01110001", "BE1", 0.4763, base) - b = _ins(db, "01110002", "BE1", 0.4763, twin) # twin: same serial+pvs, 40s apart - c = _ins(db, "01110003", "BE1", 0.4763, far) # same pvs but >window away - d = _ins(db, "01110004", "BE1", 0.9999, near) # near time but different pvs - ids = {r["id"] for r in db.find_twins(a, window_seconds=300)} - assert ids == {b} + hist_pm = _ins(db, "01110001", "BE1", 0.4763, _ts(19, 31, 17), "Histogram") + wave = _ins(db, "01110002", "BE1", 0.4763, _ts(20, 46, 44), "Waveform") + _ins(db, "01110003", "BE1", 0.0100, _ts(7, 0, 0, day=26), "Histogram") # bounds the interval + assert {r["id"] for r in db.find_twins(hist_pm)} == {wave} + assert {r["id"] for r in db.find_twins(wave)} == {hist_pm} + + +def test_same_type_not_twinned(tmp_path): + # Two waveforms, same serial + PVS, seconds apart → NOT twins (cross-type only). + db = SeismoDb(tmp_path / "s.db") + a = _ins(db, "01110001", "BE1", 0.4763, _ts(20, 19, 5), "Waveform") + _ins(db, "01110002", "BE1", 0.4763, _ts(20, 19, 45), "Waveform") + assert db.find_twins(a) == [] + + +def test_waveform_matches_only_the_containing_interval(tmp_path): + # Two overnight intervals with the same PVS; a waveform in the SECOND interval + # must twin with that histogram, never the first — even though PVS matches both. + db = SeismoDb(tmp_path / "s.db") + h1 = _ins(db, "01110001", "BE1", 0.4763, _ts(19, 0, 0, day=25), "Histogram") + h2 = _ins(db, "01110002", "BE1", 0.4763, _ts(7, 0, 0, day=26), "Histogram") + w = _ins(db, "01110003", "BE1", 0.4763, _ts(8, 0, 0, day=26), "Waveform") + assert {r["id"] for r in db.find_twins(w)} == {h2} + assert w not in {r["id"] for r in db.find_twins(h1)} + + +def test_different_pvs_not_twinned(tmp_path): + db = SeismoDb(tmp_path / "s.db") + h = _ins(db, "01110001", "BE1", 0.4763, _ts(19, 0, 0), "Histogram") + _ins(db, "01110002", "BE1", 0.9999, _ts(20, 0, 0), "Waveform") # different PVS + assert db.find_twins(h) == [] + + +def test_open_ended_latest_interval(tmp_path): + # A waveform after the latest histogram (nothing bounds the interval) still twins. + db = SeismoDb(tmp_path / "s.db") + h = _ins(db, "01110001", "BE1", 0.4763, _ts(19, 0, 0), "Histogram") + w = _ins(db, "01110002", "BE1", 0.4763, _ts(23, 30, 0), "Waveform") + assert {r["id"] for r in db.find_twins(h)} == {w} + + +def test_missing_fields_returns_empty(tmp_path): + db = SeismoDb(tmp_path / "s.db") + assert db.find_twins("nonexistent-id") == [] diff --git a/tests/test_health_version.py b/tests/test_health_version.py new file mode 100644 index 0000000..f4be03c --- /dev/null +++ b/tests/test_health_version.py @@ -0,0 +1,17 @@ +"""The /health version must track the release, not a stale literal. + +terra-view's SFM Admin page displays whatever `/health` reports. It was +hardcoded to "0.1.0" and never bumped, so the page showed 0.1.0 while the +service was actually 0.26.0. These guard against that regression — and run +without httpx (they call the endpoint function directly, no TestClient). +""" +from minimateplus.event_file_io import TOOL_VERSION +from sfm.server import app, health + + +def test_health_reports_current_tool_version(): + assert health()["version"] == TOOL_VERSION + + +def test_openapi_version_matches_tool_version(): + assert app.version == TOOL_VERSION diff --git a/tests/test_histogram_codec.py b/tests/test_histogram_codec.py index 558a0c9..1d79f0f 100644 --- a/tests/test_histogram_codec.py +++ b/tests/test_histogram_codec.py @@ -643,3 +643,38 @@ def test_multi_interval_matches_blastware_ascii_exactly(): assert hz is None elif not cell.startswith("<"): assert hz is not None and abs(hz - float(cell)) <= max(0.55, float(cell) * 0.02) + + +def test_partial_final_block_is_not_disqualified_by_missing_third_header(): + """A body can exceed two strides yet hold only two real blocks. + + Regression for BE18193 `T193L0XM.CI0H` — 51 intervals at 2 s = one full + 30-interval block plus a 21-interval remainder, in a body long enough to + demand a third block header at ``2 * stride`` that does not exist. The + third-block confirmation used to be mandatory whenever the body was long + enough, so the correct stride was discarded and the file decoded to + nothing. A missing third header means end-of-stream, not disqualification; + the block-counter check is the decisive anti-false-positive test. + """ + full = [(1, 1, 2, 2, 3, 3, 4, 4)] * 30 + partial = [(5, 5, 6, 6, 7, 7, 8, 8)] * 21 + body = (_mk_multi_block(full, ctr=256) + + _mk_multi_block(partial, ctr=257) + + b"\xff" * 700) # trailing padding past 2 * stride + stride = 12 + 20 * 30 + assert 2 * stride + 6 <= len(body), "padding must reach past two strides" + # the whole point: a third header is absent, and that must not disqualify + assert detect_multi_interval_stride(body) == stride + recs = walk_multi_interval_blocks(body) + assert len(recs) == 51 + assert recs[0]["t_peak"] == 1 + assert recs[-1]["t_peak"] == 5 + + +def test_third_block_still_rejects_a_mismatched_counter(): + """The corroboration must still bite when a third block IS present.""" + ivs = [(1, 1, 2, 2, 3, 3, 4, 4)] * 4 + body = (_mk_multi_block(ivs, ctr=256) + + _mk_multi_block(ivs, ctr=257) + + _mk_multi_block(ivs, ctr=999)) # counter jumps — not consecutive + assert detect_multi_interval_stride(body) != 12 + 20 * 4 diff --git a/tests/test_twin_propagation.py b/tests/test_twin_propagation.py index c3f7386..6d08d17 100644 --- a/tests/test_twin_propagation.py +++ b/tests/test_twin_propagation.py @@ -3,31 +3,34 @@ from sfm.database import SeismoDb from minimateplus.models import Event, Timestamp -def _ins(db, key, serial, pvs, ts): +def _ins(db, key, serial, pvs, ts, record_type="Waveform"): ev = Event(index=0) ev._waveform_key = bytes.fromhex(key) ev.timestamp = ts - # peak_vector_sum comes from peak_values; simplest: insert then UPDATE pvs directly db.insert_events([ev], serial=serial) row = [r for r in db.query_events(serial=serial) if r["waveform_key"] == key][0] with sqlite3.connect(db.db_path) as c: - c.execute("UPDATE events SET peak_vector_sum=? WHERE id=?", (pvs, row["id"])) + c.execute("UPDATE events SET peak_vector_sum=?, record_type=? WHERE id=?", + (pvs, record_type, row["id"])) return row["id"] -def test_propagate_copies_flags_to_twins(tmp_path): +def _ts(hour, minute, second=0, day=25): + return Timestamp(raw=b"", flag=0x10, year=2026, unknown_byte=0, month=2, day=day, + hour=hour, minute=minute, second=second) + + +def test_propagate_copies_flags_across_hours_apart_twins(tmp_path): + # Flagging the waveform FT propagates to its histogram twin 75 min earlier + # (the interval matcher pairs them; the old ±5-min window would have missed it). db = SeismoDb(tmp_path / "s.db") - base = Timestamp(raw=b"", flag=0x10, year=2026, unknown_byte=0, month=2, day=25, hour=20, minute=19, second=5) - twin = Timestamp(raw=b"", flag=0x10, year=2026, unknown_byte=0, month=2, day=25, hour=20, minute=19, second=45) - other = Timestamp(raw=b"", flag=0x10, year=2026, unknown_byte=0, month=2, day=25, hour=20, minute=19, second=44) + hist = _ins(db, "01110001", "BE1", 0.4763, _ts(19, 31, 17), "Histogram") + wave = _ins(db, "01110002", "BE1", 0.4763, _ts(20, 46, 44), "Waveform") # twin, 75 min later + other = _ins(db, "01110003", "BE1", 0.9999, _ts(20, 20, 0), "Waveform") # different pvs - primary_id = _ins(db, "01110001", "BE1", 0.4763, base) - twin_id = _ins(db, "01110002", "BE1", 0.4763, twin) # twin: same serial+pvs, 40s apart - non_twin_id = _ins(db, "01110003", "BE1", 0.9999, other) # near time but different pvs + db.update_event_review(wave, {"false_trigger": True}) + moved = db.propagate_review_to_twins(wave) - db.update_event_review(primary_id, {"false_trigger": True}) - moved = db.propagate_review_to_twins(primary_id) - - assert twin_id in moved - assert db.get_event(twin_id)["false_trigger"] == 1 - assert db.get_event(non_twin_id)["false_trigger"] == 0 + assert hist in moved + assert db.get_event(hist)["false_trigger"] == 1 + assert db.get_event(other)["false_trigger"] == 0