ApertureLab Open the appsign-in required

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.

Aperture normalisation. The TDBP kernel sums beam_weight * ch_weight * sample over every contributing ping and channel and, by default, does not divide by anything. The synthetic aperture grows linearly with range, so a far pixel integrates more pings than a near one and a diffuse seafloor gets steadily brighter across the swath. Measured on a bare-sand scene at 200 scatterers/m2 in the sample_images geometry (HISAS 1030, 200 m far range, 20 m altitude, 2.5 cm pixels) with cross-ping TVG enabled: mean |SLC| follows r^+0.48 over 30-190 m, a 2.82x spread; a v1.0 dataset scene gave ~20 at 10 m and ~49 at 199 m. The exponents below are specific to that configuration. TVG has already removed the spreading and per-ping footprint growth from the data, so with tvg_enabled=False the residual range law is different and aperture normalisation alone will not flatten the swath.

aperture_normalize=True divides each pixel by its own noise-equivalent aperture gain, sqrt(max(sum((beam_weight * ch_weight)^2), 1e-12)). That squared-weight sum is the kernel's second output and the second return value of run_beamforming. On the same scene it takes the profile to r^-0.02 (1.37x spread over 30-190 m), i.e. flat.

The square root matters: a point target gains the full coherent sum(w) across the aperture, but a diffuse (speckle) field gains only sqrt(sum(w^2)) because its per-ping contributions add incoherently. Dividing by sum(w) therefore over-corrects and inverts the ramp to r^-0.52 (measured). Normalising by the noise-equivalent gain is the unit-noise-gain convention: speckle level is held constant across the swath and a point target's contrast against it still improves as sqrt(N).

The default is off everywhere: the absolute image scale is pinned by the FW calibration gates (docs/fw_calibration_tests.md) and by the stored reference imagery, and normalising would move all of them. The sample-images dataset generator leaves it off as well (Isaac, 2026-09-21): the cross-ping TVG already owns range normalisation, and that is its role. When the pipeline streams the beamform in ping batches, the squared-weight sums are accumulated across batches and the division happens once after the last batch; per-batch division would normalise by a partial aperture, so run_beamforming rejects aperture_normalize=True together with return_gpu=True.

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.

The sector gate and the aperture window sit on the phase centre. For every ping, channel and pixel the kernel gates the pixel to the integrated sector (beam_width_deg_override x beam_width_safety, or the computed beamwidth) and weights it by beam_window (Hann by default) evaluated at the ground-plane angle between boresight and the pixel. That angle is measured from the channel's effective phase centre, the midpoint of the TX at transmit time and the RX at receive time (the bistatic-to-monostatic equivalence the micronavigation uses), because the synthetic aperture is laid down by phase centres. The kernel forms it as 0.5 * ((pixel - TX(0)) + (pixel - RX(t_rx))); boresight is still taken from the attitude at receive time. Before 2026-09-22 the angle was measured from the RX alone. On one element in stop-and-hop (TX = RX) that is the same point, and the image is bit-identical; on MUSCLE (36 channels, TX at the array centre, kept receivers 75 to 525 mm ahead of it, plus the travel during the echo) it centred the window ahead of the phase centres and skewed the flat-sand kx-ky angle curve: at 20 m range the misfit against the symmetric prediction fell from 1.93 to 0.93 dB rms with the fix (the model floor there is 0.90 dB), and the angle centroid from +0.34 to -0.03 deg. Point-target widths and PSLR moved by less than 0.1 %. The coarse per-ping culls (tile visibility and the 2x per-pixel early-out) use the middle channel's phase centre at transmit time (beam_apex_all) for the same reason; from the middle RX they clipped the trailing edge of the window. simulator/validation.py::along_track_psf_prediction models the gate at the phase centre to match. Test: beamformer/tests/test_tdbp_beam_window_phase_centre.py.

The slant-range carrier, and the baseband stage. The kernel multiplies every sample by exp(+j 2π f_c tof), so the raw sum carries exp(j 2 k_c r_slant) along range: a point target's response is a windowed sinc times that phase ramp, and a flat-bottom image's range spectrum sits at 2 k_c cos(θ_dep), not at baseband. Magnitude products are unaffected; anything that treats the complex SLC as a band-limited baseband signal (a 2-D spectrum, a sinc upsample by spectral zero-padding, a complex resample) is not, because the carrier aliases to a range-dependent position inside the pixel band and, where that band straddles Nyquist, an upsampled cut rings at the pixel period (seen at 35 m on the validation fixture, 2026-09-16).

Since 2026-09-16 the pipeline removes it as its last stage (BeamformSettings.baseband, default True; YAML beamform.baseband): sas/baseband.py computes, per pixel, the minimum two-way delay from the platform track to the pixel (the closest-approach ping, so heave, sway and curvature are all in it) and multiplies the finished image by exp(-j 2π f_c τ_ref). The reference depends on the pixel, not the ping, so it commutes with the aperture sum and with the streaming batches, and it is applied once after both are complete. Interferometric phase between two images on the same grid is unchanged (both lose the same reference). The convention is recorded as baseband in /SLC/tdbp attrs and SlcImage.metadata; baseband: false reproduces the historical output. simulator/validation.py::demodulate_slc is the nominal-geometry approximation for an SLC whose track is not to hand.

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
carrier_rotation_hz None Teaching/diagnostic only: the frequency of the kernel's exp(+j 2π f tof) rotation; None uses the file's fc (the only value that focuses). Wavelength, beam, micronav and the baseband stage keep the true fc. website/baseband.html uses 0, -fc and fc (1 + eps); test beamformer/tests/test_carrier_rotation.py
tvg_enabled True Cross-ping TVG on/off
tvg_smoothing_taps 101 Boxcar smoothing of TVG curve
aperture_normalize False divide each pixel by its noise-equivalent aperture gain, sqrt(sum (beam*channel weight)^2)
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.

DRC PNG format. Every DRC PNG written by save_slc_port_image (SlcImage.save_port_image / save_labels_image, and so the scene runner's <prefix>_<bf>_drc.png) with no watermark and no annotations is a palette PNG (PIL mode P): one byte per pixel holding the colormap index, with SAS_COLORMAP as the PNG palette. Resolved to RGB it is pixel-identical to the 3-channel colormapped image, at about a third of the file size. The scene YAML's output.drc_format values rgb (default) and indexed both write this PNG; geotiff writes the same indices to a georeferenced TIFF. When a watermark or annotation boxes (the _labels.png image) are drawn on top, the image is drawn in RGB and re-encoded losslessly: a palette PNG when it has at most 256 distinct colours, the RGB PNG otherwise (the watermark's antialiased logo and text usually add more than 256). Colormap-coloured pixels keep their colormap index in either palette case. A reader that needs RGB calls Image.open(path).convert("RGB"); np.asarray(Image.open(path)) on a palette PNG returns the (H, W) index array. The aperturelab tEXt metadata chunk is written on every path, and the scene runner's COCO record gives the configured format (drc.format) and what was written (drc.pixels: palette or rgb). The matplotlib figure (_port.png) and the log-magnitude PNG are unchanged.

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.
← Seafloor realism Boulder fields →