Documentation

Beamformer pipeline

From a raw cube to a focused image: TVG, matched filter, micronavigation, back-projection, and display.

End-to-end SAS processing chain used by this repository: HDF5 -> focused SLC -> DRC'd PNG. Source of truth is beamformer/sas/pipeline.py::beamform_hdf5; wrapped by beamformer/sas/slc.py::beamform_tdbp.

Hard rule: always call beamform_tdbp

Do not write a custom pipeline. This rule exists because a hand-rolled LP+MF+BF script silently dropped TVG normalization and invalidated every simulator-vs-real comparison until the shared pipeline was restored. Any new experiment adds a flag to BeamformSettings; it does not fork the pipeline.

from beamformer.sas.slc import beamform_tdbp

img = beamform_tdbp(
    'scene066.h5',
    pixel_spacing=0.0125,                     # 1.25 cm; None -> c/(4*BW)
    extra={'range_max': 100.0, 'use_micronav': False},
)
img.save_port_image('out.png')

Data flow

HDF5 (POSSM schema)
  -> load raw cube (n_pings, n_channels, n_samples) complex64
  -> GPU lowpass FIR (zero-phase, cutoff = BW/2, 65 taps)
  -> GPU matched filter (FFT-based, from baseband LFM replica)
  -> Cross-ping TVG normalization                 [optional]
  -> 8x FFT range upsampling                       [BF_RANGE_UPSAMPLE env]
  -> Glint suppression (tanh soft-knee)            [optional]
  -> Delay estimation + micronavigation solver     [optional]
  -> GPU TDBP (CuPy kernel, n_tof_iters=1)
  -> SLC (complex, along-track x ground-range)
  -> DRC + SAS colormap  + PNG outputs

Preprocessing stages

Lowpass filter. create_lowpass_filter(BW/2, fs) (Hamming-windowed firwin, 65 taps). Run as zero-phase forward-backward on GPU (_fir_filtfilt_gpu). Zero-phase prevents group-delay distortion of the chirp; the cutoff drops out-of-band noise before MF so the matched filter operates on the chirp support only.

Matched filter. create_matched_filter(BW, T_p, fs, window='hanning') generates the time-reversed conjugated baseband LFM replica and applies a single Hanning taper (the replica itself is un-windowed; a prior bug double-windowed it and widened the main lobe to ~2x Rayleigh). Applied as batched FFT convolution on the GPU cube; the n_taps - 1 group delay is trimmed from the front.

Cross-ping TVG. After MF, compute tvg_raw = median(|data|, axis=(pings, channels)), producing a 1-D curve of shape (n_t,). Median is taken across all pings x channels rather than per-ping so that a bright target visible in a subset of pings cannot bias the normalizer upward (per-ping median produces a dark halo around objects; cross-ping does not). Divide the cube by this curve to flatten the spherical-spreading envelope.

The raw curve is smoothed by a boxcar moving average of width BeamformSettings.tvg_smoothing_taps (default 101; ~5% of the range gate) with edge-padding, producing tvg_pre. The divisor is tvg_pre + 1e-10. Taps of 0 or 1 disable smoothing.

When debug=True, the pipeline writes <output_dir>/<output_prefix>_tvg_debug.{npz,png} containing the raw median, the smoothed divisor, and the post-TVG residual (should be flat ~0 dB). Useful to verify TVG did not scrub real range-structure out of the scene.

Range upsampling (8x). Pulse compression is at fs = 75 kHz (10 mm range bins). The TDBP kernel interpolates between range samples linearly; at raw fs this leaves ~-7 dB sub-sample amplitude ripple because the chirp spectrum extends to 80% of Nyquist. FFT zero-pad upsampling by 8x pushes the interpolation error below -40 dB (clean). Controlled by environment variable BF_RANGE_UPSAMPLE (default 8). The effective fs passed to the BF kernel is 8x the HDF5 value.

Glint suppression (optional, on by default). tanh soft-knee above a configurable percentile (default 99th) of magnitude, preserving phase. Keeps a single ultra-bright return from dominating the SLC histogram and collapsing DRC contrast elsewhere.

Short inputs. apply_lowpass_filter caps scipy's filtfilt padding at the input length. The side-scan front end lowpasses the matched-filter replica itself (pulse_length x fs samples: 160 for a HISAS-preset scene), which scipy's default 3 x taps padding refused (2026-09-08).

Micronavigation (optional)

