Research Blog
Welcome to my Research Blog.
This is mostly meant to document what I am working on for myself, and to communicate with my colleagues. It is likely filled with errors!
This project is maintained by ndrakos
These notes are on the (FITS) World Coordinate Systems standard. An overview of representing celestial coordinate systems is given here.
Each galaxy has an RA and Dec centred on (0,0) spanning from -0.5 to 0.5 degrees. There are 0.11 arcsec per pixel.
I followed this example.
For our problem:
deg_per_pix = 0.11/3600
num_pixels = 1.0/deg_per_pix
w = wcs.WCS(naxis=2)
w.wcs.ctype = ["RA---TAN", "DEC--TAN"] # projection type
w.wcs.crpix = [0,0] # (X,Y) reference pixel
w.wcs.crval = [-0.5,-0.5] # Sky coord of reference pixel in degrees
w.wcs.cdelt = np.array([deg_per_pix, deg_per_pix]) # increment (X,Y) in degrees
Start with some pixel coordinates: pixcrd = np.array([[0, 0], [num_pixels,num_pixels], [num_pixels/2, num_pixels/2], [0,num_pixels]])
.
To get the corresponding world coordinates: world = w.all_pix2world(pixcrd, 0)
array([[-4.99969443e-01, -4.99969444e-01],
[-3.59500185e+02, 4.99852916e-01],
[-1.17812552e-06, 1.78607360e-05],
[-4.99969448e-01, 4.99929026e-01]])
This all looks fine. Note that -359.5 is the same as 0.5 degrees. Then, to convert back to pixels, we can use pixcrd = w.all_world2pix(world, 0)
.
When writing the FITS file, we can include the WCS as follows:
header = w.to_header()
hdu = fits.PrimaryHDU(header=header)