transform

Documentation

Overview

Package transform computes forward and inverse FITS World Coordinate System mappings between image pixel coordinates and celestial spherical coordinates.

The coordinate pipeline follows Greisen & Calabretta 2002 Paper I §2.1:

1. Pixel offset from reference point: (p_i - CRPIX_i)

2. Linear transform to intermediate world coordinates (applied as either a CD matrix or a PC matrix with CDELT scaling): q_i = sum_j CD_{ij} * (p_j - CRPIX_j) or q_i = CDELT_i * sum_j PC_{ij} * (p_j - CRPIX_j)

3. Projection from intermediate world coordinates to native spherical coordinates (phi, theta) — this is the projection-specific step.

4. Spherical rotation from native to celestial coordinates (alpha, delta) using LONPOLE / LATPOLE and the reference (CRVAL_lon, CRVAL_lat).

Supported projections

Every projection from Paper II is either fully implemented or recognized and returning ErrUnsupportedProjection — never a silent wrong answer.

Zenithal family (§5.1): AZP SZP TAN STG SIN ARC ZPN ZEA Cylindrical (§5.2): CYP CEA CAR MER Pseudo-cylindrical (§5.3): SFL PAR MOL AIT Conic (§5.4): COP COE COD COO Polyconic (§5.5): BON PCO Quad-cube (§5.6): TSC CSC QSC HEALPix (CR 2007): HPX XPH

Index

• Variables
• func EclipticToEquatorial(lambda, beta, jd float64) (alpha, delta float64)
• func EquatorialToEcliptic(alpha, delta, jd float64) (lambda, beta float64)
• func EquatorialToGalactic(alpha, delta float64) (l, b float64)
• func EquatorialToSupergalactic(alpha, delta float64) (sgl, sgb float64)
• func GalacticToEquatorial(l, b float64) (alpha, delta float64)
• func GalacticToSupergalactic(l, b float64) (sgl, sgb float64)
• func PrecessFK5(alpha, delta, equinox1, equinox2 float64) (a2, d2 float64)
• func SupergalacticToEquatorial(sgl, sgb float64) (alpha, delta float64)
• func SupergalacticToGalactic(sgl, sgb float64) (l, b float64)
• type Projection
  • func Select(code string, pv map[wcs.PVKey]float64, latAxis int) (Projection, error)
• type SIP
  • func ParseSIP(h *header.Header) (*SIP, error)
  • func (s *SIP) Forward(u, v float64) (uOut, vOut float64)
  • func (s *SIP) Inverse(uPrime, vPrime float64) (u, v float64, ok bool)
• type SIPPoly
• type TNX
  • func ParseTNX(h *header.Header) (*TNX, error)
  • func (t *TNX) Forward(xi, eta float64) (xiOut, etaOut float64)
  • func (t *TNX) Inverse(xiPrime, etaPrime float64) (xi, eta float64, ok bool)
• type TPV
  • func ParseTPV(w *wcs.Header) *TPV
  • func (t *TPV) Forward(xi, eta float64) (xiOut, etaOut float64)
  • func (t *TPV) Inverse(xiPrime, etaPrime float64) (xi, eta float64, ok bool)
• type Transform
  • func New(h *wcs.Header) (*Transform, error)
  • func (t *Transform) Code() string
  • func (t *Transform) PixelToSky(p1, p2 float64) (alpha, delta float64, err error)
  • func (t *Transform) SkyToPixel(alpha, delta float64) (p1, p2 float64, err error)

Variables

var ErrUnknownProjection = errors.New("wcs/transform: unknown projection code")

ErrUnknownProjection is returned by Select when the 3-letter code is not recognized at all (i.e. it is not in the Paper II projection list).

var ErrUnsupportedProjection = errors.New("wcs/transform: unsupported projection")

ErrUnsupportedProjection is a forward-compatibility sentinel returned by Select for FITS projection codes that a future standard revision might define but this package does not support. Currently unreachable — all 27 standard projections plus HPX/XPH are implemented.

Functions

func EclipticToEquatorial

func EclipticToEquatorial(lambda, beta, jd float64) (alpha, delta float64)

EclipticToEquatorial is the inverse of EquatorialToEcliptic.

func EquatorialToEcliptic

func EquatorialToEcliptic(alpha, delta, jd float64) (lambda, beta float64)

EquatorialToEcliptic converts (alpha, delta) in ICRS degrees to ecliptic (lambda, beta) at the given Julian Date (typically derived from DATE-OBS or MJD-OBS).

