A microphone array direction-of-arrival pipeline in Python: array geometry and beam patterns, a physically-grounded simulation chain, three DOA estimators benchmarked against each other, and a real-time streaming localizer.
Everything runs on synthetic arrays and synthetic audio out of the box. There is no proprietary geometry, no vendor SDK and no hardware dependency — the array is defined by a list of coordinates, so a different array is a one-line change and a real capture is a WAV file with the right channel order.
pip install -r requirements.txt
python examples/quickstart.py| Module | Role |
|---|---|
doa/arrays.py |
Capsule coordinates, the far-field delay model, analytic beam patterns, six hypothetical array presets |
doa/signal_chain.py |
Source generation, fractional delay, array simulation, delay-and-sum, and the two slow reference localizers |
doa/fast_localize.py |
GCC-accelerated localizer with a precomputed lag table, multi-source peak picking, and the real-time stream |
doa/evaluate.py |
Wraparound-safe error metrics, grid benchmarking, operating-region reporting, plots |
Scripts wrap each of these as a CLI: make_dataset.py, run_benchmark.py,
plot_results.py, live_doa.py.
| Method | Domain | Complexity | Reverberation | Noise | Real-time |
|---|---|---|---|---|---|
| SRP | time | O(G·M·L) | weak | good | no |
| SRP-PHAT | frequency (phase) | O(M²·L·log L + G·M²·L) | strong | weak | no |
| GCC-accelerated | frequency → lag lookup | O(M²·L·log L + M²·G) | strong | weak | yes |
SRP steers the beamformer to every candidate angle and measures output power. Simple and noise-tolerant, but a strong reflection contributes real power at its own arrival direction, so it degrades in a live room.
SRP-PHAT normalises each pair's cross-spectrum to unit magnitude so only phase survives. That is what makes it reverb-robust: a loud reflection no longer outweighs the direct path just by being loud. The trade is the opposite of what people usually assume — whitening amplifies bins containing nothing but noise, so PHAT is reverb-robust, not noise-robust.
GCC-accelerated is the same computation as SRP-PHAT viewed from the other side of the correlation theorem. Instead of rebuilding the steering response at every angle, cross-correlate each capsule pair once and then read each angle's score directly out of those curves at the lag that angle implies. The per-angle term collapses from a complex exponential per pair to a single array lookup, which is the entire speedup. Precomputing the pair→angle→lag table at startup is what makes it real-time rather than merely fast.
The estimates agree to within about a degree, but not bit-for-bit, even though the math says they should be equivalent. The residual is most likely how the GCC curve is sampled between lag bins, but that hasn't been run down.
Delay sign. geom.delays(az, el) returns tau_i = (p_i · u)/c, where u
points from the array toward the source. A capsule with a larger projection
onto u hears the wavefront earlier, so a positive tau_i is an advance and
the simulation applies -tau_i. Every localizer must predict delays with the
same sign.
Get this wrong and the failure is distinctive: a uniform sign flip on all
delays is exactly u → -u, so the estimate lands precisely 180° from the
truth rather than becoming noisy. In the cross-spectrum this shows up as
conj(X_i)·X_j versus X_i·conj(X_j). A clean antipodal error is always a
convention bug, never noise and never a real ambiguity. test_doa.py has an
explicit guard for it.
Which ambiguities are real. A planar array resolves the full 360° of
azimuth uniquely — there is no front/back ambiguity to explain away. The one
inherent ambiguity a planar array has is the elevation mirror: a source
above the plane and its reflection below produce identical delays. The
ula8 preset is included as the contrasting case, since a linear array does
have a genuine front/back azimuth ambiguity, and tetra4 is the 3D geometry
that has neither.
Fractional delay is not optional. At 48 kHz one sample is ~20.8 µs, or ~7.1 mm of acoustic path — comparable to the inter-capsule delays a compact array actually sees. Rounding to whole samples quantises away most of the geometry the beamformer depends on, so delays are applied as a frequency-domain phase shift, exact for any τ.
Low-frequency omnidirectionality is not spatial aliasing. An array becomes
non-directional at low frequency because its aperture is small relative to the
wavelength. Grating lobes appear at high frequency when capsule spacing
exceeds λ/2. geom.aliasing_freq() reports the second. They are separate
phenomena and get conflated constantly.
Azimuth accuracy on a planar array degrades as a source approaches boresight (straight up, out of the array plane), because the inter-capsule delays all collapse toward zero and azimuth stops being observable. That is geometry, not an estimator defect, and it does not improve with a better algorithm.
Quoting one accuracy figure averaged across all elevations therefore
understates a method that is fine everywhere it can physically work. The
benchmark instead reports inside and outside a stated operating region — the
default boundary is EL ≤ 66°, and plot_error_vs_elevation is where that
number comes from. On the default 16-capsule dual-ring array at 20 dB SNR:
method region n mae_deg rmse_deg p95_deg within_1deg_pct
GCC in-region 88 0.024 0.035 0.088 100.000
GCC out-of-region 32 28.162 58.450 134.997 75.000
RMSE is the conventional DOA reporting metric because it penalises the occasional gross miss that a mean absolute error hides. The "within N degrees" fractions are worth carrying alongside it for a mixed audience.
from doa import arrays, signal_chain as sc, fast_localize as fl
geom = arrays.load_preset("dual16") # or ula8, uca8, grid16, spiral16, tetra4
t = sc.time_axis(0.5)
src = sc.speech_like(t)
sigs = sc.mic_signals(src, az=137.5, el=15.0, geom=geom,
snr_db=20, reflections=4)
g = fl.Geometry(geom, fs=48000) # precompute the lag table once
est, curve, angles = g.localize(sigs)
peaks = g.localize_multi(sigs, n_sources=2)Custom geometry — pass coordinates in metres, that's the whole interface:
geom = arrays.ArrayGeometry([[0.04, 0, 0], [-0.04, 0, 0],
[0, 0.04, 0], [0, -0.04, 0]], name="quad")Benchmark and plots:
python scripts/run_benchmark.py --methods gcc srp srp_phat --out results/benchmark.csv
python scripts/plot_results.py results/benchmark.csv --outdir results/figsLabelled test set on disk, then score it — this is the path real recordings take, since the filenames carry the ground truth:
python scripts/make_dataset.py --out data/synthetic --reflections 6 --snr 10
python scripts/run_benchmark.py --data-dir data/synthetic --methods gccLive:
python scripts/live_doa.py --list
python scripts/live_doa.py --device 7 --duration 30 --log run.csv--log writes per-hop estimates to CSV and the raw multichannel audio to a
.wav of the same name, sharing hop boundaries, so any logged angle traces
back to the audio that produced it at sample index round(elapsed_s * fs).
The array is just coordinates and the audio is just a WAV, so real captures need three things to line up:
- Channel order must match the row order of
geom.positions. A permuted map produces a stable, confident, wrong answer. - Sample alignment across channels. If your interface enumerates as
several smaller endpoints rather than one multichannel device, they are not
guaranteed sample-locked. On Windows an ASIO driver can merge them; that
path needs
sounddeviceinstalled via pip (not conda) andSD_ENABLE_ASIOset before the first import of sounddevice anywhere in the process.fast_localizesets it at module load, which is the only place it works. - A known-angle calibration capture. Put a source at a measured bearing and confirm the estimate lands there before trusting anything else. This single check catches both a permuted channel map and a flipped delay sign.
pytest -q # 42 tests- Extend the search to a 2D azimuth × elevation grid; the lag table already supports a fixed elevation, so this is mostly a shape change.
- Multi-source:
localize_multidoes greedy peak picking, which needs no extra computation and does not break on coherent sources. MUSIC is the heavier alternative and carries real constraints — it needs the source count known in advance, assumes narrowband signals, and breaks on coherent or reverberant sources without spatial smoothing that is awkward on a non-uniform 2D array. Exhaust peak picking first. - Move the DSP out of the audio callback into a worker thread if you ever see over/underflow status messages.
- Track the sub-degree GCC / SRP-PHAT discrepancy noted above.
MIT. See LICENSE.