Merge remote-tracking branch 'origin/main' into main

# Conflicts:
#	sras_compute.py
#	sras_viewer.py
#	sras_workers.py
#	tools/make_test_sras.py
#	tools/test_refactor.py
This commit is contained in:
Thomas Ales [M S E]
2026-08-10 18:58:21 -05:00
34 changed files with 10768 additions and 5738 deletions
+140 -8
View File
@@ -18,9 +18,12 @@ import numpy as np
sys.path.insert(0, str(Path(__file__).resolve().parent.parent))
# Packed straight from sras_format.py's own constants (not a hand-copy) so
# this generator can never silently drift from what the real parser expects.
from sras_format import GEO_FMT_V6, HDR_FMT as HDR_FMT_LEGACY, HDR_FMT_V6 # noqa: E402
# Single source of truth for the byte layout: the reader's own constants.
# The byte *assembly* below stays independent, so a writer bug can't be
# masked by a matching reader bug.
from sras_format import HDR_FMT as HDR_FMT_LEGACY # noqa: E402
from sras_format import GEO_FMT_V6, HDR_FMT_V6 # noqa: E402
from sras_compute import _rotation_matrix as _rot # noqa: E402
# Per-angle (n_rows, n_frames) — deliberately different per angle so ragged
# geometry handling is actually exercised.
@@ -46,12 +49,12 @@ def _preamble(ymult_v: float, yoff_adc: float, yzero_v: float) -> bytes:
def build(n_angles: int, seed: int, samples_per_frame: int,
geometry: list[tuple[int, int]] | None = None) -> tuple[bytes, dict]:
geometry: list[tuple[int, int]] | None = None,
bps: int = 1) -> tuple[bytes, dict]:
rng = np.random.default_rng(seed)
src_geom = geometry or _GEOMETRY
geom = [src_geom[a % len(src_geom)] for a in range(n_angles)]
n_ch = 3
bps = 1
angles_deg = np.linspace(0.0, 60.0, n_angles, dtype=np.float32)
# Distinct calibration per channel so a swapped-channel bug is visible.
@@ -104,7 +107,9 @@ def build(n_angles: int, seed: int, samples_per_frame: int,
block[r, 1, f] = np.int8((a * 7 + r * 3 + f) % 100 - 50)
block[r, 2, f] = np.int8((a * 5 + r * 11 + f * 2) % 120 - 60)
waveforms.append(block)
out += block.tobytes()
# bps=2 stores the same values big-endian int16, exercising the
# reader's >i2 memmap path.
out += (block.astype(">i2") if bps == 2 else block).tobytes()
meta = {
"n_angles": n_angles,
@@ -123,12 +128,139 @@ def build(n_angles: int, seed: int, samples_per_frame: int,
def write(path: Path, n_angles: int = 3, seed: int = 0,
samples_per_frame: int = 64,
geometry: list[tuple[int, int]] | None = None) -> dict:
payload, meta = build(n_angles, seed, samples_per_frame, geometry)
geometry: list[tuple[int, int]] | None = None,
bps: int = 1) -> dict:
payload, meta = build(n_angles, seed, samples_per_frame, geometry, bps=bps)
path.write_bytes(payload)
return meta
# ---------------------------------------------------------------------------
# Rotating-sample scan: one shape, imaged at several known rotations
# ---------------------------------------------------------------------------
#
# The scan the angle-alignment path actually has to solve: every angle images
# the *same* sample at a different known rotation and offset, and a correct
# alignment stacks them all back into one shape. Two properties are
# deliberately hostile:
#
# * every angle gets a different window size and a different, meaningless
# stage x_start / y0 — alignment must ignore per-angle stage coordinates
# entirely, so any code that reads them will visibly fail here;
# * the pixel grid is strongly anisotropic (5 µm along x, 50 µm along y),
# like the real instrument, so any registration that rotates raw indices
# instead of millimetres shears the image and cannot converge.
_ROT_DX_MM = 0.005 # x pitch, from velocity/laser_freq below
_ROT_DY_MM = 0.05 # row spacing
_ROT_BG_MV = 4.0
_ROT_FG_MV = 160.0
# How far the sample sits from the rotation axis. Non-zero on purpose: on the
# real instrument every angle's scan window is centred on the rotation axis
# while the sample is not, so each scan sees the sample somewhere else along a
# circle. That offset is exactly what a wrong rotation pivot turns into a ring
# of scans instead of a stack, so a centred test sample would hide the bug.
_ROT_SAMPLE_OFFSET_MM = (0.55, 0.40)
def _sample_shape_mv(u: np.ndarray, v: np.ndarray) -> np.ndarray:
"""An asymmetric test sample in its own mm frame, chirally distinct at
every rotation (no 180° ambiguity) and with structure at several radii so
rotation is well determined."""
u = u - _ROT_SAMPLE_OFFSET_MM[0]
v = v - _ROT_SAMPLE_OFFSET_MM[1]
img = np.full(u.shape, _ROT_BG_MV, dtype=np.float32)
img[((u / 0.85) ** 2 + (v / 0.40) ** 2) <= 1.0] = _ROT_FG_MV # bar
img[(np.abs(u - 0.55) <= 0.22) & (np.abs(v - 0.62) <= 0.22)] = _ROT_FG_MV # nub
img[((u + 0.75) ** 2 + (v + 0.30) ** 2) <= 0.20 ** 2] = _ROT_FG_MV # dot
return img
def write_rotating(path: Path, n_angles: int = 5, samples_per_frame: int = 4,
seed: int = 0) -> dict:
"""Write a v6 file whose CH4 DC image is one sample seen at n_angles known
rotations, and return the ground truth each angle should register to.
``truth[a] = (rotation_deg, (shift_x_mm, shift_y_mm))`` is the rigid map
from angle *a*'s local mm (origin at its own array center) to angle 0's —
exactly what ``register_angle_to_reference`` is supposed to recover.
"""
rng = np.random.default_rng(seed)
n_ch, bps = 3, 1
cal = [(1.5625e-3, -87.04, 0.0), (2.0e-3, -60.0, 1.0e-3), (2.5e-3, -40.0, -2.0e-3)]
ymult_mv, yoff, yzero_mv = cal[2][0] * 1000, cal[2][1], cal[2][2] * 1000
stage_angles, geom, x_starts, y_starts, thetas, offsets = [], [], [], [], [], []
for a in range(n_angles):
stage = -37.0 * a # what the rotation stage reports
stage_angles.append(stage)
# The true image rotation is the negative of the stage's reported
# angle: the stage's positive sense is the opposite of math-positive
# (x toward y) in scan mm. Nothing may depend on knowing that — the
# registration search tries both signs.
thetas.append(-stage)
offsets.append((0.0, 0.0) if a == 0
else (float(rng.uniform(-0.3, 0.3)), float(rng.uniform(-0.3, 0.3))))
# A different window per angle, all centred on the same array center —
# the real instrument grows each angle's axis-aligned bounding box to
# cover the rotated ROI. Sized so the off-axis sample stays inside every
# window at every angle, keeping the expected result unambiguous.
geom.append((88 + 8 * a, 780 + 60 * a))
# Meaningless per-angle stage positions: correct alignment never reads
# them, so scattering them proves it.
x_starts.append(float(20.0 + rng.uniform(-6.0, 6.0)))
y_starts.append(float(30.0 + rng.uniform(-6.0, 6.0)))
out = bytearray()
out += struct.pack(
HDR_FMT_V6, b"SRAS", 6, n_angles,
x_starts[0], y_starts[0], 1.0, 1.0, _ROT_DY_MM,
_VELOCITY_MM_S, _VELOCITY_MM_S / _ROT_DX_MM, # velocity/freq -> 5 µm pitch
samples_per_frame, _SAMPLE_RATE_HZ, bps, n_ch,
)
out += np.array(stage_angles, dtype=">f4").tobytes()
for a, (n_rows, n_frames) in enumerate(geom):
out += struct.pack(GEO_FMT_V6, x_starts[a], 1.0, n_frames, n_rows)
for a, (n_rows, _) in enumerate(geom):
out += (y_starts[a] + np.arange(n_rows) * _ROT_DY_MM).astype(">f4").tobytes()
for ymult_v, yoff_a, yzero_v in cal:
p = _preamble(ymult_v, yoff_a, yzero_v)
out += struct.pack(">H", len(p)) + p
background = rng.integers(-8, 9, size=samples_per_frame, dtype=np.int8)
out += struct.pack(">I", samples_per_frame) + background.tobytes()
truth, dc4_images = {}, []
for a, (n_rows, n_frames) in enumerate(geom):
# Local mm of every pixel, measured from this angle's own array center.
lx = (np.arange(n_frames) - (n_frames - 1) / 2.0) * _ROT_DX_MM
ly = (np.arange(n_rows) - (n_rows - 1) / 2.0) * _ROT_DY_MM
gx, gy = np.meshgrid(lx, ly)
# local = R(theta) @ sample + offset, so sample = R(theta)^T @ (local - offset)
rel = np.stack([gx - offsets[a][0], gy - offsets[a][1]], axis=-1)
s = rel @ _rot(thetas[a]) # == rel @ R^T.T == R^T @ rel
dc4 = _sample_shape_mv(s[..., 0], s[..., 1])
dc4_images.append(dc4)
inv = _rot(-thetas[a])
truth[a] = (-thetas[a],
tuple(float(v) for v in -(inv @ np.array(offsets[a]))))
adc4 = np.clip(np.round((dc4 - yzero_mv) / ymult_mv + yoff), -128, 127).astype(np.int8)
block = np.zeros((n_rows, n_ch, n_frames, samples_per_frame), dtype=np.int8)
block[:, 2] = adc4[:, :, None] # CH4 carries the sample
block[:, 1] = 10 # CH3 flat
block[:, 0] = rng.integers(-40, 41, size=(n_rows, n_frames, samples_per_frame),
dtype=np.int8) # CH1 noise
out += block.tobytes()
path.write_bytes(bytes(out))
return {"n_angles": n_angles, "geometry": geom, "stage_angles_deg": stage_angles,
"truth": truth, "dc4_mv": dc4_images, "x_starts": x_starts,
"y_starts": y_starts, "dx_mm": _ROT_DX_MM, "dy_mm": _ROT_DY_MM}
def write_legacy(path: Path, version: int = 4, n_angles: int = 2,
n_rows: int = 4, n_frames: int = 10,
samples_per_frame: int = 32, seed: int = 0) -> dict: