Affected Repositories

• PyAutoLens (primary)
• PyAutoArray (Phase 2a: autoarray/structures/triangles/shape.py only)

Branch Survey

worktree_check_conflict reports a read-only research claim on PyAutoArray main (numba-vs-jax-sparse, no branch, no worktree) and two heart-ack lines it misparses as claims; no worktree holds either repo. The numpy-deflections-p1 claim that held both repos this morning has been released (its PRs merged).

Suggested branch: feature/image-source-mappings-p2 (both repos, one worktree ~/Code/PyAutoLabs-wt/image-source-mappings-p2/)

Implementation Steps

Phase 2a — PyAutoArray (PR-A)

1. autoarray/structures/triangles/shape.py: add abstract Shape.contains(points) -> np.ndarray[bool] and Shape.boundary(n=100) -> (n, 2) closed loop. points uses the same axis order as triangles[..., 0] / [..., 1] — the order mask() already reads. Document the convention once on the Shape base class, including that PointSolver.solve builds Point(*source_plane_coordinate) from a (y, x) tuple, so the attribute named x holds the first (y) coordinate.
   • Point.contains raises NotImplementedError pointing at Circle / PointSolver; Point.boundary → single row.
   • Circle: radial test; parametric boundary.
   • Triangle: barycentric test against the points (the formula already in triangle_contains_mask), with a test that pins the vertex axis order — triangle_contains_mask unpacks vertices as y1, x1 = self.a while the centroid it tests is (idx0, idx1), which is a reflected containment test unless vertices are (y, x); resolve by test, fix if wrong.
   • Polygon: any-of its fan triangles (exact for convex polygons — say so in the docstring); boundary = vertices + first repeated.
   • Square: bounds with the top / bottom numeric ordering made explicit (mask assumes top < bottom); boundary = four corners closed.
2. New test_autoarray/structures/triangles/test_shape.py: per shape, contains on hand-picked inside/outside points; boundary is closed and its edge midpoints test as inside (or on the edge, by a tolerance); Circle.contains agrees with Circle.mask on a triangle whose centroid is the test point.
3. pytest test_autoarray/structures/triangles -q green; open PR-A with label pending-release.