func EquatorialToGalactic

func EquatorialToGalactic(alpha, delta float64) (l, b float64)

EquatorialToGalactic converts (alpha, delta) in ICRS degrees to galactic (l, b) in degrees.

func EquatorialToSupergalactic

func EquatorialToSupergalactic(alpha, delta float64) (sgl, sgb float64)

EquatorialToSupergalactic chains ICRS → galactic → supergalactic.

func GalacticToEquatorial

func GalacticToEquatorial(l, b float64) (alpha, delta float64)

GalacticToEquatorial is the inverse.

func GalacticToSupergalactic

func GalacticToSupergalactic(l, b float64) (sgl, sgb float64)

GalacticToSupergalactic converts galactic (l, b) to supergalactic (SGL, SGB).

func PrecessFK5

func PrecessFK5(alpha, delta, equinox1, equinox2 float64) (a2, d2 float64)

PrecessFK5 precesses equatorial coordinates from equinox equinox1 to equinox2, both given as Besselian/Julian years (e.g. 2000.0 for J2000). Uses the IAU 1976 precession model via the standard p-matrix formulas.

This is an approximation adequate for arcsecond-level precision over ~1 century; higher-precision work should use the IAU 2006 precession or delegate to an ephemeris library.

func SupergalacticToEquatorial

func SupergalacticToEquatorial(sgl, sgb float64) (alpha, delta float64)

SupergalacticToEquatorial is the inverse.

func SupergalacticToGalactic

func SupergalacticToGalactic(sgl, sgb float64) (l, b float64)

SupergalacticToGalactic is the inverse.

Types

type Projection

type Projection interface {
	Code() string
	Theta0() float64
	Forward(phi, theta float64) (x, y float64, ok bool)
	Inverse(x, y float64) (phi, theta float64, ok bool)
}

Projection is the interface implemented by every spherical projection.

All angles are in RADIANS. Implementations work with FITS "native" longitude phi and "native" latitude theta on the unit sphere.

Forward: given (phi, theta), return intermediate world coordinates (x, y) — also in radians. ok=false indicates the input falls outside the projection's valid domain (e.g. the far hemisphere for SIN).

Inverse: given intermediate (x, y), return the native (phi, theta). ok=false indicates the point is not on the projection's image.

Theta0 returns the native latitude of the projection's reference (fiducial) point, in radians. Zenithal projections have theta_0 = pi/2; most others have theta_0 = 0. Used by the transform constructor to apply the correct LONPOLE / LATPOLE defaults per Paper II §2.4.

func Select

func Select(code string, pv map[wcs.PVKey]float64, latAxis int) (Projection, error)

Select returns the Projection for a 3-letter FITS projection code, configured from the projection-parameter map pv (PVi_m values keyed by (axis, index) — conventionally axis 2 for the latitude axis). The latAxis argument is the 1-based latitude axis index; most callers pass 2 (the common convention where CRVAL2 is declination).

Returns ErrUnknownProjection for codes not in the FITS standard list. All 27 standard projections plus HPX/XPH are implemented, so ErrUnsupportedProjection is currently unreachable in normal use and exists only as a forward-compatibility sentinel.

type SIP

type SIP struct {
	A  *SIPPoly
	B  *SIPPoly
	AP *SIPPoly // may be nil (solved iteratively if absent)
	BP *SIPPoly // may be nil
}

SIP holds a parsed SIP distortion: forward (A, B) and optional inverse (AP, BP) polynomial pairs.

func ParseSIP

func ParseSIP(h *header.Header) (*SIP, error)

ParseSIP extracts SIP coefficients from the given image header. Returns nil, nil if the header does not declare SIP (i.e. CTYPE suffix is not -SIP, or A_ORDER / B_ORDER is absent). An error indicates a malformed SIP section.

func (*SIP) Forward

func (s *SIP) Forward(u, v float64) (uOut, vOut float64)

Forward applies the SIP forward correction to a pixel offset (u, v) from CRPIX. Returns the corrected (u', v') that should be fed into the linear CD step.

func (*SIP) Inverse

func (s *SIP) Inverse(uPrime, vPrime float64) (u, v float64, ok bool)

Inverse applies the SIP inverse correction. If AP/BP are present, it evaluates them directly. Otherwise it iterates Newton's method against the forward polynomial using the analytic Jacobian.

The forward map is F(u, v) = (u + A(u, v), v + B(u, v)), so the Jacobian is

J = I + [∂A/∂u  ∂A/∂v]
        [∂B/∂u  ∂B/∂v]

