feat(offset): detector v3 — pre-trigger floor with a constant-floor test
Brian's method, and better than v2's whole-record median: the pre-trigger
window is definitionally quiet (the buffer captured before the trigger fired),
whereas a median is merely robust to the event. Requiring the floor to hold
across pre-trigger / middle / end rejects transients that a median cannot.
per channel: pre/mid/end medians, spread = max - min
offset when |pre| >= floor AND spread <= 0.02 in/s
real fault >= 3 consecutive flagged events on that channel
The empirical noise floor justifies the threshold and validates the decoder:
across 19,244 non-flagged channel-events the pre-trigger floor is 62.7% exactly
0.000, 94.5% within +/-1 quantisation unit, median +0.0000, mean -0.0008. There
is no systematic zero-point bias — an independent confirmation of the
32000-count geo scale.
The result is threshold-insensitive across a 2x range (0.020 to 0.040 in/s),
which is what separates a real signal from a tuned one:
FINAL: 5 of 45 units (11%) — BE9558, BE11529, BE12599, BE13117, BE18438
BE11007 and BE10895 drop out; the spread test identifies them as transients
rather than pedestals. v1's 11% headline was right by luck — it included
BE11007 and named the wrong channel on most units.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01HgTe8CamXAHcAmaQ6QNcog
This commit is contained in:
@@ -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()
|
||||
Reference in New Issue
Block a user