serversdownandClaude Opus 4.8
2902ab373e
channel_spectrum(samples, sps) → the single-sided amplitude spectrum Blastware's FFT Report draws, and dominant_frequency() picks its peak in the 2–250 Hz band. Reverse-engineered against 7 BE12844 (MiniMate Plus) events with Blastware FFT reports as ground truth. Recipe: DC-remove, NO window (a window smears the peak and worsens the match), zero-pad to 4096 (→ 0.25 Hz bins at 1024 sps — the resolution every reported dominant frequency lands on), single-sided 2/N amplitude. Reproduces Blastware's dominant frequency to the exact bin on all 28 channels and the amplitude to report precision. This is the missing piece for both the USBM RI8507 compliance chart (its scatter is these (freq, amp) points vs the limit curve) and the FFT view. Pure numpy, series-agnostic (feed it in/s samples from either decoder). The 7 events land in tests/fixtures as the oracle (force-added past the fixtures gitignore, matching 5-11-26 / decode-re-5-8-26). Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01YDXjZCr4RqT2U3QvMDhgzf
86 lines
3.7 KiB
Python
86 lines
3.7 KiB
Python
"""Blastware-compatible channel FFT (waveform_fft).
|
|
|
|
Reverse-engineered 2026-09-14 against 7 BE12844 (MiniMate Plus) events, each with
|
|
a Blastware FFT report as ground truth. The recipe (DC-remove, no window,
|
|
zero-pad to 4096 → 0.25 Hz bins, single-sided 2/N amplitude) reproduces
|
|
Blastware's dominant frequency to the exact bin on all 28 channels and the
|
|
amplitude to report precision.
|
|
"""
|
|
from pathlib import Path
|
|
|
|
import numpy as np
|
|
|
|
from waveform_fft import channel_spectrum, dominant_frequency
|
|
from minimateplus.waveform_codec import decode_waveform_v2
|
|
|
|
FIXDIR = Path(__file__).parent / "fixtures" / "fft-oracle-2026-09-14"
|
|
GEO_LSB = 0.005 # 1 decode unit = 16 ADC counts = 0.005 in/s (series-3 Normal range)
|
|
|
|
# Blastware FFT-report ground truth: file → {channel: (dominant_hz, amplitude_ips)}.
|
|
# amplitude is None where the channel is at the noise floor (report amp 0.000/0.001)
|
|
# — the dominant frequency still matches exactly, but the amplitude isn't meaningful.
|
|
ORACLE = {
|
|
"N844LPGH.VV0W": {"Tran": (27.00, 0.018), "Vert": (26.75, 0.009), "Long": (26.50, 0.021), "MicL": (2.000, None)},
|
|
"N844LPPR.3S0W": {"Tran": (30.75, None), "Vert": (46.75, None), "Long": (26.75, None), "MicL": (49.50, None)},
|
|
"N844LQHB.ZT0W": {"Tran": (19.75, 0.040), "Vert": (26.50, 0.018), "Long": (26.50, 0.083), "MicL": (2.750, None)},
|
|
"N844LQUE.T50W": {"Tran": (21.50, 0.080), "Vert": (14.25, 0.028), "Long": (28.50, 0.046), "MicL": (5.750, None)},
|
|
"N844LR8W.790W": {"Tran": (31.00, None), "Vert": (31.00, None), "Long": (34.00, None), "MicL": (66.25, None)},
|
|
"N844LRCO.G60W": {"Tran": (32.25, 0.009), "Vert": (32.00, 0.005), "Long": (32.00, 0.008), "MicL": (32.00, None)},
|
|
"N844LRCW.F30W": {"Tran": (21.25, 0.010), "Vert": (42.25, 0.002), "Long": (21.25, 0.014), "MicL": (21.25, None)},
|
|
}
|
|
|
|
|
|
def test_pure_sine_frequency_and_amplitude():
|
|
# A pure sine at a bin-centre frequency (128 cycles over 4096 samples) has no
|
|
# leakage, so the single-sided 2/N normalisation returns the amplitude exactly.
|
|
sps, n, f0, amp = 1024.0, 4096, 32.0, 0.5
|
|
x = amp * np.sin(2 * np.pi * f0 * np.arange(n) / sps)
|
|
freqs, amps = channel_spectrum(x, sps=sps, nfft=4096)
|
|
fpk, apk = dominant_frequency(freqs, amps)
|
|
assert fpk == 32.0
|
|
assert abs(apk - amp) < 1e-3
|
|
|
|
|
|
def test_bin_resolution_is_quarter_hz():
|
|
freqs, _ = channel_spectrum(np.zeros(3328), sps=1024.0, nfft=4096)
|
|
assert abs((freqs[1] - freqs[0]) - 0.25) < 1e-9
|
|
|
|
|
|
def test_empty_input():
|
|
freqs, amps = channel_spectrum([])
|
|
assert len(freqs) == 0 and len(amps) == 0
|
|
|
|
|
|
def _spectra(fname):
|
|
raw = (FIXDIR / fname).read_bytes()
|
|
dec = decode_waveform_v2(raw[raw.find(b"STRT") + 21:])
|
|
out = {}
|
|
for ch, samples in dec.items():
|
|
ips = np.asarray(samples, float) * GEO_LSB
|
|
out[ch] = channel_spectrum(ips, sps=1024.0)
|
|
return out
|
|
|
|
|
|
def test_dominant_frequency_matches_blastware_exactly():
|
|
misses = []
|
|
for fname, chans in ORACLE.items():
|
|
spectra = _spectra(fname)
|
|
for ch, (want_hz, _) in chans.items():
|
|
got_hz, _ = dominant_frequency(*spectra[ch])
|
|
if abs(got_hz - want_hz) > 0.25:
|
|
misses.append(f"{fname}:{ch} got {got_hz} want {want_hz}")
|
|
assert not misses, "dominant-frequency mismatches:\n" + "\n".join(misses)
|
|
|
|
|
|
def test_amplitude_matches_blastware():
|
|
misses = []
|
|
for fname, chans in ORACLE.items():
|
|
spectra = _spectra(fname)
|
|
for ch, (_, want_amp) in chans.items():
|
|
if want_amp is None:
|
|
continue
|
|
_, got_amp = dominant_frequency(*spectra[ch])
|
|
if abs(got_amp - want_amp) > 0.0015:
|
|
misses.append(f"{fname}:{ch} got {got_amp:.4f} want {want_amp:.3f}")
|
|
assert not misses, "amplitude mismatches:\n" + "\n".join(misses)
|