#
# Read in the data...
#
data "public/ip4sm"
define decs  read 1 1
define decdd read 1 2
define decmm read 1 3
define file read  1 4

#
# Calculate declination in radians...
#
define dec (($decdd + $decmm/60)*(3.142/180))
if ('$decs' == '-'){
  define dec (-1*$dec)
  define decdd (-1*$decdd) # needed to get label right...
}

#
# Define latitude...
#
define lat (-29.25)
define l ($lat*(3.14/180))

#
# Define vectors to store HA/PA
#
set dimen(pa) = 76
set dimen(ha) = 76

#
# Loop over HA, and calculate PA...
#
do ha=0,7.5,.1{
  define h ($ha*15*(3.14/180))
  define spang ((sin($h)*cos($l))/ \
  sqrt(1.0-(sin($l)*sin($dec)+cos($l)*cos($dec)*cos($h))**2))
  define pa (asin($spang)*(180/3.14))
  if ($dec > $l){
    define pa (180 - $pa)
  }
  if ((cos($h)  - tan($l)/tan($dec) < 0.) && ($dec < $l)) {
    define  pa (180 - $pa)
  }
  set ha[($ha*10)] = $ha
  set pa[($ha*10)] = $pa
}


#
# Draw graph...
#
device postencap "public/"$file
location 4000 32000 4000 30000
erase
limits 0 8 0 180
expand 0.99
box
xlabel HOUR ANGLE WEST OF MERIDIAN (HOURS)
ylabel OPTIMAL POSITION ANGLE (DEGREES)
connect ha pa

expand 0.8
define dec ($dec*(180/3.142))

relocate (18000 1000)
putlabel 5 (Note: Corresponding position angles for objects east \
of the meridian are negative.)

relocate (5000 29000)
#putlabel 6 Declination=$decdd\u0$decmm'
putlabel 6 Declination=$decdd  $decmm'
relocate (5000 28200)
putlabel 6 Geo. Latitude=$lat degrees