function gps(latitude, longitude, timezone)
{
  this.latitude = latitude
  this.longitude = longitude
  this.timezone = timezone
}

function astro(where, when)
{
  this.where = where
  this.when = when
  this.sunrise = {time:0, azimuth:0}
  this.sunset = {time:0, azimuth:0}
  this.moonrise = {time:0, azimuth:0}
  this.moonset = {time:0, azimuth:0}
}
  
function SunUp(data)
{
  // function to calculate times of sunrise and sunset for given
  // longitude, latitude, time zone, year, month, day

  // This program by Roger W. Sinnott calculates the times of sunrise
  // and sunset on any date, accurate to the minute within several
  // centuries of the present.  It correctly describes what happens in the 
  // arctic and antarctic regions, where the Sun may not rise or set on
  // a given date.  Enter north latitudes positive, west longitudes
  // negative.  For the time zone, enter the number of hours west of
  // Greenwich (e.g., 5 for EST, 4 for EDT).  The calculation is
  // discussed in Sky & Telescope for August 1994, page 84.

  // define constants
  p1 = Math.PI; p2 = p1+p1
  dr = p1/180; k1 = 15*dr*1.0027379

  with (data) {  
    b5 = where.latitude; l5 = where.longitude/360; z0 = where.timezone/24
    y = when.getFullYear(); month = when.getMonth()+1; day = when.getDate()
  }

  // convert date to Julian days
  with (Math) {
    d1 = floor(day); f = day-d1-0.5
    j = -floor(7 * (floor((month+9)/12) + y) / 4)
    if (y >= 1583) 
    {
      mm = month-9; am = Math.abs(mm); s = (mm < 0) ? -1 : ((mm == 0) ? 0 : 1)
      j3 = floor(y + s*floor(am/7))
      j3 = -floor((floor(j3/100)+1)*3/4)
      j += j3 + 2
    }
    j = j + floor(275*month/9) + d1
    if (f < 0) {
      f++; j--
    }
  }
  t = j + (367*y - 730518) + f
  tt = t / 36525 + 1    // tt = centuries from 1900

  // LST at 0h zone time
  t0 = t / 36525
  s = 24110.5 + 8640184.813*t0
  s += 86636.6*z0 + 86400*l5
  s = s/86400; s = s - Math.floor(s)
  t0 = s*360*dr
  t += z0

  // get sun's position

  // fundamental arguments (Van Flandern & Pulkkinen, 1979)
  A = new Array(2); D = new Array(2)

  for ( i = 0; i < 2; i++)
  {
    with (Math) {
      l = .779072+.00273790931*t
      l = (l-floor(l))*p2
      g = .993126+.0027377785*t
      g = (g-floor(g))*p2
      v = .39785*sin(l)
      v = v-.01000*sin(l-g)
      v = v+.00333*sin(l+g)
      v = v-.00021*tt*sin(l)
      u = 1-.03349*cos(g)
      u = u-.00014*cos(2*l)
      u = u+.00008*cos(l)
      w = -.00010-.04129*sin(2*l)
      w = w+.03211*sin(g)
      w = w+.00104*sin(2*l-g)
      w = w-.00035*sin(2*l+g)
      w = w-.00008*tt*sin(g)
    
      // Compute sun's RA and Dec
      s = w / sqrt(u - v*v)
      a5 = l + atan(s / sqrt(1 - s*s))
      s = v / sqrt(u)
      d5 = atan(s / sqrt(1 - s*s))
      r5 = 1.00021 * sqrt(u)
      A[i] = a5; D[i] = d5
      t = t+1
    }
  }
  if (A[1] < A[0]) A[1] += p2
  z1 = dr * 90.833      // Zenith dist.
  s = Math.sin(b5 * dr); c = Math.cos(b5 * dr)
  z = Math.cos(z1); m8 = w8 = 0  
  a0 = a2 = A[0]; d0 = d2 = D[0]
  da = (A[1] - a0)/24; dd = (D[1] - d0)/24
  for (c0 = 0; c0 < 24; c0++)
  {
    a2 += da; d2 += dd
    
    // test an hour for an event
    l0 = t0 + c0*k1; l2 = l0 + k1
    h0 = l0 - a0; h2 = l2 - a2
    h1 = (h2 + h0) / 2  // hour angle,
    d1 = (d2 + d0) / 2  // declination at half hour
    if (c0 == 0) v0 = s*Math.sin(d0) + c*Math.cos(d0)*Math.cos(h0) - z
    v2 = s*Math.sin(d2) + c*Math.cos(d2)*Math.cos(h2) - z
    if ((v0<0 && v2>0) || (v0>0 && v2<0) || (v0==0 && v2!=0))
    {
      v1 = s * Math.sin(d1) + c * Math.cos(d1) * Math.cos(h1) - z
      a = 2*v2 - 4*v1 + 2*v0; b = 4*v1 - 3*v0 - v2
      d = b*b - 4*a*v0
      if (d >= 0)
      {
        d = Math.sqrt(d)
        e = (-b+d)/(2*a)
        if (e > 1 || e < 0) e = (-b-d)/(2*a)
        t3 = c0 + e + 1/120     // round off
        h3 = Math.floor(t3); h3 += Math.floor((t3 - h3)*60)/100
        h7 = h0 + e * (h2 - h0)
        n7 = - Math.cos(d1) * Math.sin(h7)
        d7 = c*Math.sin(d1) - s*Math.cos(d1)*Math.cos(h7)
        az = Math.atan(n7/d7)/dr
        if (d7 < 0) az += 180
        if (az < 0) az += 360
        if (az > 360) az -= 360
        az = Math.floor(az+0.5)
        if (v0 < 0 && v2 > 0) {
          m8 = 1; data.sunrise.time = h3; data.sunrise.azimuth = az
        }
        else if (v0 > 0 && v2 < 0) {
          w8 = 1; data.sunset.time = h3; data.sunset.azimuth = az
        }
      }
    }
    a0 = a2; d0 = d2; v0 = v2;
  }

  // if no sunrise or sunset found then special case
  if (m8 != 0 || w8 != 0) {
    if (m8 == 0)
      data.sunrise.time = "No sunrise"
    else if (w8 == 0)
      data.sunset.time = "No sunset"
  }
  else {
    if (v2 < 0)
      data.sunrise.time = data.sunset.time = "Sun down"
    else
      data.sunrise.time = data.sunset.time = "Sun up"
  }
}  