Both partial derivatives are computed in a single walk over the A/B coefficient tables, making each iteration cost the same as a forward evaluation (vs. 3× cost for the old numerical version).

type SIPPoly

type SIPPoly struct {
	Order  int
	Coeffs [][]float64 // [i][j], 0-indexed; entries beyond Order are zero
}

SIPPoly holds a SIP polynomial as a 2D coefficient grid indexed by (i, j), where Coeffs[i][j] is the coefficient of u^i * v^j. Order is the maximum total degree (i + j <= Order).

type TNX

type TNX struct {
	Lon tnxSurface
	Lat tnxSurface
}

TNX holds a parsed TNX distortion for both axes.

func ParseTNX

func ParseTNX(h *header.Header) (*TNX, error)

ParseTNX reads TNX coefficients from a *header.Header. Returns nil if no TNX keywords are present (WAT1_* or WAT2_* absent).

func (*TNX) Forward

func (t *TNX) Forward(xi, eta float64) (xiOut, etaOut float64)

Forward applies the TNX correction at intermediate-world coordinates (xi, eta). Output is (xi + lng_correction, eta + lat_correction).

func (*TNX) Inverse

func (t *TNX) Inverse(xiPrime, etaPrime float64) (xi, eta float64, ok bool)

Inverse iterates Newton's method against Forward, using the analytic Jacobian computed via evalWithDeriv. The forward map is F(xi, eta) = (xi + Lon(xi, eta), eta + Lat(xi, eta)), so the Jacobian is

J = I + [∂Lon/∂xi  ∂Lon/∂eta]
        [∂Lat/∂xi  ∂Lat/∂eta]

type TPV

type TPV struct {
	Axis1 [40]float64
	Axis2 [40]float64
	Lon   int // 1-based longitude axis index (conventionally 1)
	Lat   int // 1-based latitude axis index (conventionally 2)
}

TPV holds a parsed TPV distortion: one coefficient array per axis.

func ParseTPV

func ParseTPV(w *wcs.Header) *TPV

ParseTPV extracts TPV coefficients from a parsed WCS header. Returns nil if no TPV coefficients are present (i.e. no PV1_* or PV2_* values).

func (*TPV) Forward

func (t *TPV) Forward(xi, eta float64) (xiOut, etaOut float64)

Forward applies the TPV polynomial to intermediate world coordinates. Input (xi, eta) is the pre-distortion position; output is the distortion-corrected position, both in degrees.

func (*TPV) Inverse

func (t *TPV) Inverse(xiPrime, etaPrime float64) (xi, eta float64, ok bool)

Inverse applies the inverse TPV correction via Newton iteration against the forward polynomial, using the analytic Jacobian. Each iteration evaluates both the residual and the 2×2 Jacobian in one pass per axis, so total cost is 2 polynomial walks per iteration (vs. 6 for the previous numerical-Jacobian version). The analytic Jacobian is also numerically stable at any distance from the origin, unlike finite differences which lose precision at the h ~ 1e-7 scale.

type Transform

type Transform struct {
	// contains filtered or unexported fields
}

Transform maps between image pixel coordinates and celestial spherical coordinates for a 2-D celestial WCS.

Higher-dimensional WCS with more than one celestial axis pair (rare in practice) is not supported.

Pixel coordinates are 0-BASED to match Go slice indexing: pixel (0, 0) is the center of the first image cell on disk, which corresponds to FITS pixel (1, 1). The library applies the +1 offset internally.

Angles in the public API (alpha, delta) are in DEGREES. Internally the pipeline works in radians.

func New

func New(h *wcs.Header) (*Transform, error)

New builds a Transform from a parsed *wcs.Header. Returns an error if the header is not a 2-D celestial WCS, if the linear transform is singular, or if the projection code is not supported (see ErrUnsupportedProjection in projection.go).

func (*Transform) Code

func (t *Transform) Code() string

Code returns the projection code in use (e.g. "TAN").

func (*Transform) PixelToSky

func (t *Transform) PixelToSky(p1, p2 float64) (alpha, delta float64, err error)

PixelToSky maps 0-based pixel coordinates (p1, p2) to celestial coordinates (alpha, delta) in degrees. alpha is normalized to [0, 360).

func (*Transform) SkyToPixel

func (t *Transform) SkyToPixel(alpha, delta float64) (p1, p2 float64, err error)

SkyToPixel is the inverse of PixelToSky. Returns 0-based pixel coordinates.