WCS Transformations

class WCS:
    """
    World Coordinate System transformation object.
    
    Parameters:
    - header: FITS header containing WCS keywords
    - fobj: file object or HDU
    - key: WCS key (single character, default ' ')
    - keysel: sequence of flags for selecting WCS keys
    - colsel: sequence of flags for selecting WCS axes
    - fix: apply standard fixes to WCS
    - translate_units: translate non-standard units
    - naxis: number of axes
    - relax: degree of WCS standard relaxation
    - hdukw: header-data unit keywords
    """
    def __init__(self, header=None, fobj=None, key=' ', keysel=None, colsel=None, 
                 fix=True, translate_units='', naxis=None, relax=True, hdukw=None): ...
    
    def pixel_to_world(self, *pixel_arrays):
        """
        Convert pixel coordinates to world coordinates.
        
        Parameters:
        - *pixel_arrays: pixel coordinate arrays
        
        Returns:
        tuple: world coordinate arrays or SkyCoord objects
        """
    
    def world_to_pixel(self, *world_arrays):
        """
        Convert world coordinates to pixel coordinates.
        
        Parameters:
        - *world_arrays: world coordinate arrays or SkyCoord objects
        
        Returns:
        tuple: pixel coordinate arrays
        """
    
    def all_pix2world(self, pixcrd, origin):
        """
        Convert pixel to world coordinates (lower-level interface).
        
        Parameters:
        - pixcrd: pixel coordinates
        - origin: pixel coordinate origin (0 or 1)
        
        Returns:
        ndarray: world coordinates
        """
    
    def all_world2pix(self, world, origin):
        """
        Convert world to pixel coordinates (lower-level interface).
        
        Parameters:
        - world: world coordinates
        - origin: pixel coordinate origin (0 or 1)
        
        Returns:
        ndarray: pixel coordinates
        """
    
    def footprint_to_file(self, filename='footprint.reg', color='green', width=2):
        """Write image footprint to DS9 region file."""
    
    def calc_footprint(self, header=None, undistort=True, axes=None, center=True):
        """Calculate image footprint in world coordinates."""
    
    @property
    def pixel_shape(self):
        """Shape of pixel grid."""
    
    @property
    def pixel_bounds(self):
        """Bounds of pixel coordinates."""
    
    @property
    def world_axis_names(self):
        """Names of world coordinate axes."""
    
    @property
    def world_axis_units(self):
        """Units of world coordinate axes."""
    
    @property
    def world_axis_physical_types(self):
        """Physical types of world coordinate axes."""
    
    @property
    def pixel_axis_names(self):
        """Names of pixel coordinate axes."""
    
    @property
    def array_shape(self):
        """Shape of associated data array."""
    
    def copy(self):
        """Create a copy of the WCS object."""
    
    def deepcopy(self):
        """Create a deep copy of the WCS object."""
    
    def to_header(self, relax=None, key=None):
        """Convert WCS to FITS header."""

def find_all_wcs(header, relax=True, keysel=None):
    """
    Find all WCS transformations in a FITS header.
    
    Parameters:
    - header: FITS header
    - relax: degree of standard relaxation
    - keysel: key selection criteria
    
    Returns:
    list: list of WCS objects
    """

def wcs_to_celestial_frame(wcs):
    """
    Convert WCS to astropy coordinate frame.
    
    Parameters:
    - wcs: WCS object
    
    Returns:
    BaseCoordinateFrame: corresponding coordinate frame
    """

def celestial_frame_to_wcs(frame, projection='TAN'):
    """
    Convert astropy coordinate frame to WCS.
    
    Parameters:
    - frame: coordinate frame
    - projection: WCS projection type
    
    Returns:
    WCS: corresponding WCS object
    """

def proj_plane_pixel_scales(wcs):
    """
    Calculate pixel scales in projection plane.
    
    Parameters:
    - wcs: WCS object
    
    Returns:
    ndarray: pixel scales
    """

def proj_plane_pixel_area(wcs):
    """
    Calculate pixel area in projection plane.
    
    Parameters:
    - wcs: WCS object
    
    Returns:
    Quantity: pixel area
    """

