Skip to content

HSM: measure PSF/star shapes directly in sky coords (use_sky_coords) - #812

Open
cailmdaley wants to merge 6 commits into
developfrom
feat/hsm-sky-coords
Open

HSM: measure PSF/star shapes directly in sky coords (use_sky_coords)#812
cailmdaley wants to merge 6 commits into
developfrom
feat/hsm-sky-coords

Conversation

@cailmdaley

@cailmdaley cailmdaley commented Jul 12, 2026

Copy link
Copy Markdown
Contributor

Closes #775.

ShapePipe measured PSF/star HSM adaptive moments in the pixel frame, then rotated the ellipticities into sky coordinates in a separate post-processing step (scripts/python/convert_psf_pix2world.py, a WCS-Jacobian decomposition). galsim's FindAdaptiveMom now does this at the point of measurement via use_sky_coords=True (available since galsim 2.5; we pin ≥2.8) — per @aguinot's suggestion on #741. This PR moves the WCS knowledge into the HSM call and drops the redundant rotation.

What changed

  • Direct sky-frame measurement. Every PSF/star FindAdaptiveMom call site — PSFEx single-exposure and multi-epoch (psfex_interp.py), MCCD interpolation (mccd_interpolation_script.py, shapepipe_auxiliary_mccd.py) — now attaches the local WCS (linearised at the object's image position) to the stamp and passes use_sky_coords=True. Multi-epoch paths use the per-CCD WCS log already in scope; single-exposure paths read the WCS from the galcat's LDAC_IMHEAD header.
  • Rotation removed where redundant; script renamed. convert_psf_pix2world.py no longer rotates the in-repo E1/E2/SIGMA_*_HSM shapes. Per review, it is renamed to collate_star_cat.py to reflect its remaining purpose: collate the per-exposure PSF validation catalogues into star catalogues (positions, MCCD CCD id, merge).
  • Patch-less v2.0 support. The script takes -V (default v2.0), following the get_ccds_with_psf convention from Coverage nexp masks #797. the patch-less default reads from the single run root and writes validation_psf_conv-<idx>.fits (no patch token; downstream combine_runs.bash / config_Ms_psfex_conv.ini globs match both schemes). -V v1.6 -P … selects the unchanged v1.x behaviour. Invalid versions raise. The collate section of pipeline_canfar.md shows both invocations.
  • Exception kept deliberately: the MCCD fit-validation PSF_MOM_LIST/STAR_MOM_LIST columns come from the external mccd library, which still measures in the pixel frame — that one branch keeps its rotation (with a docstring caveat) until mccd itself adopts use_sky_coords.
  • Equivalence proven in-tree. tests/module/test_hsm_sky_coords.py: a synthetic Gaussian stamp under a nontrivial local WCS (0.187″/pix, 28° rotation, shear), measured pixel-frame + old transform_shape math vs use_sky_coords=True — g1/g2/sigma agree to rtol 1e-6. Plus test_collate_star_cat.py pinning the v1.x/v2.0 path and filename logic. All 14 tests pass in the develop container.

Notes for review

  • Sigma unit change (intended): moments_sigma now comes back in world units (arcsec) rather than pixels, so HSM_T_* = 2σ² is arcsec² — which is what the old pipeline produced after the rotation step, so end-of-pipeline catalogues are unchanged. Downstream consumers (merge_starcat, ngmix T-guess) checked: no pixel-scale assumption.
  • No double rotation: grep-verified nothing downstream applies a second Jacobian/getDecomposition to the HSM ellipticities; merge_starcat passes them through.
  • A subtle off-by-one (galsim's FITS 1-indexing vs astropy origin-0 positions) in the local-WCS evaluation point was caught in internal review and fixed (f153ea28).

An independent adversarial review of the full diff ran before opening this PR; both findings it raised (the external-mccd branch above, and the indexing fix) are addressed.

— Claude on behalf of Cail

🤖 Generated with Claude Code

https://claude.ai/code/session_01Gao7Pi5QDBgnTuagV7HT2w

cailmdaley and others added 5 commits July 12, 2026 23:34
Attach the local (Jacobian) WCS at each object position to the PSF/star
stamps and pass use_sky_coords=True to galsim FindAdaptiveMom, so PSFEx and
MCCD HSM ellipticities and sizes come out already in world coordinates. This
makes the downstream WCS-Jacobian shape rotation (convert_psf_pix2world)
redundant.

Threaded at every in-repo HSM PSF/star call site:
- psfex_interp: multi-epoch (per-CCD WCS log), single-exposure validation and
  classic (WCS from the galcat FITS_LDAC image header via _galcat_wcs)
- mccd_interpolation_script: multi-epoch (per-CCD WCS log)
- shapepipe_auxiliary_mccd: interpolation (WCS from galcat LDAC header)

local_wcs_list() builds the per-object local WCS: use_sky_coords evaluates the
stamp WCS at true_center, so it must be the full WCS linearised at the object's
image position, not at the stamp centre.

SIGMA_*_HSM (hence HSM_T_* = 2 sigma^2) is now in sky units (arcsec) rather
than pixels; this is the intended, physically correct size.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01Gao7Pi5QDBgnTuagV7HT2w
HSM ellipticities and sizes are measured directly in world coordinates
upstream, so the WCS-Jacobian shape rotation (transform_shape,
getDecomposition) is redundant. Remove it from both the PSFEx and MCCD
branches; the script keeps its remaining job of collating positions (X/Y/
RA/DEC), assigning the MCCD focal-plane CCD id, and merging per-exposure
catalogues into per-patch star catalogues. Drops the galsim dependency.

Caveat documented in the module docstring: the MCCD PSF_MOM_LIST/STAR_MOM_LIST
columns come from the external mccd fit-validation code, still measured in the
pixel frame; they are now passed through unrotated. Update the CANFAR doc to
describe the step as collation, not pix->world shape conversion.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01Gao7Pi5QDBgnTuagV7HT2w
Equivalence test for use_sky_coords=True: draw an elliptical Gaussian on a
stamp with a nontrivial local WCS (rotation + shear + 0.187 arcsec/pix scale,
positive determinant), and assert that measuring in the pixel frame then
applying the removed convert_psf_pix2world.transform_shape rotation (copied in
verbatim as reference) agrees with measuring directly in sky coordinates, to
numerical precision, for e1/e2 and sigma. A second test pins the deliberate
unit change: the sky-frame sigma is the pixel sigma times the WCS scale.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01Gao7Pi5QDBgnTuagV7HT2w
…ranch

The MCCD PSF_MOM_LIST/STAR_MOM_LIST columns come from the external mccd
fit-validation code (mccd.auxiliary_fun.mccd_validation), which still
measures HSM adaptive moments in the pixel frame. Stripping the WCS-Jacobian
rotation there (71d502f) would have handed MergeStarCatMCCD pixel-frame
ellipticities. Restore transform_shape + the per-CCD getDecomposition rotation
for the MCCD MOM_LIST branch only (and the galsim import it needs); the in-repo
PSFEx and MCCD-interpolation paths stay unrotated, as they now measure in world
coordinates upstream. Update the module docstring caveat: this branch keeps the
rotation until mccd itself adopts use_sky_coords.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01Gao7Pi5QDBgnTuagV7HT2w
local_wcs_list wraps galsim.AstropyWCS, which follows the 1-indexed FITS pixel
convention, but the multi-epoch paths fed it positions from all_world2pix(...,
0) (origin 0) — evaluating the local Jacobian one pixel off. Compute a separate
origin=1 position array (gal_pos_wcs) for local_wcs_list; the origin-0 gal_pos
stays for the PSF-model interpolation (interpsfex / interp_MCCD), which is its
only other use.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01Gao7Pi5QDBgnTuagV7HT2w
@martinkilbinger

Copy link
Copy Markdown
Contributor

We should revisit convert_psf_pix2world.py. In v2.0 we will no longer define patches. Renaming or re-scoping this script should reflect its main purpose if it is still needed: to merge the star catalogues. If pixel to WCS is needed for MCCD, this should just be implemented in the MCCD library, we can check in with @tobias-liaudat .

PSF/star ellipticities and sizes are now measured directly in sky coordinates
during PSF interpolation (galsim `FindAdaptiveMom(use_sky_coords=True)`), so
`convert_psf_pix2world.py` only collates and passes the shapes through.

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Update to reflect removal of the patch concept in v2.0

cailmdaley added a commit that referenced this pull request Jul 16, 2026
Martin (PR #812 review): the patch concept is removed in v2.0. For
-V v2.0 the handler reads summary/missing_job_32_all.txt and
exp_numbers.txt relative to the run root instead of iterating P1..P9;
-p restores an explicit P1..P<n> layout when needed.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01MqH9wjbJicVmaQdhZ4kBYr
@cailmdaley cailmdaley mentioned this pull request Jul 16, 2026
12 tasks
…add patch-less v2.0

Re-scope the script to its actual remaining purpose: collating per-exposure PSF
validation catalogues into star catalogues. The WCS-Jacobian shape rotation is
gone (HSM moments are measured in sky coordinates upstream), so "convert" no
longer describes it. The external-mccd PSF_MOM_LIST/STAR_MOM_LIST rotation
branch stays until mccd adopts use_sky_coords.

Add a -V/--version option mirroring get_ccds_with_psf. v1.x is unchanged
(loop over -P patches, <base>/P<patch>/output, validation_psf_conv-<patch>-<idx>).
v2.0 is patch-less: single <base>/output root, patch token dropped from the
output filename (validation_psf_conv-<idx>). The downstream combine_runs.bash
glob (validation_psf_conv-*) matches both.

Path/filename logic extracted into pure collate_paths/output_filename helpers
with unit tests; docs collate section updated with both invocations. The script
command name follows the file (Dockerfile symlinks by basename).

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01WqNFhq9uGpDpprYeuGeb8k
@cailmdaley
cailmdaley force-pushed the feat/hsm-sky-coords branch from 79fafe0 to baea03a Compare July 20, 2026 19:51
@cailmdaley

Copy link
Copy Markdown
Contributor Author

Thanks Martin — both points addressed.

  • convert_psf_pix2world.py is renamed to collate_star_cat.py, re-scoped to its remaining purpose: collating the per-exposure PSF validation catalogues into star catalogues. No shape rotation happens there any more, with one exception: the external-mccd PSF_MOM_LIST/STAR_MOM_LIST columns are still measured in the pixel frame upstream, so that one branch keeps the rotation until mccd adopts use_sky_coords — agreed that belongs in the MCCD library; we can raise it with Tobias separately.
  • The script now takes -V following the same convention as get_ccds_with_psf in Coverage nexp masks #797, defaulting to the patch-less v2.0 (single run root); v1.x runs pass -V v1.6 -P … and behave as before. The collate section of pipeline_canfar.md shows both.

— Fable on behalf of Cail

cailmdaley added a commit that referenced this pull request Jul 20, 2026
Align with the collate_star_cat convention (PR #812): v2.0 (patch-less) is
now the default catalogue version. Older, patched layouts are selected
explicitly with -V. Help text and the pipeline_canfar.md invocations updated;
the v2.0 call drops the flag, v1.x shows -V v1.6. Tests set the version
explicitly, so none pinned the old default.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01WqNFhq9uGpDpprYeuGeb8k
@aguinot

aguinot commented Jul 21, 2026

Copy link
Copy Markdown
Contributor

Hello,
While you're looking into this, it might be good to add higher order moments. It has been very useful in other surveys.

@cailmdaley

Copy link
Copy Markdown
Contributor Author

Good point — agreed, and it turns out this was attempted before: #698 (fourth-order moments, PSFHOME-style) was closed unmerged and is now stale. It edits exactly the two functions this PR rewrites. We've reopened #697/#698 to track it and will implement as a follow-up on top of this PR: the spin-2 fourth-moment combinations measured in the sky frame (which removes the pixel→sky rotation the old code needed), plus galsim's spin-0 moments_rho4, with tests for both.

— Fable on behalf of Cail

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

HSM: measure PSF/star shapes directly in sky coords (use_sky_coords), drop the WCS re-transform

3 participants