extract_spectrum_optimal¶

Optimal (Horne 1986) aperture extraction of ctx.image.

Faithful Horne 1986 recipe (Table I) : sky from ctx.extras / the flanking rows / none ; spatial profile from low-order polynomial fits along the dispersion of the normalised data (D − S)/f per detector row (eq. 14, weights f²/V, σ-clipped, P ≥ 0, unit sum per column) ; model variance V = (RON/gain)² + |f·P + S|/gain (eq. 13) ; optimal sum f = Σ M P (D − S)/V ÷ Σ M P²/V with variance 1/Σ M P²/V (eqs. 8–9) ; cosmic rays rejected one pixel per column per iteration (§II.E) ; the profile → extraction sequence runs twice, the second time normalised by the optimal spectrum and without the rejected pixels (§II.F). v2.0.0 (audit A4/A9) : v1 built the variance from the data (|D|/gain + RON²) and the profile from a 21-column running median of the raw cutouts, whose noise is correlated with the pixel noise — measured bias on a σ = 1.5 px Gaussian trace at gain 1, RON 3 e⁻, 65 536 columns : −23.1 % at 5 e⁻/column, −13.7 % at 20, −4.8 % at 100, −0.8 % at 1000 (v2 : −0.3 % ± 0.6 %, +0.6 %, +0.1 %, 0.0 %, the last three identical to the boxcar control on the same frames) ; v1 added RON² in e⁻² to a variance in ADU², so the reported σ was ×2.05 the empirical scatter at gain 2.5 / RON 6 e⁻, 100 e⁻/column (v2 : ×1.00, 0.995–1.03 over 20–1000 e⁻) ; v1 rejected every pixel above the threshold at once. The aperture is now the symmetric window round(y) ± trace_half_width (2·hw + 1 rows) instead of [int(y) − hw, int(y) + hw) — the old window sat 0.5–1.5 px below the trace, so the captured fraction rippled with the fractional trace position (0.947 at row 20.0, 0.860 at 20.9 for σ = 1.5 px, hw = 3 ; now 0.983, 0.976 at the worst half-pixel offset, 0.982 at 20.9). Profile normalisation : a column is normalised by its own flux when its S/N ≥ 10, else by the running median of the spectrum over profile_smooth_columns — Horne's per-column normaliser attenuates the profile at low S/N (errors-in-variables ; +17 % at 5 e⁻/column, profile peak 0.195 instead of 0.266), which biases the sum ; the running median keeps the estimator unbiased down to a few e⁻ per column. Single-pixel spikes along the dispersion are excluded from the profile fit only. The flanking sky fit is unweighted : Horne's 1/V weights (eq. 12) bias a Poisson mean by ≈ −1 count per pixel, i.e. +(2·hw + 1) ADU per column in the flux. Uncertainty : √(1/Σ M P²/V), sky taken as noise-free (Horne §II.A). Edge columns keep the aperture rows that fit the detector. Strongly tilted or curved traces (drift of several rows) need a higher profile_poly_order or a prior geometric rectification (correct_slant_affine / correct_smile_polynomial) — the profile of a detector row is fitted along the dispersion in fixed rows, as in the paper, not in trace-relative coordinates (Marsh 1989). When to choose this vs the lighter alternatives : extract_spectrum_sum is the simplest pure sum (no trace fit) ; extract_spectrum_boxcar is the edge-safe boxcar with optional Gaussian profile weighting but no variance model, no uncertainty propagation and no cosmic rejection — pick it for high-SNR sources or extended emission. Pick extract_spectrum_optimal when SNR matters most (faint, read-noise-limited targets), when residual cosmics may remain after clip_cosmic_rays, or when the downstream pipeline needs per-pixel Spectrum1D.uncertainty (e.g. radial-velocity error propagation). The pixel axis is 0-based (numpy.arange(npix)).

Parameters¶

Use it¶

from spectro_kernel import run_algorithm

output = run_algorithm("extract_spectrum_optimal", ctx, {
    "trace_method": "argmax",
    "trace_poly_order": 2,
    "trace_y_pixel_range": 15,
    "trace_peak_height": 100.0,
    "trace_peak_distance": 50,
    "trace_half_width": 12,
    "gain": 1.0,
    "readnoise": 6.0,
    "sky_source": "auto",
    "sky_extras_key": "sky_spectrum",
    "sky_offset": 4,
    "sky_half_width": 10,
    "sky_poly_order": 1,
    "profile_poly_order": 2,
    "profile_smooth_columns": 401,
    "reject_sigma": 5.0,
    "iterations": 3
})

{
  "tool": "extract_spectrum_optimal",
  "arguments": {
    "session_id": "<session_id>",
    "params": {
      "trace_method": "argmax",
      "trace_poly_order": 2,
      "trace_y_pixel_range": 15,
      "trace_peak_height": 100.0,
      "trace_peak_distance": 50,
      "trace_half_width": 12,
      "gain": 1.0,
      "readnoise": 6.0,
      "sky_source": "auto",
      "sky_extras_key": "sky_spectrum",
      "sky_offset": 4,
      "sky_half_width": 10,
      "sky_poly_order": 1,
      "profile_poly_order": 2,
      "profile_smooth_columns": 401,
      "reject_sigma": 5.0,
      "iterations": 3
    }
  }
}

References¶

• Horne 1986, PASP 98, 609 — optimal extraction algorithm for CCD spectroscopy.
• Marsh 1989, PASP 101, 1032 — empirical profile for tilted traces.
• Tody 1986, SPIE 627, 733 — IRAF apall heritage.

Related algorithms¶

• detect_trace - Detect the spectral trace on ctx.image and expose it on ctx.extras["trace"].
• extract_sky_lateral_bands - Median-combine two off-trace sky bands into a 1D pixel-axis reference.
• extract_spectrum_boxcar - Trace (easyspec) + pure-numpy aperture extraction on ctx.image.
• extract_spectrum_easyspec - Trace and extract a 1D spectrum from a 2D frame via easyspec.
• extract_spectrum_sum - Extract a 1D spectrum from ctx.image by summing across the trace.