Phase 2 — PyAutoLens (PR-B), validation before features
4. test_autolens/point/triangles/test_shape_solver.py (beside the existing test_solver*.py, test_extended.py, test_use_jax.py; fixtures in that directory's conftest.py): (a) SIS + small Circle on-axis → ring; off-axis inside θ_E → two images whose summed find_magnification matches the analytic SIS |θ/(θ−θ_E)| to ~1% as radius → 0; (b) solver regions vs the grid-trace oracle (tracer.traced_grid_2d_list_from → shape.contains → slim bool → aa.inversion.mappings.image_regions_from_slim_mask) agree within one pixel for a source spanning ≥5 pixels; (c) steps() kept area converges monotonically; (d) Triangle / Polygon / Square behave like Circle; (e) plane_redshift on a three-plane tracer; (f) use_jax=True and numpy give identical kept triangles. Revive the commented-out test in test_extended.py.
5. autolens/point/solver/shape_solver.py: fix what step 4 exposes, each with a regression test. Known audit items (2026-08-19, draft/feature/autolens/area_magnification_leggos.md): find_magnification hardcodes xp=np and ignores self._xp; kept-triangle area over-covers by up to one triangle edge per boundary with no convergence statement; JAX triangle area cannot survive jax.jit; no per-image split. Defects that live in autoarray/structures/triangles/* go into PR-A instead.
6. ShapeSolver.image_regions_from(tracer, shape, plane_redshift=None, xp=None) -> List[aa.ImageRegion]: solve_triangles → connected components over shared vertices → per component the triangle polygons (arcsec, (y, x)) and pixel coverage on the solver's grid (slim_indexes + mask), built on the Phase 1 ImageRegion so brightest_coordinate_from / centroid_from / flux_from are shared. ShapeSolver.mapping_from(...) -> aa.Mapping with source_contours = [shape.boundary()]. find_magnification gains a per-image split from the same grouping.
7. New autolens/lens/mappings.py: source_mapping_from(tracer, grid, shape, plane_index=-1, plane_redshift=None, pixel_scale_precision=None); traced_region_from(tracer, region_grid, plane_index=-1) (source-plane convex hull); mappings_from_fit(fit, plane_index=-1, shape=None, **clump_kwargs) — engine A when tracer.planes[plane_index].has(cls=aa.Pixelization) (the predicate at fit_imaging_plots.py:106), else engine B with Circle at the first light profile centre and half_light_radius; multiple_image_positions_from(fit, plane_index=-1, use_centroid=False) -> aa.Grid2DIrregular from fit.model_images_of_planes_list[plane_index]; multiple_image_pixel_coordinates_from(...) via Geometry2D.pixel_coordinates_wcs_2d_from; magnifications_from(fit, plane_index=-1). Numpy-only, never jitted. No RA/Dec.
8. autolens/imaging/plot/fit_imaging_plots.py::subplot_mappings(fit, plane_index=-1, ...): 2×2 — lens-subtracted data + regions + critical curves (_compute_critical_curves_from_fit); model image + regions; source plane via _plot_source_plane + clumps/shape + caustics; unzoomed source. Export from autolens/plot/__init__.py (the fit_imaging_plots import block); wire the existing plots.yaml: subplot_mappings key in autolens/imaging/model/visualizer.py.
9. Tests: test_autolens/lens/test_mappings.py (Isothermal + Circle inside the caustic → ≥2 regions; outside → 1; parametric vs pixelized dispatch); subplot_mappings on fit_imaging_x2_plane_7x7 and fit_imaging_x2_plane_inversion_7x7 in the fit-imaging plot tests.
10. Verification: pytest test_autolens/lens/test_mappings.py test_autolens/point/triangles test_autolens/imaging/plot -q; render a parametric fit's subplot_mappings and eyeball it; multiple_image_positions_from on a simulated quad within a pixel of PointSolver positions on a full-res dataset (never under PYAUTO_SMALL_DATASETS, whose PointSolver.solve short-circuit returns a fixed pair). PR-B body carries ## API Changes and a ## ShapeSolver audit section (each finding: fixed / documented / deferred).

Key Files

• PyAutoArray/autoarray/structures/triangles/shape.py — Phase 2a contains / boundary
• PyAutoArray/test_autoarray/structures/triangles/test_shape.py — new
• PyAutoLens/autolens/point/solver/shape_solver.py — audit fixes, image_regions_from, mapping_from, per-image magnification
• PyAutoLens/test_autolens/point/triangles/test_shape_solver.py — validation suite (new); test_extended.py revived
• PyAutoLens/autolens/lens/mappings.py — new fit-level API
• PyAutoLens/autolens/imaging/plot/fit_imaging_plots.py — subplot_mappings
• PyAutoLens/autolens/plot/__init__.py, PyAutoLens/autolens/imaging/model/visualizer.py — export + plots.yaml wiring
• PyAutoLens/test_autolens/lens/test_mappings.py — new

ShapeSolver-based source-to-image region mappings, fit-level subplot_mappings, brightest multiple-image positions, and the ShapeSolver validation suite in PyAutoLens

Type: feature
Target: PyAutoLens
Repos:

• PyAutoLens
• PyAutoArray
  Difficulty: large
  Autonomy: supervised
  Priority: high
  Status: formalised
  Consequence: judge
  Review-minutes: 25
  Unattended: needs-slicing
  Epic: image-source-mappings
  Phase: 2
  Filed: 2026-09-02

ShapeSolver-based source-to-image region mappings, fit-level subplot_mappings, brightest multiple-image positions, and the ShapeSolver validation suite in PyAutoLens.

Phase 2 of the image-source-mappings epic — ledger
draft/feature/autoarray/image_source_mappings_epic.md. Depends on Phase 1
(complete/2026/09/image-source-mappings-p1.md — merged 2026-09-02 as PyAutoArray#517, not yet released)
being merged AND released: this phase imports
Mapping / ImageRegion and the regions= overlay from the released PyAutoArray. Phase 3 (workspace + tutorials)
opens only after this one releases.

Engine B — ShapeSolver as the source→image engine

User decision (mid-plan revision, 2026-09-02): for non-pixelized sources, reuse as much ShapeSolver code as
possible so the same code serves the shape-tracing visuals AND the precise numerical calculations (magnification,
source size) — and this phase becomes the validation ShapeSolver never had. The simpler grid-trace membership approach
is kept only as a test oracle (pixel-resolution cross-check), never as API.

ShapeSolver.for_grid(grid, pixel_scale_precision, ...) → solve_triangles(tracer, shape, plane_redshift) returns the
kept image-plane triangles whose traced images contain the source Shape. From them: (i) group triangles into connected
components by shared vertices (the separate multiple images), (ii) each component's polygons = its triangles (drawn as
filled triangles, union outline optional), (iii) pixel coverage = image pixels containing any kept-triangle vertex or
centroid → slim_indexes/mask on the same ImageRegion object as engine A, so brightest_coordinate_from etc. work
identically. Sub-pixel sources still produce regions because the triangles carry the geometry, not the pixel grid.
steps() exposes every refinement level for the guide's "how the solver homes in" figure (Phase 3).
find_magnification gives the per-image / total magnification numerically. Default shape for a parametric source:
Circle at the profile centre, radius = half_light_radius (autogalaxy/profiles/light/abstract.py:150).
PointSolver.solve stays the point-level overlay of the centre. Test oracle only: trace the image grid,
shape.contains(points) → slim bool → image_regions_from_slim_mask; the ShapeSolver regions must agree with it to
within a pixel for sources larger than a few pixels.

Verified facts (do not re-derive)

3. Shape classes (autoarray/structures/triangles/shape.py) only have mask(triangles) (which reads (x, y) order);
   there is no point-containment method.
4. Array2D.from_fits retains the FITS header (uniform_2d.py:900), so dataset.data.header carries WCS keywords when the
   source FITS had them — enough for the guide-level astropy recipe (Phase 3).
5. PointSolver.solve short-circuits under PYAUTO_SMALL_DATASETS to a fixed pair — never assert on its values in smoke.

(Numbering matches the ledger's list; facts 1, 2, 4 and 6 are Phase 1's.)

Cross-reference — the existing ShapeSolver audit

draft/feature/autolens/area_magnification_leggos.md (cluster-strong-lensing phase 6) carries the ShapeSolver audit
verdict of 2026-08-19: effectively unmaintained — one method (find_magnification = kept-triangle image area /
source shape.area), its only test commented out since ~2024-09, use_jax=True silently ignored (find_magnification
hardcodes xp=np and never consults self._xp), the JAX triangle area cannot survive jax.jit (__len__ returns a
traced count_nonzero), no per-image split, no magnification-threshold filter, no error estimate, and image area
over-covers by up to one triangle edge per boundary with no convergence check; sole consumer is one unused demo in the
workspace point_source/fit.py. That prompt frames it as "rehabilitate or retire" against a LEGGOS per-pixel-inversion
primary API. Do not merge the two prompts. This phase settles the rehabilitate half from the mappings side — the
per-image split, the xp dispatch, the revived tests and the convergence statement are all listed below — and the
LEGGOS prompt keeps its own scope (area magnification of an arc via the μ-map). Read that audit before touching
shape_solver.py and reconcile any finding that has since changed.

Files

• PyAutoArray/autoarray/structures/triangles/shape.py — add Shape.contains(points_yx) -> bool array (Circle radial,
  Triangle barycentric, Polygon even-odd, Square bounds, Point raises pointing at Circle/PointSolver) and
  Shape.boundary(n=100) -> (N,2) yx. Takes (y, x); the docstring notes the legacy (x, y) order of
  mask(triangles). Small PyAutoArray change — it did not ship in the Phase 1 PR (feat: redesign subplot_fit panels — add Source Plane (Mid Zoom) #517 leaves shape.py
  untouched), so it is the tiny Phase 2a PyAutoArray PR, and it must merge and release before the
  PyAutoLens PR opens.
• autolens/point/solver/shape_solver.py — ShapeSolver audit + reuse: read solve_triangles, steps,
  find_magnification, the triangle containers (autoarray/structures/triangles/{array,coordinate_array}*.py:
  containing_indices, neighborhood, up_sample, for_indexes) and the PYAUTO_SMALL_DATASETS / JAX paths; fix
  what the validation below exposes (bugs fixed in place, each with a regression test — never worked around in the
  mappings layer). Add ShapeSolver.image_regions_from(tracer, shape, plane_redshift=None) -> List[ImageRegion]
  (triangle grouping by shared vertices → per-image triangle polygons + pixel coverage on the solver's grid mask) and
  ShapeSolver.mapping_from(tracer, shape, ...) -> Mapping (source contour = shape.boundary()).
• new autolens/lens/mappings.py:
  • source_mapping_from(tracer, grid, shape, plane_index=-1, plane_redshift=None, pixel_scale_precision=None) -> Mapping
    — builds ShapeSolver.for_grid and delegates (thin).
  • traced_region_from(tracer, region_grid, plane_index=-1) -> np.ndarray (source-plane convex hull of a traced region).
  • mappings_from_fit(fit, plane_index=-1, shape=None, **clump_kwargs) -> List[Mapping] — engine A when the plane
    has(cls=aa.Pixelization) (predicate at fit_imaging_plots.py:~105), engine B (ShapeSolver) otherwise (Circle at
    the first light profile centre, half_light_radius).
  • multiple_image_positions_from(fit, plane_index=-1, use_centroid=False) -> aa.Grid2DIrregular — per image region,
    brightest pixel (or flux-weighted centroid) of fit.model_images_of_planes_list[plane_index]; plus
    multiple_image_pixel_coordinates_from(...) via Geometry2D.pixel_coordinates_wcs_2d_from.
    No RA/Dec in the library — WCS/RA-Dec conversion is guide-level only (user decision, 2026-09-02): the library
    returns arcsec + pixel coordinates and the Phase 3 guide shows the astropy conversion inline.
  • magnifications_from(fit, plane_index=-1) — per-image flux over source-plane flux.
• autolens/imaging/plot/fit_imaging_plots.py::subplot_mappings(fit, plane_index=-1, ...) — lens-subtracted data +
  regions + critical curves (_compute_critical_curves_from_fit, :25); model image + regions; source plane via
  _plot_source_plane (:68) + clumps/shape + caustics; unzoomed source. Export from autolens/plot/__init__.py; wire
  plots.yaml subplot_mappings into the imaging visualizer (autolens/imaging/model/visualizer.py) so the key goes live.

Tests

• test_autolens/lens/test_mappings.py — Isothermal + Circle inside the caustic → ≥2 regions; outside → 1; parametric
  vs pixelized dispatch.
• test_fit_imaging_plots.py — subplot_mappings on fit_imaging_x2_plane_7x7 and the inversion fixture.
• test_autoarray/structures/triangles/test_shape.py — contains / boundary per shape.

ShapeSolver validation suite (test_autolens/point/solver/test_shape_solver.py, extending what exists):

(a) SIS lens, small Circle source on-axis → Einstein ring; off-axis inside the Einstein radius → two images whose
find_magnification matches the analytic SIS magnification |θ/(θ−θ_E)| summed over images to ~1% as the circle
radius → 0; (b) ShapeSolver regions vs the grid-trace oracle agree to within one pixel for a source spanning ≥5 pixels;
(c) steps() refinement is monotone (kept area converges); (d) Triangle/Polygon/Square shapes behave like Circle;
(e) plane_redshift on a three-plane tracer; (f) the JAX (use_jax=True) and numpy paths give identical kept triangles.

Failures here are ShapeSolver bugs to fix in this phase, reported in the PR body under a "ShapeSolver audit" heading.

Verification

pytest test_autolens/lens/test_mappings.py test_autolens/point/solver test_autolens/imaging/plot -q; a parametric fit
produces the same figure with a Circle in the source plane and its ShapeSolver triangles filled in the image plane;
multiple_image_positions_from on a simulated quad returns 4 positions within a pixel of the PointSolver positions of
the centre (full-res dataset, not smoke — verified fact 7: PointSolver.solve short-circuits under
PYAUTO_SMALL_DATASETS); the ShapeSolver validation suite passes and its audit findings are listed in the PR body.