#!/bin/sh -- # perl
# \
eval 'exec perl -S $0 ${1+"$@"}'
    if 0;
use POSIX;

# uy4.phases.pl
#
# Lance A. M. Benner
# Created:   1999 April 2  (for 1992 SK)
# Modified:  2003 September 24
#            2004 February 23 (3908 Nyx in 2004; significant improvements)
#            2005 January 12   2004 MN4
#            2007 August 29    4660 Nereus: Arecibo

# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 

# Compute the rotation phases for NEA Nereus.  Normally
# this script uses the UTC rise and set times and two minute
# intervals to compute phases as radial line segments suitable for plotting with
# tkgraph (and probably grace too, although that hasn't been tested), but 
# in this case I used the actual receive start-stop times on specific dates
# when we observed the asteroid.

# The user has to provide the UTC rise and set times in doy.dddd, the number of 
# days to check, and the orbital and rotation periods.

# The script outputs results in radial line segment format suitable for plotting
# in tkgraph with polar plots.  When reading the data into tkgraph, be sure to 
# select "r theta" as the input format.

# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 

# Zero phase epoch: 2002 January 09 00:00:00 UTC (DOY 009)
$zeroepoch = 9.0;                 # Zero phase epoch in days-of-year (2002 January 9)
$period = 15.1/24.0;              # rotation period in days

# UTC Start and stop epochs.  Format: "DOY hh mm ss" in UTC
$start[1] = "009 18 46 04";          # 2002 Jan 09 Arecibo
$stop[1]  = "009 20 35 31";

$start[2] = "010 18 48 46";          # 2002 Jan 10 Arecibo
$stop[2]  = "010 20 21 26";

$start[3] = "011 18 27 32";          # 2002 Jan 11 Arecibo
$stop[3]  = "011 20 04 19";

$start[4] = "012 18 21 29";          # 2002 Jan 12 Arecibo
$stop[4]  = "012 19 44 50";

$start[5] = "013 18 11 43";          # 2002 Jan 13 Arecibo
$stop[5]  = "013 18 41 26";

# Output files
$outfile[1] = "nereus.jan09.phases.a.dat";
$outfile[2] = "nereus.jan10.phases.a.dat";
$outfile[3] = "nereus.jan11.phases.a.dat";
$outfile[4] = "nereus.jan12.phases.a.dat";
$outfile[5] = "nereus.jan13.phases.a.dat";

$daymin = 1;
$daymax = 5;

# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 

$dt = (60.0)/(86400.0);             # Time step in units of days (= 1 minute)
print "  \n";
for ($day = $daymin; $day <= $daymax; $day++)
{
	print "$outfile[$day]\n";

	@begin = split(/ /,$start[$day]);
	@end   = split(/ /,$stop[$day]);

	print "start = DOY $begin[0]  $begin[1]:$begin[2]:$begin[3] UTC\n";
	print "stop  = DOY $end[0]  $end[1]:$end[2]:$end[3] UTC\n";

	$rcsta = $begin[0] + (($begin[1] + $begin[2]/60.0 + $begin[3]/3600.0)/24.0);
	$rcend = $end[0]   + (($end[1]   + $end[2]/60.0   + $end[3]/3600.0)/24.0);

	open(FILE, ">$outfile[$day]") || die "Can't open $outfile[$day]: $!\n";	
#	$time = $start[$day];                          # Set to the beginning of the track
#	while($time <= $stop[$day])                    # Loop over each track in 1 minute intervals

	$time = $rcsta;
	while($time <= $rcend )                        # Loop over each track in 1 minute intervals
	{
	    $rotations = ($time - $zeroepoch)/$period; # Number of rotations (decimal)
	    $temp1 = int $rotations;                   # Number of complete rotations (integer)
	    $phase = ($rotations - $temp1)*(360);      # Rotation phase in degrees
#	    print "$time  $rotations  $temp1  $phase\n";
	    print FILE " 8  $phase\n";                 # 8 specifies the lower radius limit (segments)
	    print FILE "12  $phase\n";                 # 12     "        upper          "
	    $time = $time + $dt;
	}

	print "\n";
#	print "$time\n";

	close(FILE);
	print "\n";
}

print "\n";

# --------------------------------------------------------------------