code: public static double alongTrackDistanceRad(double lat1, double lon1, double lat2, double lon2, double latQ, double lonQ, double dist1Q, double ctd) {
      // FIXME: optimize the sign computation!
      int sign = Math.abs(bearingRad(lat1, lon1, lat2, lon2) - bearingRad(lat1, lon1, latQ, lonQ)) < HALFPI ? +1 : -1;
      return sign * acos(cos(dist1Q) / cos(ctd));
      // TODO: for short distances, use this instead?
      // asin(sqrt( (sin(dist_1Q))^2 - (sin(XTD))^2 )/cos(XTD))
    }
docs: The along track distance, is the distance from S to Q along the track S to
  E.
  
  ATD=acos(cos(dist_SQ)/cos(XTD))
  
  TODO: optimize: can we do a faster sign computation?
  
  @param lat1 Latitude of starting point in radians.
  @param lon1 Longitude of starting point in radians.
  @param lat2 Latitude of destination point in radians.
  @param lon2 Longitude of destination point in radians.
  @param latQ Latitude of query point in radians.
  @param lonQ Longitude of query point in radians.
  @param dist1Q Distance S to Q in radians.
  @param ctd Cross-track-distance in radians.
  @return Along-track distance in radians. May be negative - this gives the
  side.