#!/usr/bin/pyhton
# get_photometric.py
# Jure Skvarc, May 2000
# For the given list of images, get all photometric standards stars,
# read star lists and find matching stars
# Show those stars on screen

import stars
from FITS import *
import fitsblink
from astrotrig import *
import re
import os
from string import *
import fileinput

# Access PHOTOM star catalog
phot = stars.Catalog('PHOT', '/home/jure/astro/tass/disk15/diary')


#  Initialize fitsblink
# Fitsblink must be started with option -s to know it is in script mode
#arguments = ["fitsblink", "-s"] #  Define the command line tokens
#argc = len(arguments)          
#args = fitsblink.string_array(argc + 1) #  Create an array of chars
#  Put arguments into array
#for i in range(argc):
#    fitsblink.string_set(args, i, arguments[i])

#fitsblink.string_set(args, argc, "")
# Start fitsblink and initialize graphic interface
#fitsblink.main0(argc, args)


imagenames = '^H.*fts'
imagere = re.compile(imagenames)
files = os.listdir('.')
for name in files:
    if imagere.search(name) != None:
    	if os.path.isfile(name):
	    print name
	    #fitsblink.load_image(name, 0)
	    #  Read keywords
	    keywords = getFITSkeywords(name)
	    width = float(keywords['CDELT1'][0]) * float(keywords['NAXIS1'][0])
	    height = float(keywords['CDELT2'][0]) * float(keywords['NAXIS2'][0])
	    ra = float(keywords['CRVAL1'][0])
	    dec = float(keywords['CRVAL2'][0])
	    fi = float(keywords['CROTA1'][0])
	    minmag = 0
	    maxmag = 15
	    stars = phot.get_stars(ra, dec, width = width, height = height, \
				   fi = fi,\
				   w = int(keywords['NAXIS1'][0]),\
				   h = int(keywords['NAXIS2'][0]),\
				   minmag = minmag, maxmag = maxmag)
	    for a in stars:
		a.ra = a.ra * 0.01745329252
		a.dec = a.dec * 0.01745329252
	        #fitsblink.Append_marker(0, a.x, a.y, 5.0, a.id, "red", "circle")
		#print a
            #  Read the original star list
	    datname = split(name, '.')[0] + '.dat'
	    for line in fileinput.input(datname):
		#  Skip the first line
		if fileinput.lineno() == 1: continue
		(intensity, ra, dec) = map(float, split(line)[3:6])
		#  Find if any star from photometric catalog matches
		ra = ra * 0.01745329252
		dec = dec * 0.01745329252
		for a in stars:
		    if a.i != 99.0 and a.i != -99.0 and angular_distance(ra, dec, a.ra, a.dec) < 0.0000096:
		    	print '%20s %8d %7.3f %7.3f %8.4f %8.4f %8.4f %8.4f' %\
			( a.id, intensity, a.i, a.si, a.ra / 0.01745329252, \
	ra /  0.01745329252, a.dec / 0.01745329252, dec / 0.01745329252)
	    #fitsblink.Set_image(0)
	    #fitsblink.Draw_image(0)
            #while 1:
	    #	fitsblink.loop()