feat: WFC3 (UVIS + IR) reduction support (phase 2)

autoreduce/instruments/__init__.py

@@ -7,3 +7,5 @@
 
 from .adapter import InstrumentAdapter, get, register, registered_keys
 from .acs_wfc import ACS_WFC
+from .wfc3_uvis import WFC3_UVIS
+from .wfc3_ir import WFC3_IR

autoreduce/instruments/adapter.py

@@ -20,7 +20,10 @@ class InstrumentAdapter:
     crds_reference_subpath: str  # where CRDS syncs this instrument's files
     supports_cte_correction: bool
     default_drizzle_kwargs: Dict[str, object]
-    saturation_dn: float  # conservative full-well / saturation level
+    saturation_dn: float  # conservative full-well / saturation level, electrons
+    # The adapter's recommendation for TargetSpec.final_scale (which remains
+    # the user-facing dial); documents sensible sampling for this detector.
+    recommended_final_scale: float = 0.05
 
     def scale_ratio(self, final_scale: float) -> float:
         """s = output scale / native scale, as used by the Casertano factor."""

autoreduce/instruments/wfc3_ir.py

@@ -0,0 +1,34 @@
+"""
+WFC3/IR — adapter #3 (roadmap phase 2), the genuinely different path.
+
+No CTE correction exists for the IR channel, so the calibrated product is
+``_flt`` (already in e-/s from up-the-ramp fitting, which also rejects most
+cosmic rays per read — AstroDrizzle's driz_cr then handles the residue when
+multiple exposures exist). Native scale ~0.128″/pix under-samples the PSF, so
+dithered programs conventionally drizzle to a finer grid; 0.065″/pix is the
+adapter's recommendation (half-native, within the 0.06–0.08 range common in
+deep-field practice) — `TargetSpec.final_scale` remains the user dial and
+star-poor or poorly-dithered data may prefer coarser values.
+"""
+
+from .adapter import InstrumentAdapter, register
+
+WFC3_IR = register(
+    InstrumentAdapter(
+        key="wfc3_ir",
+        mast_instrument_name="WFC3/IR",
+        native_scale=0.128,
+        calibrated_suffix="FLT",
+        reference_env_key="iref",
+        crds_reference_subpath="references/hst/wfc3",
+        supports_cte_correction=False,
+        default_drizzle_kwargs={
+            "skymethod": "globalmin+match",
+            "final_wht_type": "IVM",
+            "final_units": "cps",
+            "final_rot": 0.0,
+        },
+        saturation_dn=78_000.0,
+        recommended_final_scale=0.065,
+    )
+)

autoreduce/instruments/wfc3_uvis.py

@@ -0,0 +1,30 @@
+"""
+WFC3/UVIS — adapter #2 (roadmap phase 2).
+
+The ACS-like path: CTE-corrected ``_flc`` exposures, ``iref`` references,
+cps/IVM/north-up outputs. Native plate scale 0.0396″/pix — the published
+lensing anchor is the Bayer et al. (arXiv:1803.05952) F390W reduction of
+SDSS J0252+0039 at exactly this output scale with pixfrac 1.0.
+"""
+
+from .adapter import InstrumentAdapter, register
+
+WFC3_UVIS = register(
+    InstrumentAdapter(
+        key="wfc3_uvis",
+        mast_instrument_name="WFC3/UVIS",
+        native_scale=0.0396,
+        calibrated_suffix="FLC",
+        reference_env_key="iref",
+        crds_reference_subpath="references/hst/wfc3",
+        supports_cte_correction=True,
+        default_drizzle_kwargs={
+            "skymethod": "globalmin+match",
+            "final_wht_type": "IVM",
+            "final_units": "cps",
+            "final_rot": 0.0,
+        },
+        saturation_dn=63_000.0,
+        recommended_final_scale=0.0396,
+    )
+)

autoreduce/pipeline.py