function MoonUp(data)
{
  // this program computes the times of moonrise and moonset
  // anywhere in the world. Form Sky and Telescope, July 1988
  
  // define constants
  p1 = Math.PI; p2 = p1+p1
  dr = p1/180; k1 = 15*dr*1.0027379

  with (data) {  
    b5 = where.latitude; l5 = where.longitude/360; z0 = where.timezone/24
    y = when.getFullYear(); month = when.getMonth()+1; day = when.getDate()
  }

  // convert date to Julian days
  with (Math) {
    d1 = floor(day); f = day-d1-0.5
    j = -floor(7 * (floor((month+9)/12) + y) / 4)
    if (y >= 1583) 
    {
      mm = month-9; am = Math.abs(mm); s = (mm < 0) ? -1 : ((mm == 0) ? 0 : 1)
      j3 = floor(y + s*floor(am/7))
      j3 = -floor((floor(j3/100)+1)*3/4)
      j += j3 + 2
    }
    j = j + floor(275*month/9) + d1
    if (f < 0) {
      f++; j--
    }
  }
  t = j + (367*y - 730518) + f
  tt = t / 36525 + 1    // tt = centuries from 1900

  // LST at 0h zone time
  t0 = t / 36525
  s = 24110.5 + 8640184.813*t0
  s += 86636.6*z0 + 86400*l5
  s = s/86400; s = s - Math.floor(s)
  t0 = s*360*dr

  t += z0

  // position loop
  mm = new Array(3)
  for (i = 0; i < 3; i++) mm[i] = new Array(3)
  for (i = 0; i < 3; i++) {
    // fundamental arguments
    with (Math) {
      l=0.606434+0.03660110129*t
      m=0.374897+0.03629164709*t
      f=0.259091+0.03674819520*t
      d=0.827362+0.03386319198*t
      n=0.347343-0.00014709391*t
      g=0.993126+0.00273777850*t
      l=(l-floor(l))*p2; m=(m-floor(m))*p2
      f=(f-floor(f))*p2; d=(d-floor(d))*p2
      n=(n-floor(n))*p2; g=(g-floor(g))*p2
      d2 = d+d; f2 = f+f; m2 = m+m; n2 = n+n
      v=0.39558*sin(f+n)
      v+=0.08200*sin(f)
      v+=0.03257*sin(m-f-n)
      v+=0.01092*sin(m+f+n)
      v+=0.00666*sin(m-f)
      v-=0.00644*sin(m+f-d2+n)
      v-=0.00331*sin(f-d2+n)
      v-=0.00304*sin(f-d2)
      v-=0.00240*sin(m-f-d2-n)
      v+=0.00226*sin(m+f)
      v-=0.00108*sin(m+f-d2)
      v-=0.00079*sin(f-n)
      v+=0.00078*sin(f+d2+n)
      u=1-0.10828*cos(m)
      u-=0.01880*cos(m-d2)
      u-=0.01479*cos(d2)
      u+=0.00181*cos(m2-d2)
      u-=0.00147*cos(m2)
      u-=0.00105*cos(d2-g)
      u-=0.00075*cos(m-d2+g)
      w=0.10478*sin(m)
      w-=0.04105*sin(f2+n2)
      w-=0.02130*sin(m-d2)
      w-=0.01779*sin(f2+n)
      w+=0.01774*sin(n)
      w+=0.00987*sin(d2)
      w-=0.00338*sin(m-f2-n2)
      w-=0.00309*sin(g)
      w-=0.00190*sin(f2)
      w-=0.00144*sin(m+n)
      w-=0.00144*sin(m-f2-n)
      w-=0.00113*sin(m+f2+n2)
      w-=0.00094*sin(m-d2+g)
      w-=0.00092*sin(m2-d2)
      
      // compute ra, dec, dist
      s=w/sqrt(u-v*v)
      a5=l+atan(s/sqrt(1-s*s))
      s=v/sqrt(u); d5=atan(s/sqrt(1-s*s))
      r5=60.40974*sqrt(u)
    }
    mm[i][0]=a5; mm[i][1]=d5; mm[i][2]=r5
    t += 0.5
  }
  if (mm[1][0] <= mm[0][0]) mm[1][0] += p2
  if (mm[2][0] <= mm[1][0]) mm[2][0] += p2
  z1=dr*(90.567-41.685/mm[1][2])
  s = Math.sin(b5*dr); c = Math.cos(b5*dr)
  z = Math.cos(z1); m8 = w8 = 0
  a0 = mm[0][0]; d0 = mm[0][1]
  for (c0 = 0; c0 < 24; c0++)
  {
    p = (c0+1)/24
    f0 = mm[0][0]; f1 = mm[1][0]; f2 = mm[2][0]
    a = f1-f0; b = f2-f1-a; a2 = f0+p*(a+a+b*(p+p-1))   // 3-point interp
    f0 = mm[0][1]; f1 = mm[1][1]; f2 = mm[2][1]
    a = f1-f0; b = f2-f1-a; d2 = f0+p*(a+a+b*(p+p-1))   // 3-point interp
    // test an hour for an event
    l0=t0+c0*k1; l2=l0+k1
    if (a2 < a0) a2 += p2
    h0 = l0-a0; h2 = l2-a2
    h1 = (h2+h0)/2      // hour angle
    d1 = (d2+d0)/2      // dec
    if (c0 == 0) v0 = s*Math.sin(d0)+c*Math.cos(d0)*Math.cos(h0)-z
    v2 = s*Math.sin(d2)+c*Math.cos(d2)*Math.cos(h2)-z
    if ((v0<0 && v2>0) || (v0>0 && v2<0) || (v0==0 && v2!=0))
    {
      v1=s*Math.sin(d1)+c*Math.cos(d1)*Math.cos(h1)-z
      v1=4*v1
      a=v2+v2-v1+v0+v0; b=v1-3*v0-v2
      d = b*b-4*a*v0
      if (d >= 0)
      {
        d = Math.sqrt(d)
        e = (-b+d)/(2*a)
        if (e > 1 || e < 0) e = (-b-d)/(2*a)
        t3 = c0 + e + 1/120     // round off
        h3 = Math.floor(t3); h3 += Math.floor((t3 - h3)*60)/100
        h7 = h0 + e * (h2 - h0)
        n7 = - Math.cos(d1) * Math.sin(h7)
        d7 = c*Math.sin(d1) - s*Math.cos(d1)*Math.cos(h7)
        az = Math.atan(n7/d7)/dr
        if (d7 < 0) az += 180
        if (az < 0) az += 360
        if (az > 360) az -= 360
        az = Math.floor(az+0.5)
        if (v0 < 0 && v2 > 0) {
          m8 = 1; data.moonrise.time=h3; data.moonrise.azimuth=az
        } else if (v0 > 0 && v2 < 0) {
          w8 = 1; data.moonset.time=h3; data.moonset.azimuth=az
        }
      }
    }
    a0=a2; d0=d2; v0=v2
  }
  
  // if no moonrise or moonset found then special case
  if (m8 != 0 || w8 != 0) {
    if (m8 == 0)
      data.moonrise.time = "No moonrise"
    else if (w8 == 0)
      data.moonset.time = "No moonset"
  }
  else {
    if (v2 < 0)
      data.moonrise.time = data.moonset.time = "Moon down"
    else
      data.moonrise.time = data.moonset.time = "Moon up"
  }
}


