docs: flagship 'PyAutoLens for Lenstool users' example — SMACS J0723 vs the RELICS Lenstool model

[Jammy2211]

Overview

The flagship adoption example: take the first JWST cluster, SMACS J0723.3−7327, and do a full Lenstool-style analysis in PyAutoLens — same dPIE profiles (via the merged from_lenstool / dPIEMassLenstool API, PyAutoGalaxy#487), same reference-anchored member scaling relation (autolens_workspace#238), same positional constraints — then compare directly against the published RELICS Lenstool model products on MAST. Written for a reader whose mental model is a Lenstool .par file. A real prospective user is available for beta iteration once a draft exists.

Autonomy: supervised --auto, in-session authorization (user directive 2026-07-09: "merge PRs and then do 8_lenstool_user_example").

Reference model (research complete, from Sharon et al. RELICS paper + MAST probe)

• RELICS Lenstool HST model (arXiv:2208.08483; MAST doi:10.17909/T9SP45): 1 cluster-scale PIEMD (7 free params, r_cut fixed 1000 kpc) + BCG PIEMD (pos/e/θ fixed from F814W SExtractor; r_core/r_cut/σ free) + 143 red-sequence members on the Limousin+05 scaling relation (mag0 = 19.12, r_core* = 0.15 kpc fixed, σ*/r_cut* free). 8 systems / 25 images (Table 1 extracted — positions in hand), photo-z (system 1 fixed z=2.2, others free U(0.5,5)). Image-plane RMS 0.58″. z_lens = 0.388, Ω_m=0.3/Ω_Λ=0.7/h=0.7 (the Lenstool-typical cosmology — exactly the documented from_lenstool cosmology knob).
• Public products confirmed on MAST (archive.stsci.edu/hlsps/relics/smacs0723m73/): best-fit κ/γ/ψ/deflection maps + μ(z=1,2,4,9) + 100 MCMC map samples (models/lenstool/, range/); ACS+WFC3IR photometric catalogues with BPZ photo-z (catalogs/); 30/60mas drizzled images (images/). No .par/params file published → the comparison target is the maps + image-plane RMS, which is the honest "same analysis, compare outputs" framing anyway.
• Mahler et al. 2023 (arXiv:2207.07101, Lenstool + MUSE) provides spectroscopic redshifts for several systems — to be extracted in Phase 1 and used where available (documented against the RELICS photo-z choices).

Plan (phased per the Brain routing)

Phase 1 — data prep (scripts/cluster/lenstool/data.py or equivalent): download the F606W/F814W 60mas cutout, the RELICS photometric catalogue, and the Lenstool κ + deflection maps from MAST (runtime download with cached local copies; nothing heavy committed). Build point_datasets.csv from the paper's Table 1 (25 images, 8 systems; spec-z where Mahler+23 provides, else RELICS photo-z treatment) and scaling_galaxies.csv from red-sequence member selection on the catalogue (F606W−F814W CMD, magnitudes → relative luminosities vs mag0=19.12).

Phase 2 — the example (scripts/cluster/lenstool/modeling.py + README): compose the RELICS-equivalent model — cluster-scale dPIEMassLenstool (fitting directly in σ/r_core/r_cut, priors in Lenstool units), BCG dPIEMassLenstool with fixed geometry, member tier on the reference-anchored relation with r_cut scaling — fit with source-plane χ² (Lenstool's default), Lenstool-typical cosmology. Prose maps every concept to its .par equivalent (potentiel section ↔ model component; v_disp/ellipticite/angle_pos ↔ constructor args incl. the √(3/2) σ trap; image catalogue ↔ point_datasets.csv; source-plane vs image-plane χ²). Tutorial prose = judgment tier.

Phase 3 — comparison + wrap (scripts/cluster/lenstool/comparison.py): our max-likelihood tracer's convergence + deflection fields vs the MAST Lenstool maps (pixel-level, with the 100-sample range as the tolerance band); image-plane RMS vs the published 0.58″; per-plane critical curves via the merged cluster plotters (PyAutoLens#578). Honest section on convention differences and where parity is/isn't expected (LTM-free, photo-z freedom, member-selection differences). Then beta-user iteration.

Affected Repositories

• autolens_workspace (primary; worktree ~/Code/PyAutoLabs-wt/lenstool-example/, branch feature/lenstool-example)

Original Prompt

Click to expand starting prompt

$(cat PyAutoMind/docs/cluster/8_lenstool_users_example.md)

[Jammy2211]

Phase 1 research complete — reference model upgraded from RELICS to Mahler et al. 2023 (arXiv:2207.07101), whose public repo (github.com/guillaumemahler/SMACS0723-mahler2022, ICLv2/) ships the complete Lenstool workflow: best.par (149 potentials: 5 named halos + 144 members, every one with x/y offsets from ref (110.826989, −73.454723), ellipticity, angle_pos, core/cut radius in arcsec AND kpc, v_disp — i.e. exact from_lenstool inputs), input.par (potfile: type 81, corekpc 0.15 fixed, mag0 19.12, σ* U(50,500), r_cut* U(1,100) kpc, vdslope 4 ⇒ σ∝L^¼, slope 4 ⇒ r_cut∝L^½ — identical to our shipped scaling convention), arcs.dat (21 systems incl. spec-z 1.449/1.3779/1.9914/1.425/5.17 and the famous z≈11.8 model source), galcat.cat (144 members with mags), bayes.dat (chain). Cosmology FlatLCDM H0=70/Ωm=0.3. Published image-plane RMS 0.32″ (χ²/ν=1.0, 46 free params).

This makes the example strictly better than planned: (a) reconstruct the published best-fit verbatim via a small .par parser + from_lenstool and verify PyAutoLens forward-predicts the observed images at the published RMS — a true cross-code validation; (b) refit in PyAutoLens with dPIEMassLenstool priors mirroring input.par and compare parameter posteriors table-to-table. RELICS MAST maps kept as an optional appendix. Files are runtime-downloaded with attribution (repo is GPL-licensed; nothing committed).

Now implementing: scripts/cluster/lenstool/ — data.py (download+parse → CSVs), modeling.py (reconstruct → verify → refit), README.md (the '.par ↔ PyAutoLens' mapping for Lenstool users).

[Jammy2211]

Shipped to PR-open: #240 (pending-release ✓).

The headline result — cross-code parity on a published model. Reading all 149 dPIE potentials of Mahler et al.'s best.par straight into al.mp.dPIEMass.from_lenstool and ray-tracing the 60 observed multiple images through the 22-plane tracer collapses every system to a compact source-plane group:

median source-plane rms = 0.068"   (all 21 systems ≤ 0.29")
point_16 z=1.097: 0.015"  ...  point_13 z=2.996: 0.270"

consistent with the published image-plane RMS of 0.32" through typical magnifications. This validates the full convention chain in one number: .par parsing, the West-positive relative frame, emass→ell_comps, the √(3/2) fiducal-σ convention, kpc↔arcsec at the run's own cosmology (H0=70/Ωm=0.3), and multi-plane β factors.

A real trap found and documented: PyAutoLens multi-plane tracers normalize profile deflections to the tracer's final plane — from_lenstool must take the highest source redshift (here z=11.76) as redshift_source. Normalizing to z=2 mis-scales every deflection by ~15% and smears the source groups to ~2" (we measured exactly that before fixing it). Now prominent in modeling.py + README.

Refit composition validated: 72 free parameters (30 mass — matching input.par's limits verbatim in Lenstool units via dPIEMassLenstool — + 42 source positions; the 16 photometric source-z's held at published values). The potfile member tier is two shared priors (σ* U(50,500), r_cut* U(1,100) kpc) with the fixed L^0.25/L^0.5 exponents — identical to the convention shipped in #238. Full search is a production job, gated (LENSTOOL_EXAMPLE_RUN_FIT=1), as are the forward-solve RMS check (RUN_IMAGE_PLANE) and the per-plane critical-curve figure (MAKE_FIGURES; ~10 min/plane through 149 profiles in numpy — a JAX-side speedup is a natural follow-up).

Next: merge when happy with the README voice — then this is the artifact to put in front of the beta tester for the back-and-forth iteration loop.

[Jammy2211]

PR #240 merged — the 'PyAutoLens for Lenstool users' example is live on main. The LensTool-parity thread (tasks 6 → 7 → 8, PyAutoGalaxy#485 / autolens_workspace#237 / here) is complete: Lenstool-native construction + fitting API, community-convention scaling relations, and a published-model walkthrough validating the whole chain at 0.068" source-plane parity. Next: beta-tester iteration on this example (tracked in the cluster tracker), and the remaining home-straight items 5_profiling + 10_solver_over_under_prediction. Closing.