/** This is a test of altitude-azimuth to XYZ Cartesian coordinate conversions
    both forwared and backward. */

// Define and test conversions between cartesian (x, y, z) and
// alt-azimuth (radius, altitude, azimuth) coordinates.

// Altitude is the angle up from the horizon.

// Azimuth is the angle clockiwse from north with north being 0 degrees and
// east being 90 degrees (like standard compass bearing.s)
// The x axis is positive toward east.
// The y axis is positive toward north.
// z is up.

altAzToXYZ[r, alt, az] :=
{
   x = r * cos[alt] sin[az]
   y = r * cos[alt] cos[az]
   z = r * sin[alt]
   return [x,y,z]
}

xyzToAltAz[x, y, z] :=
{
   r   = sqrt[x^2 + y^2 + z^2]

   if r != 0 deg
      alt = arcsin[z / r]
   else
      alt = 0 deg
   
   az  = arctan[x, y]
   
   return [r, alt, az]
}

test[x, y, z, err=1e-14] :=
{
   [r, alt, az] = xyzToAltAz[x, y, z]
   [xp, yp, zp] = altAzToXYZ[r, alt, az]

   dx = x-xp
   dy = y-yp
   dz = z-zp

   if abs[dx] > err or abs[dy] > err or abs[dz] > err 
      println["$x $y $z : $dx $dy $dz : $r " + (alt->deg) + " " + (az->deg)]
}

p = new range[-2, 2, 1/2]
//p = [-1, 1]
multifor [sx, sy, sz] = [p, p, p]
{
   test[sx, sy, sz]
}