Enabled by use_micronav=True. Estimates sway (body-y) and heave (body-z) velocities per ping from the redundant-phase-center (RPC) principle: pings k and k+1 share a subset of phase-center locations, so overlapping channels should see identical echoes up to a micronav residual delay.

Steps: 1. Find overlapping channel pairs for all adjacent ping pairs (overlap.find_overlapping_channel_pairs). 2. Slide 128-sample windows (64-sample overlap) across the 20-90 m slant range gate (hard-coded in pipeline.py; stays fixed regardless of imaging range). 3. GPU batched cross-correlation with 8x upsample + phase refinement (estimate_delay_batch_gpu), producing coarse/fine delays and NCC. 4. Filter by NCC >= 0.66. 5. Least-squares solve for 2*N_PINGS unknowns (yvel, zvel per ping) with Huber loss, sqrt(NCC) weights, and jerk-continuity smoothness penalty (micronav_solver.compute_residuals_numba). 6. Re-run TDBP with micronav-estimated body velocities substituted for the POSSM ground-truth sway/heave.

POSSM truth is always kept for comparison (slc_possm); micronav SLC is an additional output (slc_micronav).

GPU TDBP kernel

beamformer.run_beamforming -> CuPy CUDA kernel BEAMFORM_KERNEL_TILED (see beamformer/sas/beamformer.py). Tile-based; each block owns a (tile_at x tile_r) image patch. Precomputes per-ping TX position, phase center, rotation matrix, beam-center angle, and a dense (n_pings, n_tof_samples, ...) trajectory for intra-ping motion compensation. Accumulates per-channel contributions coherently with linear range interpolation plus a synthetic-aperture beam-angle window.

n_tof_iters = 1 always. This is the motion-compensation sub-iteration count; the CLAUDE.md rule is non-negotiable for production imaging. n_tof_iters_override=0 is reserved for stop-and-hop apples-to-apples comparison against the FW simulator.

BeamformSettings (user-facing knobs)

Field Default Purpose
pixel_spacing 0.0125 m/pixel (isotropic)
along_track_length None Image length (m); None -> track extent + buffer
along_track_buffer 2.0 Buffer on each end when length is None
range_min, range_max 0.0, 100.0 Slant-range window (m)
chirp_bandwidth 60e3 Hz; HDF5 does not store this
pulse_length 4e-3 s; HDF5 does not store this
use_micronav False Enable RPC delay-est + LS solver
glint_suppression True tanh soft-knee above glint_percentile
glint_percentile 99.0 Threshold percentile
beam_width_safety 1.5 Multiplier on computed beamwidth
beam_window 'hann' SA window: hann/uniform/hamming/taylor[:sll[:nbar]]
first_channel 'auto' 'auto' drops leading channels per ping; or int
n_tof_iters_override None Keep None -> 1
tvg_enabled True Cross-ping TVG on/off
tvg_smoothing_taps 101 Boxcar smoothing of TVG curve
output_dir, output_prefix 'sim_output', 'sim'
debug False Dump TVG diagnostic npz+png
save_images, verbose True, True

Beamwidth is computed, never hardcoded: min(lambda/pixel_spacing, 0.886 * lambda/channel_spacing) * beam_width_safety.

Output artifacts

Written to <output_dir>/<output_prefix>_*:

  • *_drc_possm.png - Schlick-toned, SAS_COLORMAP. Primary visual output.
  • *_logmag_60dB_possm.png - 60 dB log-mag PNG (no tone mapping).
  • *_drc_micronav.png, *_logmag_60dB_micronav.png - if use_micronav.
  • *_tvg_debug.{npz,png} - if debug=True.

beamform_tdbp also returns an SlcImage whose save_port_image() / save_port_image_labeled() render the canonical starboard layout (track on left edge, along-track bottom-to-top, range left-to-right). See coordinate-systems.md for the layout rationale.

MSHDF export (SlcImage.save_mshdf / beamformer/sas/mshdf_writer.py) writes the focused complex SLC as a valid MCM Sensor HDF v2.1 file: a single /Sonar/Sensor1/Data1 with PingData shaped [1, n_along, n_range] of MSHDF_COMPLEXNUMBER, plus a full /Platform log and committed named datatypes. SAS-only; complex phase is preserved (unlike the magnitude-only XTF export).

Invariants

  • n_tof_iters = 1 (CLAUDE.md hard rule).
  • Beamwidth computed from wavelength and channel_spacing.
  • Initial position P_0 = [0, 0, -altitude[0]] (seafloor at Z=0).
  • Match-filtered data required for cross-correlation in micronav.
  • PNG outputs are not committed to git.