def skycoord_to_pixel(coords, wcs, origin=0, mode='all'):
    """
    Convert SkyCoord to pixel coordinates.
    
    Parameters:
    - coords: SkyCoord object
    - wcs: WCS object
    - origin: pixel coordinate origin
    - mode: conversion mode
    
    Returns:
    tuple: pixel coordinates
    """

def pixel_to_skycoord(xp, yp, wcs, origin=0, mode='all', cls=None):
    """
    Convert pixel coordinates to SkyCoord.
    
    Parameters:
    - xp, yp: pixel coordinates
    - wcs: WCS object
    - origin: pixel coordinate origin
    - mode: conversion mode
    - cls: SkyCoord subclass
    
    Returns:
    SkyCoord: world coordinates
    """

from astropy.wcs import WCS
from astropy.io import fits
import numpy as np

# Load WCS from FITS header
with fits.open('image.fits') as hdul:
    wcs = WCS(hdul[0].header)
    data = hdul[0].data

# Convert pixel to world coordinates
x_pix, y_pix = 100, 200
world_coords = wcs.pixel_to_world(x_pix, y_pix)
print(f"Pixel ({x_pix}, {y_pix}) -> World {world_coords}")

# Convert world to pixel coordinates
pixel_coords = wcs.world_to_pixel(world_coords)
print(f"World {world_coords} -> Pixel {pixel_coords}")

# Convert arrays of coordinates
x_pixels = np.array([10, 50, 100, 150])
y_pixels = np.array([20, 60, 120, 180])

# Pixel to world (returns SkyCoord for celestial coordinates)
world_coords = wcs.pixel_to_world(x_pixels, y_pixels)
print(f"RA: {world_coords.ra}")
print(f"Dec: {world_coords.dec}")

# Lower-level interface for more control
world_array = wcs.all_pix2world(np.column_stack([x_pixels, y_pixels]), 0)
ra_deg = world_array[:, 0]
dec_deg = world_array[:, 1]

from astropy.coordinates import SkyCoord
from astropy.wcs.utils import skycoord_to_pixel, pixel_to_skycoord
import astropy.units as u

# Create SkyCoord objects
coords = SkyCoord(ra=[10.1, 20.2, 30.3]*u.degree, 
                  dec=[40.4, 50.5, 60.6]*u.degree)

# Convert to pixel coordinates
x_pix, y_pix = skycoord_to_pixel(coords, wcs)
print(f"Pixel coordinates: x={x_pix}, y={y_pix}")

# Convert back to SkyCoord
coords_back = pixel_to_skycoord(x_pix, y_pix, wcs)
print(f"Round-trip check: {coords.separation(coords_back).to(u.arcsec)}")

# Examine WCS properties
print(f"Pixel shape: {wcs.pixel_shape}")
print(f"World axis names: {wcs.world_axis_names}")
print(f"World axis units: {wcs.world_axis_units}")
print(f"Physical types: {wcs.world_axis_physical_types}")

# Calculate pixel scales
from astropy.wcs.utils import proj_plane_pixel_scales, proj_plane_pixel_area
pixel_scales = proj_plane_pixel_scales(wcs)
pixel_area = proj_plane_pixel_area(wcs)

print(f"Pixel scales: {pixel_scales * 3600} arcsec/pixel")
print(f"Pixel area: {pixel_area.to(u.arcsec**2)}")

from astropy.coordinates import ICRS
from astropy.wcs.utils import celestial_frame_to_wcs

# Create WCS from coordinate frame
frame = ICRS()
custom_wcs = celestial_frame_to_wcs(frame, projection='TAN')

# Set reference pixel and coordinate
custom_wcs.wcs.crpix = [512, 512]  # Reference pixel
custom_wcs.wcs.crval = [180.0, 45.0]  # Reference coordinate (RA, Dec)
custom_wcs.wcs.cdelt = [-0.1/3600, 0.1/3600]  # Pixel scale (degrees)
custom_wcs.wcs.ctype = ['RA---TAN', 'DEC--TAN']  # Coordinate types

# Use custom WCS
test_coord = custom_wcs.pixel_to_world(512, 512)
print(f"Reference coordinate: {test_coord}")