WCS not appearing in generated FITS files

Question

I am trying to generate a FITS file with random stars. For some reason, when I open the file in SAO DS9, the mouse-over doesn't display the on-sky coordinates, but instead only shows the pixel locations to which I am mousing over. Below I'm showing the script I am using in Python to generate the FITS file, and what I am seeing in DS9:

import numpy as np
from pathlib import Path
import matplotlib.pyplot as plt
import pandas as pd
import os
from astropy.wcs import WCS
from astropy import units as u
from astropy.io import fits
from astropy.visualization import AsinhStretch, LogStretch
from astropy.visualization.mpl_normalize import ImageNormalize
from convenience_functions import show_image
from photutils.datasets import make_gaussian_sources_image, make_random_gaussians_table
from astropy.nddata import Cutout2D



def add_wcs(cntrpix,cntrval,cdelt,ctype, radesys):
    #zhdrbb.remove('ZODIL')
    #CREATE WCS KEYWORDS
    w = WCS(naxis=2)
    w.wcs.crpix = cntrpix
    w.wcs.crval = cntrval
    w.wcs.cdelt = cdelt
    w.wcs.ctype = ctype
    w.wcs.radesys = radesys
    w.wcs.equinox = 2000
    w.wcs.set_pv([(2, 1, 45.0)])
    w.wcs.cunit = [u.deg, u.deg]

    header = w.to_header()
    return header

def stars(image, number, mag_min, mag_max):
    """
    Add some stars to the image.
    """

    # Most of the code below is a direct copy/paste from
    # https://photutils.readthedocs.io/en/stable/_modules/photutils/datasets/make.html#make_100gaussians_image

    flux_min = 10**mag_min
    flux_max = 10**mag_max
    flux_range = [flux_min, flux_max]

    y_max, x_max = image.shape
    xmean_range = [0.01 * x_max, 0.99 * x_max]
    ymean_range = [0.01 * y_max, 0.99 * y_max]
    xstddev_range = [4, 4]
    ystddev_range = [4, 4]
    params = dict([('flux', flux_range),
                  ('x_mean', xmean_range),
                  ('y_mean', ymean_range),
                  ('x_stddev', xstddev_range),
                  ('y_stddev', ystddev_range),
                  ('theta', [0, 2*np.pi])])
    
    sources = make_random_gaussians_table(number, params, seed=int(random.random()*1e4))

    star_im = make_gaussian_sources_image(image.shape, sources)

    return sources, star_im

RA_0, RA_F = 30, 70
DEC_0, DEC_F = -20,20
MAG_0, MAG_F = 3, 8
N_STARS = 250

D_RA = np.abs(RA_0 - RA_F)
D_DEC = np.abs(DEC_0 - DEC_F)


GRID_LENGTH = int(2e3) # N PIXELS PER SIDE

# WCS STUFF
CRPIX = (int(GRID_LENGTH / 2), int(GRID_LENGTH / 2))
CRVAL = (RA_0 + D_RA / 2., DEC_0 + D_DEC / 2.)

cdelt_x = D_RA / GRID_LENGTH
cdelt_y = D_DEC / GRID_LENGTH
CDELT = (cdelt_x, cdelt_y)
CTYPE = ["RA-TAN", "DEC-TAN"]
RADESYS = 'ICRS'
base_pix_full_image = np.zeros((GRID_LENGTH,GRID_LENGTH))
star_data, star_cat_img = stars(base_pix_full_image, N_STARS, MAG_0, MAG_F)
# WCS
wcs_header = add_wcs(CRPIX, CRVAL, CDELT, CTYPE, RADESYS)

hdr = fits.Header()
hdu_cat = fits.PrimaryHDU(header=wcs_header,data=star_cat_img)

# WRITE FITS PIXEL FILE
hdu_cat.writeto(fits_pixel_file, overwrite=True)

From the image, the mouseover on the image shows x and y pointing in pixel coordinates, even though the header has WCS information in it.

I feel like I'm missing something very simple and just need a fresh set of eyes to troubleshoot.

Answer 1

It seems that you need to set CTYPE as CTYPE = ["RA---TAN", "DEC--TAN"].

Here is a cite from FITS Standard.

CTYPEi – [string; indexed; default: '␣' (i.e. a linear, undefined axis)]. Type for the Intermediate-coordinate Axis i. Any coordinate type that is not covered by this Standard or an officially recognized FITS convention shall be taken to be linear. All non-linear coordinate system names must be expressed in ‘4–3’ form: the first four characters specify the coordinate type, the fifth character is a hyphen (‘-’), and the remaining three characters specify an algorithm code for computing the world coordinate value. Coordinate types with names of fewer than four characters are padded on the right with hyphens, and algorithm codes with fewer than three characters are padded on the right with blanks11. Algorithm codes should be three characters.