import sys
import math
import pyfits
import time

def getNativeCoordinate( (alpha, delta), (alpha0, delta0) ):
  da = alpha - alpha0
  sina = math.sin( math.radians(da) )
  sind = math.sin( math.radians(delta) )
  sin0 = math.sin( math.radians(delta0) )
  cosa = math.cos( math.radians(da) )
  cosd = math.cos( math.radians(delta) )
  cos0 = math.cos( math.radians(delta0) )

  phi = math.atan2( cosd*sina, cos0*cosd*cosa+sin0*sind )
  theta = math.asin( -sin0*cosd*cosa + cos0*sind )

  return ( math.degrees(phi), math.degrees(theta) )

def getAngle(theta,phi):
  cost=math.cos(math.radians(theta))
  sint=math.sin(math.radians(theta))
  cosp=math.cos(math.radians(phi))
  sinp=math.sin(math.radians(phi))
  ang=math.atan2(math.sqrt((cost*sinp)**2+sint**2),cost*cosp)
  return ang

def MET2date (met, formatString = '%Y/%j:%H:%M:%S'):
  MET_OFFSET=978307200.
  return time.strftime(formatString, time.gmtime(met + MET_OFFSET))


if __name__=='__main__':

    ft2file1=sys.argv[1]
    ra=float(sys.argv[2])
    dec=float(sys.argv[3])
    ang=float(sys.argv[4])
    rockang=50.
    
    # retrive data from the FT2
    hdulist=pyfits.open(ft2file1)
    
    SC_data1=hdulist['SC_DATA'].data
    
    ra_zenith1=SC_data1.field('RA_ZENITH')
    dec_zenith1=SC_data1.field('DEC_ZENITH')
    time1=SC_data1.field('START')
    
    hdulist.close()

    MET1=[]
    RAZ1=[]
    DECZ1=[]
    AngZenith=[]
    
    for i in range(len(time1)):
      (phi,theta)=getNativeCoordinate((ra,dec),(ra_zenith1[i],dec_zenith1[i]))
      AngZenith.append(math.degrees(getAngle(phi,theta)))
      if i==1: # this is intentional, the i=0 entry is bogus
        if AngZenith[i]>ang: 
          print "//", MET2date(time1[i]), "Survey     Begin "
          print "//        offset = ",rockang, " " 
        else:
          print "//", MET2date(time1[i]), "Obs    Begin "
          print "//        RA = ",ra, " deg  "
          print "//        dec = ",dec, " deg  "
      if i > 1:
        if AngZenith[i]>ang and AngZenith[i-1]<ang:
          print "//", MET2date(time1[i]), "Obs    End "
          print " "
          print "// ---------------------------------------"
          print " "
          print "//", MET2date(time1[i]), "Survey    Begin "
          print "//            offset = ",rockang, " deg " 
        if AngZenith[i]<ang and AngZenith[i-1]>ang:
          print "//", MET2date(time1[i]), "Survey    End "
          print " "
          print "// ---------------------------------------"
          print " "
          print "//", MET2date(time1[i]), "Obs    Begin "
          print "//        RA = ",ra, " " 
          print "//        dec = ",dec, " "
          rockang*=-1.

#      MET1.append(time1[i])
#      RAZ1.append(ra_scz1[i])
#      RAZ2.append(ra_scz2[i])