UPDATE(03/10/2016): Thought I would put this here for anybody in the future to see. I have essentially put this problem to rest. I was essentially trying to "make" a world coordinate system after the image was taken. My thought was, if I could get a WCS on hundreds of images from a particular night reduction would be much easier. However, something like this is not so easily done without reference points in the image or just making one at the source when the image is taken. As of writing this my University's telescope is being upgraded to have an auto-guiding system so we should hopefully be able to implement a WCS in our image headers this way.
Hopefully there isn't to much assumed knowledge. This is purely an astronomy reduction question using astropy.
Some Background: Right now I have a whole bunch of fits images; hundreds for many different targets. The telescope I used has some nasty drift if you just sit on the target and so now I need to align all my images. I could use non-free software like PixInsight that I have heard does a beautiful job, but I am a curious undergraduate who wants to try it the free way/hard way. After extensively looking for a way to align many images with drift in IRAF it seems what I need most are star coordinates. So I came up with the idea that if I get the world coordinates for a reference star and then use that against all the other images I have I can extract pixel coordinates for the reference star in each image. This will allow me to calculate shifts and thus align the data.
On to my troubles:
Following the documentation on getting the wcs from a fits header with astropy I come into the problem following error:
WARNING: FITSFixedWarning: EPOCH = '2015.5 ' / Epoch for RA and Dec (years) A floating-point value was expected. [astropy.wcs.wcs]
The epoch is strange because a telescope can't slew using the EPOCH of 2000 so the system converts it to the current year to properly transform RA and DEC and that is what is cataloged in the header, I believe. If I press on from here and try to convert pixel coordinates to world coordinates then I end up with the pixel coordinates again so clearly this isn't working. My thoughts then are to manually create a WCS object for each image but I can't really follow the astropy documentation so well on this. The attributes for the wcs object confuse me when all I seem to know is the RA, and Dec from the header.
Could somebody help me sort this out?
Don't tell me to use imstar, please, either because that requires me to install wfctools to IRAF and everybody knows IRAF is a finicky beast which I would likely break.
Thanks!
PS: Hopefully somebody is still willing to help me here, but as to my excuse summer quarter ended and this question went to the side for vacation. However, I want to revive it and solve it as well as rephrase it.
Essentially, my University has a telescope that only undergrads really use because it is bottom of research grade for today's standards. So it goes on to be a good exercise for Astronomy undergrads to drive a 30 inch telescope. Ultimately, the fits headers on the images do not include a world coordinate system (WCS) and I believe this to be a problem as I feel it would be a lot easier to align images, in mass, if all the images had a WCS. So what I am trying to do is add a WCS correctly to the image headers which I don't have much of an idea on how to do.
My first attempts have been to just test making a WCS within astropy, which I will provide some preliminary code, but they have been fraught with failure because the key words to making a WCS confuse me. I have found what seems to be a good description and possibly tutorial on WCSs here: http://astroweb.iag.usp.br/~moser/notes/GAi_FITSimgs.html.
Here is my code for trying make a WCS (it is a bit pathetic):
import numpy as np
from astropy.io import fits
from astropy import wcs
def main():
w = wcs.WCS(naxis=2)
w.wcs.crpix = [256, 256] # center pixel
w.wcs.crval = [18.73337, 0.54708] # RA and dec values in hours and degrees
w.wcs.cdelt = np.array([-1.0, 1.0])
w.wcs.ctype = ["RA-TAN", "DEC-TAN"]
#print w.wcs.print_contents()
wx, wy = w.wcs_pix2world(232.85, 278.14, 1) # pixel coords of interested star
print wx, wy
pix_x, pix_y = w.wcs_world2pix(wx, wy, 1)
# for gstandard010; another images WCS
w2 = wcs.WCS(naxis=2)
w2.wcs.crpix = [256, 256]
w2.wcs.crval = [18.733427777777777, 0.5471388888888888]
w2.wcs.cdelt = np.array([-1.0, 1.0])
w2.wcs.ctype = ["RA-TAN", "DEC-TAN"]
wx2, wy2 = w.wcs_pix2world(pix_x, pix_y, 1)
print wx2, wy2
# Using first WCS wx, and wy to find the pixel coordinates in this new images.
# It only outputs the first images pixel values which is incorrect.
pix_x2, pix_y2 = w2.wcs_world2pix(wx, wy, 1)
print pix_x2, pix_y2
main()
OutPut:
41.88337 22.68708
41.88337 22.68708
232.850057778 278.139941111
I believe I need to further fill out the details of the WCS parameters, including ones that are missing, so could somebody help me figure this out?