[ppp-one]

Add optional numba backend to find_transform

Summary

This PR makes plate solving much faster when many asterism pairs must be tested. A new backend argument on find_transform and compute_wcs selects a numba-compiled loop. numba is optional. Without it, twirl uses the current NumPy loop and gives the same results.

Why

Testing one candidate asterism pair costs about 25 µs. Almost all of that is NumPy call overhead on tiny arrays (lstsq on 3–4 points, a 20×40 distance array). The number of pairs grows as about n⁸ for quads and n⁶ for triangles. When min_match is not reached, find_transform tests every pair. With 20 detections and 40 catalog stars, one solve then takes 10–20 s.

On 357 real SPECULOOS images, run through the astra pointing pipeline, 114 images (32%) took more than 1 s with the current code, and the slowest took 19 s.

Changes

• find_transform(..., backend="auto") and compute_wcs(..., backend="auto"): "numba", "numpy", or "auto" (numba when it is installed). The default is "auto", so existing calls need no change.
• twirl/_numba.py: the compiled pair loop (least-squares affine fit, projection, cross-match count, early stop). numba is imported only when this backend is used, so import twirl does not get slower.
• Candidate pairs are built as NumPy arrays, not as a Python list. With quads the list could hold millions of pairs.
• quads.reorder uses np.take_along_axis, not a Python loop. The output is identical, and this is about 2× faster on 40 stars.
• pyproject.toml: new optional extra twirl[numba], and numba is added to the dev group so CI tests both backends. uv.lock is updated (no existing pins changed).
• README: a short note on pip install "twirl[numba]".
• benchmarks/: plate_solve.py times compute_wcs on synthetic fields and checks each solution against the true WCS. complexity.py plots how the cost grows with the number of stars. Both have a --backend option.

Results

Synthetic benchmark (benchmarks/plate_solve.py, median per case):

case	before	after	speed-up
crowded, 20 stars, quads	627 ms	30 ms	21×
crowded, 20 stars, triangles	391 ms	12 ms	32×
30 stars, quads	282 ms	160 ms	1.8×
10–15 stars	6–15 ms	6–9 ms	1.0–1.7×

When early stop happens, the fixed cost of compute_wcs (about 5 ms) dominates, so the gain is small there.

357 SPECULOOS images, 20 detections and 40 Gaia stars each:

	median	slowest
quads (default), before → after	558 ms → 218 ms	19.3 s → 0.61 s
triangles (asterism=3, min_match=0.5), before → after	95 ms → 23 ms	4.0 s → 0.26 s

The solutions are identical to the current code on all 357 images, for both asterisms.

Measured on an Apple M2 Max with 8 parallel workers.

Compatibility

• No breaking change. The new argument has a default, and results are the same with or without numba.
• numba stays optional, for these reasons: it adds about 150 MB (llvmlite), it caps NumPy (numba 0.68 requires numpy<2.6), and it usually supports a new Python version some months after its release.
• The first numba call in a process compiles the loop (about 0.6 s). numba then caches it on disk if the folder is writable.

Testing

• tests/test_backends.py (50 cases):
  • reorder matches the old loop exactly
  • both backends return the same transform on simulated fields (n = 6, 12, 25; triangles and quads; with and without noise; with and without early stop)
  • compute_wcs gives the same WCS with both backends
  • an unknown backend raises ValueError
• The full suite passes: 83 tests.
• Checked end to end on the 357 SPECULOOS images against the unchanged code, as above.

Notes

• I did not touch the ruff warnings that were already on main.