@@ -118,11 +118,21 @@ def reduce_target(
 
     target_xy = WCS(header).world_to_pixel_values(spec.ra, spec.dec)
     selection = stars_mod.StarSelection()
+    # Saturation is per exposure, not per stack: a star saturates when its
+    # rate fills the well within one exposure, so the cps cap divides the
+    # full well by the longest single-exposure time — never the mosaic total.
+    max_single_exptime = max(
+        float(fits.getheader(p).get("EXPTIME", 0.0)) for p in exposures
+    )
+    if max_single_exptime <= 0.0:
+        raise ValueError("no exposure carries a positive EXPTIME header")
     stars = stars_mod.find_stars(
         sci,
         selection,
         target_xy=(float(target_xy[0]), float(target_xy[1])),
-        peak_max=selection.saturation_fraction * adapter.saturation_dn / float(exptime),
+        peak_max=selection.saturation_fraction
+        * adapter.saturation_dn
+        / max_single_exptime,
     )
     psf, psf_full, psf_diag = epsf_mod.build_epsf(
         sci, stars, spec.psf_shape, spec.psf_full_shape

docs/design/roadmap.md

@@ -13,17 +13,14 @@ geometry, calibrated-product naming (`_flc` vs `_flt` vs `_cal`), units and
 gain, recommended combine parameters, and the PSF-model source. ACS/WFC is
 adapter #1; nothing outside `instruments/` may mention a detector by name.
 
-## HST/WFC3 (IR + UVIS)
-
-- Same stage graph; new adapter. UVIS is ACS-like (CTE-corrected `_flc`,
-  drizzle); IR differs (no CTE correction, up-the-ramp sampling, different
-  DQ semantics, 0.13″ native pixels → different `final_scale` choice).
-- The [ajshajib/hst-lens](https://github.com/ajshajib/hst-lens) notebooks
-  (download / IR / UVIS) serve as a step checklist to audit against — not as
-  architecture. Their gaps (no ACS, notebook-form, unclear noise treatment)
-  are what this project exists to fix.
-- Other ACS filters (F435W, F606W…) are config, not code: the adapter already
-  parameterizes filter-dependent zero-points and PSF models.
+## HST/WFC3 (IR + UVIS) — **in progress (phase 2, PyAutoReduce#4)**
+
+- Design deltas live in [`wfc3.md`](wfc3.md); adapters `wfc3_uvis` /
+  `wfc3_ir` implemented. UVIS is ACS-like (CTE-corrected `_flc`); IR differs
+  (`_flt`, no CTE correction, up-the-ramp CR rejection, 0.128″ native →
+  recommended 0.065″ output).
+- Other ACS/WFC3 filters (F435W, F606W…) are config, not code: the adapter
+  already parameterizes the filter-dependent pieces.
 
 ## JWST (NIRCam first)
 

docs/design/wfc3.md

@@ -0,0 +1,65 @@
+# WFC3 (UVIS + IR) — per-stage deltas vs the ACS/WFC design
+
+Phase 2. `hst_acs_pipeline.md` remains the stage-by-stage spec; this page
+records only what WFC3 changes. Everything below lives in the two adapters
+(`instruments/wfc3_uvis.py`, `instruments/wfc3_ir.py`) — no stage module
+mentions a detector.
+
+## WFC3/UVIS — the ACS-like path
+
+| Stage | Delta vs ACS |
+|-------|--------------|
+| acquire | `_flc` (CTE-corrected) as ACS; references key **`iref`**, CRDS subpath `references/hst/wfc3` |
+| align / drizzle / noise / psf / package | unchanged — same recipes, same dials |
+| scale | native **0.0396″/pix**; adapter recommends output at native scale |
+| saturation | ~63 ke- full well (star selection peak cut) |
+
+**Validation anchor:** the published [Bayer et al.
+(arXiv:1803.05952)](https://arxiv.org/abs/1803.05952) F390W reduction of
+SDSS J0252+0039 — output 0.0396″/pix, **pixfrac 1.0**, noise recipe
+σ = √(N/W + σ²_sky) with σ_sky ≈ 0.002 e-/s. Our integration script
+(`scripts/reduce_j0252_wfc3.py --channel uvis`) reduces the same data with those dials and
+checks: output units e-/s, empirical σ_sky in that regime, noise-map
+consistency, WHT uniformity. (Note their correlated-noise treatment is
+blank-sky *realizations*, not the scalar R — with pixfrac 1.0 at native
+scale, R = 1.5; comparisons account for whether R is applied.)
+
+## WFC3/IR — the genuinely different path
+
+| Stage | Delta vs ACS |
+|-------|--------------|
+| acquire | **`_flt`** — no CTE correction exists for the IR channel; `iref` references |
+| drizzle | up-the-ramp fitting in `calwf3` already rejects most CRs per read; `driz_cr` still runs on multi-exposure stacks for the residue (defaults-first) — documented, revisit if IR integrations show over-flagging |
+| scale | native **0.128″/pix** under-samples the PSF; adapter recommends **0.065″/pix** for dithered data (half-native, in the 0.06–0.08 deep-field range). The dial stays user-facing; the fine-grid Casertano branch (s < p) then applies, so R is materially larger — reported per run as always |
+| units | `_flt` IR data are already e-/s (count rates); `final_units='cps'` unchanged |
+| saturation | ~78 ke- effective full well |
+| psf | same tiers; STScI focus-diverse ePSF grids exist for IR when tier 2 lands |
+
+## Coverage audit vs `ajshajib/hst-lens` (the checklist, not the architecture)
+
+Their three notebooks (Download / IR / UVIS) cover: archive download,
+per-channel calibrated products, AstroDrizzle combination, and cutouts. Ours
+adds what they lack: instrument adapters (theirs is notebook-per-channel),
+provenance (`reduction.json`), an explicit noise-map recipe with correlated-
+noise handling, tiered PSF construction with diagnostics, cache/eviction, and
+loud-failure contracts. Nothing in their steps is absent from our stage graph.
+
+## Integration finding — IR pixfrac/coverage (2026-07-08)
+
+Reducing the J0252+0039 F160W snapshot (program 11202) at the recommended
+0.065″/pix with the phase-1 default pixfrac 0.8 left **230 zero-weight
+speckle pixels inside the cutout** — the few-dither + fine-grid (s ≈ 0.51 <
+p) regime — and the finite-noise packaging guard refused to ship the dataset,
+exactly as designed. pixfrac 1.0 closes coverage at the cost of a larger
+correlated-noise factor (reported per run, as always). Rule of thumb recorded
+here: **on the IR channel, few-dither data at sub-native output scales needs
+pixfrac → 1.0 (or a coarser scale dial)**; the WHT-uniformity diagnostic and
+the finite-noise guard enforce the trade-off loudly rather than letting a
+holey dataset through.
+
+## Open items
+
+- IR integration target: discovered via MAST at run time; if no suitable
+  IR lens dataset is reachable, the leg parks as a batched question.
+- Tier-2 PSF for IR (focus-diverse ePSF grids) — with the roadmap's tier-2
+  work, not phase 2.

scripts/reduce_j0252_wfc3.py

@@ -0,0 +1,147 @@
+"""
+WFC3 integration + acceptance (issue #4): SDSS J0252+0039.
+
+--channel uvis : F390W through wfc3_uvis at the published Bayer et al.
+                 (arXiv:1803.05952) dials — 0.0396"/pix, pixfrac 1.0 — then
+                 check our products against their published numbers (units
+                 e-/s; sigma_sky ~ 0.002 e-/s; noise closure; R accounting).
+--channel ir   : discover WFC3/IR imaging at the same coordinates via MAST
+                 and reduce it at the adapter-recommended 0.065"/pix;
+                 internal validation only (no published anchor).
+
+Run:  ~/venv/PyAuto/bin/python scripts/reduce_j0252_wfc3.py --channel uvis
+Network + drizzlepac required; unit tests never import this.
+"""
+
+import argparse
+import json
+import sys
+from pathlib import Path
+
+import numpy as np
+
+REPO = Path(__file__).resolve().parent.parent
+sys.path.insert(0, str(REPO))
+
+from autoreduce import TargetSpec, reduce_target  # noqa: E402
+from autoreduce import instruments  # noqa: E402
+
+# SDSS J0252+0039: 02h52m45.21s +00d39m58.4s (legacy dataset info.json)
+RA, DEC = 43.188375, 0.666222
+CACHE_ROOT = REPO / "scripts" / "cache"
+OUTPUT_ROOT = REPO / "scripts" / "output"
+
+BAYER_SIGMA_SKY = 0.002  # e-/s, published for the F390W reduction
+
+
+def spec_for(channel: str) -> TargetSpec:
+    if channel == "uvis":
+        return TargetSpec(
+            name="j0252+0039_f390w",
+            ra=RA,
+            dec=DEC,
+            instrument="wfc3_uvis",
+            filter_name="F390W",
+            final_scale=0.0396,   # Bayer dial
+            final_pixfrac=1.0,    # Bayer dial
+        )
+    if channel == "ir":
+        filter_name = discover_ir_filter()
+        return TargetSpec(
+            name=f"j0252+0039_{filter_name.lower()}",
+            ra=RA,
+            dec=DEC,
+            instrument="wfc3_ir",
+            filter_name=filter_name,
+            final_scale=instruments.get("wfc3_ir").recommended_final_scale,
+            # Few-dither snapshot on a half-native grid: pixfrac 0.8 leaves
+            # zero-weight speckle (the finite-noise guard caught it); the
+            # full drop closes coverage at the cost of a larger R — the
+            # dial trade-off working as documented.
+            final_pixfrac=1.0,
+            # 281 px at 0.065" spans 18.3" — far more sky than the ACS-era
+            # 14" footprint, and it clips a zero-coverage detector-defect
+            # blob 8.5" from the lens (the guard refused it). Match the ACS
+            # sky footprint instead: 14" / 0.065 -> 215 px.
+            cutout_shape=(215, 215),
+        )
+    raise ValueError(channel)
+
+
+def discover_ir_filter() -> str:
+    """Find which WFC3/IR filter (if any) covers the target."""
+    from astropy.coordinates import SkyCoord
+    from astroquery.mast import Observations
+
+    from autoreduce.acquire.mast import select_observations
+
+    obs = Observations.query_criteria(
+        coordinates=SkyCoord(RA, DEC, unit="deg"),
+        radius="0.5 arcmin",
+        obs_collection="HST",
+        instrument_name="WFC3/IR",
+        dataproduct_type="image",
+    )
+    direct = select_observations(obs)
+    if not direct:
+        sys.exit(
+            "[ir] no direct WFC3/IR observations at J0252+0039 — the IR leg "
+            "needs a different target (parked question on issue #4)"
+        )
+    # HAP composite rows carry the pseudo-filter 'detection'; never reduce it.
+    filters = sorted(
+        {str(row["filters"]) for row in direct} - {"detection"}
+    )
+    if not filters:
+        sys.exit("[ir] only HAP composite products found — no real filter")
+    print(f"[ir] direct WFC3/IR observations found; filters: {filters}")
+    preferred = [f for f in filters if f == "F160W"] or filters
+    return preferred[0]
+
+
+def validate(channel: str, record: dict, out_dir: Path):
+    from astropy.io import fits
+
+    noise = fits.getdata(out_dir / "noise_map.fits").astype(float)
+    data = fits.getdata(out_dir / "data.fits").astype(float)
+    with fits.open(out_dir / "data.fits") as hdul:
+        bunit = hdul[0].header.get("BUNIT", "unknown")
+
+    r_factor = record["noise"]["correlated_noise_factor"]
+    sky_rms = record["noise"]["empirical_background_rms"]
+    summary = {
+        "channel": channel,
+        "n_exposures": record["acquire"]["n_exposures"],
+        "bunit": bunit,
+        "weight_uniformity": record["drizzle"]["weight_uniformity"],
+        "correlated_noise_factor": r_factor,
+        "empirical_sky_rms_cps": sky_rms,
+        "noise_map_min_median": [float(np.nanmin(noise)), float(np.nanmedian(noise))],
+        "psf": record["psf"],
+    }
+    if channel == "uvis":
+        # Published anchor: Bayer et al. sigma_sky ~ 0.002 e-/s at these dials.
+        # Our empirical sky RMS is pre-R; theirs enters sigma_n pre-realization.
+        summary["bayer_sigma_sky_cps"] = BAYER_SIGMA_SKY
+        summary["sky_rms_over_bayer"] = sky_rms / BAYER_SIGMA_SKY
+        # Noise-floor closure: in blank sky the map should approach R * sky_rms.
+        summary["noise_floor_over_R_times_sky"] = float(
+            np.nanpercentile(noise, 5) / (r_factor * sky_rms)
+        )
+    print(f"[{channel}] ---- validation ----")
+    print(json.dumps(summary, indent=2))
+    (out_dir / "validation_summary.json").write_text(json.dumps(summary, indent=2))
+
+
+def main():
+    parser = argparse.ArgumentParser(description=__doc__)
+    parser.add_argument("--channel", required=True, choices=["uvis", "ir"])
+    args = parser.parse_args()
+
+    spec = spec_for(args.channel)
+    record = reduce_target(spec, cache_root=CACHE_ROOT, output_root=OUTPUT_ROOT)
+    validate(args.channel, record, OUTPUT_ROOT / spec.name)
+
+
+if __name__ == "__main__":
+    main()