Today I am looking further into the ra/dec offset between the catalog and SEP extracted sources.
%matplotlib inline
import numpy as np
import sep
from astropy.io import fits
import matplotlib.pyplot as plt
from matplotlib.patches import Ellipse
from astropy import wcs
from math import *
from matplotlib.patches import Ellipse
# returns the indices of a list in descending order
def sort_reverse_index(x):
return sorted(range(len(x)), key=lambda k: x.max()-x[k])
#returns indices of list in ascending order
def sort_index(x):
return sorted(range(len(x)), key=lambda k: x[k])
#imports fits and puts into numpy array
def load_fits(filename):
hdu_list = fits.open(filename, do_not_scale_image_data=True)
tbdat = hdu_list[0].data
tbdat = tbdat.byteswap().newbyteorder()
return tbdat
#cuts the exposure images to 2000X2000, as the mosaic is just 2000X2000
def cut_exp(tbdat):
tbdat = np.delete(tbdat, list(range(2000, 2036)),0)
tbdat = np.delete(tbdat, list(range(2000, 2036)),1)
return tbdat
#function gives zero point from exposure 1. All exposures should have the same zero point.
def zero(filename):
hdu_list = fits.open(filename, do_not_scale_image_data=True)
zero_point = hdu_list[0].header['ABMAG']
return zero_point
#subtracts background from image
def subtraction(tbdat, bkg):
bkg_rms = bkg.rms() #background noise as 2d array
tbdat_sub = tbdat - bkg #subtract background
return tbdat_sub
#perform the extraction
def extraction(tbdat_sub, bkg):
objects = sep.extract(tbdat_sub, thresh = 5.0, err = bkg.globalrms)
return objects
#Gets kron radius information
def kron_info(objects, tbdat_sub):
kronrad, kronflag = sep.kron_radius(tbdat_sub, objects['x'], objects['y'], objects['a'], objects['b'], objects['theta'], r=6.0)
return kronrad, kronflag
#Finds flux, flux error, and extraction flags
def find_flux(tbdat_sub, bkg, objects, kronrad, kronflag):
flux, fluxerr, flag = sep.sum_ellipse(tbdat_sub, objects['x'], objects['y'], objects['a'], objects['b'], objects['theta'], r = (2.5*kronrad), err = bkg.globalrms, subpix=1)
flag |=kronflag #combines all flags
r_min = 1.75 #minimum diameter = 3.5
use_circle = kronrad * np.sqrt(objects['a'] * objects['b']) < r_min
cflux, cfluxerr, cflag = sep.sum_circle(tbdat_sub, objects['x'][use_circle], objects['y'][use_circle], r_min, subpix=1)
flux[use_circle] = cflux
fluxerr[use_circle] = cfluxerr
flag[use_circle] = cflag
r, rflag = sep.flux_radius(tbdat_sub, objects['x'], objects['y'], rmax = 6.0*objects['a'], frac = 0.5, normflux = flux, subpix =5)
sig = 2.0 / (2.35*r) # r from sep.flux_radius() above, with fluxfrac = 0.5
xwin, ywin, wflag = sep.winpos(tbdat_sub, objects['x'], objects['y'], sig)
return flux, fluxerr, flag, r, xwin, ywin
#convert extraction flux to AB magnitude
def magnitude(flux, zero_point):
mag = -2.5*np.log10(flux) + zero_point
return mag
#reads in catalog used to generated exposures and puts info into arrays
def read_cat(filename):
f = open(filename, "r")
line = f.readlines()[1:]
f.close()
num_cat = np.array([])
ra_cat = np.array([])
dec_cat = np.array([])
mag_cat = np.array([])
z_cat = np.array([])
s_cat = np.array([])
a_cat = np.array([])
b_cat = np.array([])
theta_cat = np.array([])
for i in range(len(line)):
num_cat = np.append(num_cat, int(line[i].split()[0]))
ra_cat = np.append(ra_cat, float(line[i].split()[1]))
dec_cat = np.append(dec_cat, float(line[i].split()[2]))
mag_cat = np.append(mag_cat, float(line[i].split()[3]))
z_cat = np.append(z_cat, float(line[i].split()[4]))
a_cat = np.append(a_cat, float(line[i].split()[5]))
b_cat = np.append(b_cat, float(line[i].split()[6]))
theta_cat = np.append(theta_cat, float(line[i].split()[7]))
s_cat = np.append(s_cat, float(line[i].split()[9]))
return num_cat, ra_cat, dec_cat, mag_cat, z_cat, s_cat, a_cat, b_cat, theta_cat
#convert pixel location of objectes extracted from image to ra and dec
def world(filename, objects):
hdu_list = fits.open(filename)
w = wcs.WCS(hdu_list[0].header)
hdu_list.close()
wrd = w.all_pix2world(objects['x'], objects['y'], 0)
ra = wrd[:][0]
dec = wrd[:][1]
return ra, dec
#convert xwin and ywin to ra and dec
def world_win(filename, xwin, ywin):
hdu_list = fits.open(filename)
w = wcs.WCS(hdu_list[0].header)
hdu_list.close()
wrd = w.all_pix2world(xwin, ywin, 0)
ra = wrd[:][0]
dec = wrd[:][1]
return ra, dec
#convert pixel location of corners of image to ra and dec
def world_area(filename, tbdat):
hdu_list = fits.open(filename)
w = wcs.WCS(hdu_list[0].header)
hdu_list.close()
world_min = w.all_pix2world(0,0,0)
world_max = w.all_pix2world(len(tbdat), len(tbdat), 0)
return world_min, world_max
fname_exp1 = 'simple_sim_cube_F090W_487_01.slp.fits'
tbdat_exp1 = load_fits(fname_exp1)
tbdat_exp1 = cut_exp(tbdat_exp1)
tbdat_exp1 = np.ascontiguousarray(tbdat_exp1) #Converts cut array to C-order memory
zero_point = zero(fname_exp1)
fname_img = 'F090W_8exp.fits'
tbdat = load_fits(fname_img)
tbdat = np.nan_to_num(tbdat)
for i in range(len(tbdat)):
for j in range(len(tbdat)):
if tbdat[i][j] == 0:
tbdat[i][j] = 0.00001
bkg = sep.Background(tbdat) #measures background
tbdat_sub = subtraction(tbdat, bkg)
objects = extraction(tbdat_sub, bkg) #extraction of mosaic image
print('number of mosaic sources', len(objects['x']))
kronrad, kronflag = kron_info(objects, tbdat_sub) #Finds Kron Radius
flux, fluxerr, flag, r, xwin, ywin= find_flux(tbdat_sub, bkg, objects, kronrad, kronflag) #Finds flux from aperature phot.
mag = magnitude(flux, zero_point)
('number of mosaic sources', 271)
Next I'll get all of the info from the catalog
num_cat, ra_cat, dec_cat, mag_cat, z, s, a_cat, b_cat, theta_cat= read_cat('candels_with_fake_mag.cat')
I'll get RA/DEC info on all of the sources extracted from the mosaic.
ra, dec = world(fname_img, objects)
Now I want to find the RA and DEC of the total area of the mosaic. I'll use this to plot the extracted sources over the mosaic image.
world_min, world_max = world_area(fname_img, tbdat)
Delete all catalog values outside of image area.
print "mosaic minimum RA/DEC", world_min[0], world_min[1]
print "mosaic maximum RA/DEC", world_max[0], world_max[1]
print "minimum RA/DEC in catalog", min(ra_cat), min(dec_cat)
print "maximum RA/DEC in catalog", max(ra_cat), max(dec_cat)
mosaic minimum RA/DEC 53.1014876537 -27.8135469046 mosaic maximum RA/DEC 53.1292161985 -27.7895521325 minimum RA/DEC in catalog 52.992519 -27.959661 maximum RA/DEC in catalog 53.270374 -27.656631
ra_cat_new = np.asarray([])
dec_cat_new = np.asarray([])
mag_cat_new = np.asarray([])
a_cat_new = np.asarray([])
b_cat_new = np.asarray([])
theta_cat_new = np.asarray([])
for i in range(len(ra_cat)):
if (ra_cat[i] >= world_min[0]) and (ra_cat[i] <= world_max[0]):
ra_cat_new = np.append(ra_cat_new, ra_cat[i])
dec_cat_new = np.append(dec_cat_new, dec_cat[i])
mag_cat_new = np.append(mag_cat_new, mag_cat[i])
a_cat_new = np.append(a_cat_new, a_cat[i])
b_cat_new = np.append(b_cat_new, b_cat[i])
theta_cat_new = np.append(theta_cat_new, theta_cat[i])
ra_cat_cut = np.asarray([])
dec_cat_cut = np.asarray([])
mag_cat_cut = np.asarray([])
a_cat_cut = np.asarray([])
b_cat_cut = np.asarray([])
theta_cat_cut = np.asarray([])
for j in range(len(ra_cat_new)):
if (dec_cat_new[j] >= world_min[1]) and (dec_cat_new[j] <= world_max[1]):
ra_cat_cut = np.append(ra_cat_cut, ra_cat_new[j])
dec_cat_cut = np.append(dec_cat_cut, dec_cat_new[j])
mag_cat_cut = np.append(mag_cat_cut, mag_cat_new[j])
a_cat_cut = np.append(a_cat_cut, a_cat_new[j])
b_cat_cut = np.append(b_cat_cut, b_cat_new[j])
theta_cat_cut = np.append(theta_cat_cut, theta_cat_new[j])
Let's plot the mosaic image and the extracted sources over.
fig=plt.figure(figsize=(16, 16))
plt.imshow(tbdat, vmin=0.093537331, vmax=0.21167476, cmap='gray' ,interpolation='none', origin='lower', extent = [world_min[0], world_max[0], world_min[1], world_max[1]])
plt.plot(ra, dec, '.', color = 'red')
[<matplotlib.lines.Line2D at 0x10fcf7750>]
The ra/dec offsets are seen in the above image.
To see some numbers let's list the 25 brightest sources given in the catalog.
#prints souces in magnitude increasing order with corresponding ra/dec (highest 50 flux)
index_cat_cut = sort_index(mag_cat_cut)
print "ABMag", "\t","RA", "\t\t", "DEC"
for j in range(0,25):
k = index_cat_cut[j]
print mag_cat_cut[k], "\t", ra_cat_cut[k], "\t", dec_cat_cut[k]
ABMag RA DEC 17.962 53.118019 -27.802065 18.292 53.103878 -27.804071 19.217 53.101707 -27.80225 19.599 53.124931 -27.810972 20.93 53.10648 -27.800976 20.968 53.121548 -27.809208 21.017 53.11866 -27.807377 21.098 53.124756 -27.799213 21.107 53.111118 -27.799646 21.188 53.125084 -27.790775 21.712 53.104717 -27.790092 21.773 53.112545 -27.80802 21.824 53.119873 -27.798738 22.123 53.123108 -27.803391 22.152 53.109596 -27.78956 22.346 53.116508 -27.806746 22.449 53.105938 -27.808401 22.591 53.127308 -27.810123 22.652 53.107235 -27.805758 22.655 53.103172 -27.806795 22.716 53.114033 -27.811705 22.773 53.111141 -27.809668 22.815 53.128803 -27.80895 22.816 53.120098 -27.808273 22.82 53.109676 -27.792929
List the brightes 25 extracted sources and their positions.
#prints souces in magnitude increasing order with corresponding ra/dec (highest 50 flux)
index = sort_index(mag)
print "ABMag", "\t\t","RA", "\t\t", "DEC"
for j in range(0,25):
k = index[j]
print mag[k], "\t", ra[k], "\t", dec[k]
ABMag RA DEC 18.4270536033 53.1179584593 -27.8021192292 18.662848017 53.103818779 -27.8041255775 20.9791018403 53.1017108602 -27.802308584 21.0591447351 53.1094483908 -27.7896864632 21.0729577287 53.118601134 -27.8074302089 21.1140913907 53.1064200452 -27.8010299362 21.1641615272 53.1246959748 -27.799265158 21.1704690828 53.1110581308 -27.7996995012 21.3634907376 53.1250219133 -27.7908270532 21.5199311699 53.1201394124 -27.7989000524 21.6837848695 53.1198166833 -27.7987918466 21.9223937575 53.1124861866 -27.8080731775 22.2948675839 53.1230501023 -27.8034433825 22.5019814114 53.1164501547 -27.8067990094 22.5103897294 53.1277223527 -27.7951315176 22.5928427253 53.1058804401 -27.8084546182 22.7126199969 53.1031140496 -27.8068486871 22.760590268 53.1071729434 -27.8058111373 22.7653230682 53.1277855597 -27.7948490914 22.782749451 53.1139740435 -27.8117594326 22.8232764127 53.1095840799 -27.7929856005 22.86717375 53.1200406331 -27.8083256515 22.9100785809 53.108784454 -27.811773562 22.9471532065 53.1117553728 -27.8052971082 23.0152936755 53.1110833401 -27.8097215357
Below I am plotting the pixel locations of the sources extracted from SEP on top of my mosaic image.
fig=plt.figure(figsize=(16, 16))
plt.imshow(tbdat, vmin=0.093537331, vmax=0.21167476, cmap='gray' ,interpolation='none', origin='lower')
for i in range(0,100):
k = index[i]
plt.plot(objects[k]['x'], objects[k]['y'], '.', color = 'red')
So this is interesting. The pixel locations are pretty spot on. Something is going awry when converting these pixel locations to sky coordinates.
I am curious to see if my header information for my mosaic image is not correct. I'll perform source extraction on exposure 1 to see if the offset is still there. I know the header information from exposure one should be correct (since they are given by Christopher who generated the exposures).
fname_exp1 = 'simple_sim_cube_F090W_487_01.slp.fits'
tbdat_exp1 = load_fits(fname_exp1)
bkg_exp1 = sep.Background(tbdat_exp1) #measures background
tbdat_exp1_sub = subtraction(tbdat_exp1, bkg_exp1)
objects_exp1 = extraction(tbdat_exp1_sub, bkg_exp1) #extraction of mosaic image
print('number of exp1 sources', len(objects_exp1['x']))
kronrad_exp1, kronflag_exp1 = kron_info(objects_exp1, tbdat_exp1_sub) #Finds Kron Radius
flux_exp1, fluxerr_exp1, flag_exp1, r_exp1, xwin_exp1, ywin_exp1= find_flux(tbdat_exp1_sub, bkg_exp1, objects_exp1, kronrad_exp1, kronflag_exp1) #Finds flux from aperature phot.
mag_exp1 = magnitude(flux_exp1, zero_point)
('number of exp1 sources', 65)
Below is a plot of the pixel locations over exposure 1.
index_exp1 = sort_index(mag_exp1)
fig=plt.figure(figsize=(16, 16))
plt.imshow(tbdat_exp1, vmin=0.093537331, vmax=0.21167476, cmap='gray' ,interpolation='none', origin='lower')
for i in range(0,50):
k = index_exp1[i]
plt.plot(objects_exp1[k]['x'], objects_exp1[k]['y'], '.', color = 'red')
All seems well. I'm going to use astropy wcs to convert to sky coordinates and plot.
ra_exp1, dec_exp1 = world(fname_exp1, objects_exp1)
world_min_exp1, world_max_exp1 = world_area(fname_exp1, tbdat_exp1)
fig=plt.figure(figsize=(16, 16))
plt.imshow(tbdat_exp1, vmin=0.093537331, vmax=0.21167476, cmap='gray' ,interpolation='none', origin='lower', extent = [world_min_exp1[0], world_max_exp1[0], world_min_exp1[1], world_max_exp1[1]])
for i in range(0,50):
k = index_exp1[i]
plt.plot(ra_exp1[k], dec_exp1[k], '.', color = 'red')
There is still an offset. I believe that the cause is not from the fits header information.
I'm going to try to use xwin and ywin positions to see if that changes things. These locations are found in SEP and are a more accurate location that is calculated by weighing pixels by their flux and a 2d Gaussian amplitude within a circular aperature. See http://sep.readthedocs.io/en/v1.0.x/api/sep.winpos.html?highlight=xwin.
First I'll overplot the xwin, ywin on the image.
fig=plt.figure(figsize=(16, 16))
plt.imshow(tbdat, vmin=0.093537331, vmax=0.21167476, cmap='gray' ,interpolation='none', origin='lower')
for i in range(0,100):
k = index[i]
plt.plot(xwin[k], ywin[k], '.', color = 'red')
All looks well. Now I'll find the RA/DEC from xwin,ywin and plot.
ra_win, dec_win = world_win(fname_img, xwin, ywin)
fig=plt.figure(figsize=(16, 16))
plt.imshow(tbdat, vmin=0.093537331, vmax=0.21167476, cmap='gray' ,interpolation='none', origin='lower', extent = [world_min[0], world_max[0], world_min[1], world_max[1]])
for i in range(0,50):
k = index[i]
plt.plot(ra_win[k], dec_win[k], '.', color = 'red')
Still distorted.
Brant and I discussed how it might be a SIP distortion. However, if this was the case, this information would be given in the header. For example:
(Header example from The SIP Convention for Representing Distortion in FITS Image Headers, D. L. Shupe and R. N. Hook)
CTYPE1 = 'RA---TAN-SIP' / RA---TAN with distortion in pixel space
CTYPE2 = 'DEC--TAN-SIP' / DEC--TAN with distortion in pixel space
CRVAL1 = 202.581507417836 / [deg] RA at CRPIX1,CRPIX2 (using Pointing Recon
CRVAL2 = 47.2465528124827 / [deg] DEC at CRPIX1,CRPIX2 (using Pointing Reco
CRPIX1 = 128. / Reference pixel along axis 1
CRPIX2 = 128. / Reference pixel along axis 2
CD1_1 = 0.000248349650353678 / Corrected CD matrix element with Pointing Recon
CD1_2 = 0.000232107213140475 / Corrected CD matrix element with Pointing Recon
CD2_1 = 0.000232418393583541 / Corrected CD matrix element with Pointing Recon
CD2_2 = -0.000246562617306562 / Corrected CD matrix element with Pointing Reco
A_ORDER = 3 / polynomial order, axis 1, detector to sky
A_0_2 = 9.0886E-06 / distortion coefficient
A_0_3 = 4.8066E-09 / distortion coefficient
A_1_1 = 4.8146E-05 / distortion coefficient
A_1_2 = -1.7096E-07 / distortion coefficient
A_2_0 = 2.82E-05 / distortion coefficient
A_2_1 = 3.3336E-08 / distortion coefficient
I feel like this could still be an issue as some sources are shifted to the left and some to the right. SIP distortion corrections will perform transformations to correct the tilting of a detector. I can force astropy to just do the SIP conversion from pixel location to focal plane location.
hdu_list = fits.open(fname_img)
w = wcs.WCS(hdu_list[0].header)
hdu_list.close()
sip = w.sip_pix2foc(objects['x'],objects['y'], 0)
x_sip = sip[:][0]
y_sip = sip[:][1]
Now I'll throw this focal pixel locations through a pix to world transformation.
hdu_list = fits.open(fname_img)
w = wcs.WCS(hdu_list[0].header)
hdu_list.close()
wrd = w.wcs_pix2world(x_sip,y_sip, 0)
ra_sip = wrd[:][0]
dec_sip = wrd[:][0]
for i in range(0,25):
print 'x', '\t\t', 'x \w SIP correction', '\t', 'y', '\t\t', 'y \w SIP correction'
print objects[i]['x'], '\t', x_sip[i], '\t', objects[i]['y'], '\t', y_sip[i]
x x \w SIP correction y y \w SIP correction 813.865785428 813.865785428 8.86240601893 8.86240601893 x x \w SIP correction y y \w SIP correction 193.480208381 193.480208381 27.9025023437 27.9025023437 x x \w SIP correction y y \w SIP correction 2027.06841865 2027.06841865 25.3599877547 25.3599877547 x x \w SIP correction y y \w SIP correction 5.48527688048 5.48527688048 34.8609121738 34.8609121738 x x \w SIP correction y y \w SIP correction 1354.81191407 1354.81191407 57.1678745183 57.1678745183 x x \w SIP correction y y \w SIP correction 2027.84540692 2027.84540692 61.7920760011 61.7920760011 x x \w SIP correction y y \w SIP correction 257.647912315 257.647912315 89.2671436756 89.2671436756 x x \w SIP correction y y \w SIP correction 1301.31350226 1301.31350226 93.5109876108 93.5109876108 x x \w SIP correction y y \w SIP correction 1352.01901889 1352.01901889 103.996218974 103.996218974 x x \w SIP correction y y \w SIP correction 1716.45122611 1716.45122611 111.048613825 111.048613825 x x \w SIP correction y y \w SIP correction 1239.72562891 1239.72562891 105.374944283 105.374944283 x x \w SIP correction y y \w SIP correction 1509.32478718 1509.32478718 116.060636236 116.060636236 x x \w SIP correction y y \w SIP correction 928.210170443 928.210170443 119.919171367 119.919171367 x x \w SIP correction y y \w SIP correction 1612.41900134 1612.41900134 124.29068092 124.29068092 x x \w SIP correction y y \w SIP correction 23.577588416 23.577588416 141.155180883 141.155180883 x x \w SIP correction y y \w SIP correction 1617.77640045 1617.77640045 171.495843865 171.495843865 x x \w SIP correction y y \w SIP correction 2027.23218264 2027.23218264 95.4716705845 95.4716705845 x x \w SIP correction y y \w SIP correction 2027.22340636 2027.22340636 167.669867658 167.669867658 x x \w SIP correction y y \w SIP correction 2027.09397788 2027.09397788 130.719057249 130.719057249 x x \w SIP correction y y \w SIP correction 993.944234165 993.944234165 182.860864412 182.860864412 x x \w SIP correction y y \w SIP correction 1792.42722707 1792.42722707 204.633053362 204.633053362 x x \w SIP correction y y \w SIP correction 775.285934263 775.285934263 202.013653097 202.013653097 x x \w SIP correction y y \w SIP correction 743.344498555 743.344498555 214.295947954 214.295947954 x x \w SIP correction y y \w SIP correction 1641.51115726 1641.51115726 226.227646072 226.227646072 x x \w SIP correction y y \w SIP correction 1862.60646036 1862.60646036 225.126132841 225.126132841
for j in range(0,25):
print 'RA', '\t\t', 'RA \w SIP correction', '\t', 'DEC', '\t\t', 'DEC \w SIP correction'
print ra[j], '\t', ra_sip[j], '\t\t', dec[j], '\t', dec_sip[j]
RA RA \w SIP correction DEC DEC \w SIP correction 53.1094452644 53.1094452644 -27.8135696054 53.1094452644 RA RA \w SIP correction DEC DEC \w SIP correction 53.1033830185 53.1033830185 -27.8133279029 53.1033830185 RA RA \w SIP correction DEC DEC \w SIP correction 53.1213078705 53.1213078705 -27.8135740015 53.1213078705 RA RA \w SIP correction DEC DEC \w SIP correction 53.1015461527 53.1015461527 -27.8132442652 53.1015461527 RA RA \w SIP correction DEC DEC \w SIP correction 53.1147402787 53.1147402787 -27.8132154019 53.1147402787 RA RA \w SIP correction DEC DEC \w SIP correction 53.1213205032 53.1213205032 -27.8132571112 53.1213205032 RA RA \w SIP correction DEC DEC \w SIP correction 53.1040188911 53.1040188911 -27.812801866 53.1040188911 RA RA \w SIP correction DEC DEC \w SIP correction 53.114222323 53.114222323 -27.8128926647 53.114222323 RA RA \w SIP correction DEC DEC \w SIP correction 53.1147194742 53.1147194742 -27.8128076221 53.1147194742 RA RA \w SIP correction DEC DEC \w SIP correction 53.1182831309 53.1182831309 -27.812790671 53.1182831309 RA RA \w SIP correction DEC DEC \w SIP correction 53.1136218889 53.1136218889 -27.8127819238 53.1136218889 RA RA \w SIP correction DEC DEC \w SIP correction 53.1162589638 53.1162589638 -27.8127218338 53.1162589638 RA RA \w SIP correction DEC DEC \w SIP correction 53.1105785484 53.1105785484 -27.8126173213 53.1105785484 RA RA \w SIP correction DEC DEC \w SIP correction 53.1172679498 53.1172679498 -27.8126627878 53.1172679498 RA RA \w SIP correction DEC DEC \w SIP correction 53.1017378851 53.1017378851 -27.8123216565 53.1017378851 RA RA \w SIP correction DEC DEC \w SIP correction 53.1173268648 53.1173268648 -27.8122527226 53.1173268648 RA RA \w SIP correction DEC DEC \w SIP correction 53.1213191646 53.1213191646 -27.8129640007 53.1213191646 RA RA \w SIP correction DEC DEC \w SIP correction 53.1213290603 53.1213290603 -27.8123358248 53.1213290603 RA RA \w SIP correction DEC DEC \w SIP correction 53.1213226865 53.1213226865 -27.812657307 53.1213226865 RA RA \w SIP correction DEC DEC \w SIP correction 53.111229914 53.111229914 -27.8120777215 53.111229914 RA RA \w SIP correction DEC DEC \w SIP correction 53.1190388256 53.1190388256 -27.811985668 53.1190388256 RA RA \w SIP correction DEC DEC \w SIP correction 53.1090949993 53.1090949993 -27.8118843358 53.1090949993 RA RA \w SIP correction DEC DEC \w SIP correction 53.108784454 53.108784454 -27.811773562 53.108784454 RA RA \w SIP correction DEC DEC \w SIP correction 53.1175664772 53.1175664772 -27.8117794086 53.1175664772 RA RA \w SIP correction DEC DEC \w SIP correction 53.1197277267 53.1197277267 -27.8118159017 53.1197277267
Well just performing the SIP transformation does not change the pixel location and thus sky coordinates. I'm not sure why I expected it to as I had previously tried using all_pix2world() that includes this transformation.