Replace ast_julday_to_date with a BSD-3-Clause implementation

The implementation is taken from ERFA, which is a relicensed (with
permission) variant of IAU SOFA and adopted to SPP. The license is
also added to local/LICENSES.

local/LICENSES/ERFA

@@ -0,0 +1,53 @@
+Copyright (C) 2013-2014, NumFOCUS Foundation.
+All rights reserved.
+
+This library is derived, with permission, from the International
+Astronomical Union's "Standards of Fundamental Astronomy" library,
+available from http://www.iausofa.org.
+
+The ERFA version is intended to retain identical
+functionality to the SOFA library, but made distinct through
+different function and file names, as set out in the SOFA license
+conditions. The SOFA original has a role as a reference standard
+for the IAU and IERS, and consequently redistribution is permitted only
+in its unaltered state. The ERFA version is not subject to this
+restriction and therefore can be included in distributions which do not
+support the concept of "read only" software.
+
+Although the intent is to replicate the SOFA API (other than replacement of
+prefix names) and results (with the exception of bugs; any that are
+discovered will be fixed), SOFA is not responsible for any errors found
+in this version of the library.
+
+If you wish to acknowledge the SOFA heritage, please acknowledge that
+you are using a library derived from SOFA, rather than SOFA itself.
+
+
+TERMS AND CONDITIONS
+
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions are met:
+
+1 Redistributions of source code must retain the above copyright
+   notice, this list of conditions and the following disclaimer.
+
+2 Redistributions in binary form must reproduce the above copyright
+   notice, this list of conditions and the following disclaimer in the
+   documentation and/or other materials provided with the distribution.
+
+3 Neither the name of the Standards Of Fundamental Astronomy Board, the
+   International Astronomical Union nor the names of its contributors
+   may be used to endorse or promote products derived from this software
+   without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
+IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
+TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
+HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
+TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
+PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
+LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
+NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
+SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

local/LICENSES/README

@@ -5,5 +5,5 @@ to third-party sources in the IRAF system:
     GPL					For the math$slalib and ecl$readline
     OpenSolaris_License-CDDL.pdf	For the XYACC parser (boot$xyacc)
     UCAR				For NCAR graphics
-
+    ERFA				For ast_julday_to_date in noao$astutil/asttools/asttimes.x
 

noao/astutil/asttools/asttimes.x

@@ -147,41 +147,96 @@ end
 
 
 # AST_JULDAY_TO_DATE -- Convert Julian date to calendar date.
-# This is taken from Numerical Receipes by Press, Flannery, Teukolsy, and
-# Vetterling.
+# This procedure is taken from `eraJd2cal` of ERFA (Essential Routines
+# for Fundamental Astronomy) and converted to SPP.  ERFA is a C
+# library containing key algorithms for astronomy, and is based on the
+# `SOFA library <http://www.iausofa.org/>`_ published by the
+# International Astronomical Union (IAU).
+#
+# ERFA is intended to replicate the functionality of SOFA (aside from
+# possible bugfixes in ERFA that have not yet been included in SOFA),
+# but is licensed under a three-clause BSD license to enable its
+# compatibility with a wide range of open source licenses. Permission
+# for this release has been obtained from the SOFA board
+#
+# ERFA can be found at https://github.com/liberfa/erfa/
+#
+#  Notes:
+#
+#  1) The earliest valid date is -68569.5 (-4900 March 1).  The
+#     largest value accepted is 1e9.
+#
+#  2) The Julian Date is apportioned in any convenient way between
+#     the arguments dj1 and dj2.  For example, JD=2450123.7 could
+#     be expressed in any of these ways, among others:
+#
+#            dj1             dj2
+#
+#         2450123.7           0.0       (JD method)
+#         2451545.0       -1421.3       (J2000 method)
+#         2400000.5       50123.2       (MJD method)
+#         2450123.5           0.2       (date & time method)
+#
+#     In this implementation, dj2 is set to 0.0.
+#
+#  3) In early eras the conversion is from the "proleptic Gregorian
+#     calendar";  no account is taken of the date(s) of adoption of
+#     the Gregorian calendar, nor is the AD/BC numbering convention
+#     observed.
+#
+#  Reference:
+#
+#     Explanatory Supplement to the Astronomical Almanac,
+#     P. Kenneth Seidelmann (ed), University Science Books (1992),
+#     Section 12.92 (p604).
+#
+#  Copyright (C) 2013-2017, NumFOCUS Foundation.
+#  Derived, with permission, from the SOFA library.
 
 procedure ast_julday_to_date (j, year, month, day, t)
 
-double	j			# Julian day
+double	j,dj2			# Julian day
+data dj2 /0./
 int	year			# Year
 int	month			# Month (1-12)
 int	day			# Day of month
 double	t			# Time for date (mean solar day)
 
-int	ja, jb, jc, jd, je
+long jd, l, n, i, k
+double d1, d2, f1, f2, f, d
+double mod(), nint()
 
 begin
-	ja = nint (j)
-	t = 24. * (j - ja + 0.5)
-
-	if (ja >= 2299161) {
-	    jb = int (((ja - 1867216) - 0.25) / 36524.25)
-	    ja = ja + 1 + jb - int (jb / 4)
+	# Copy the date, big then small, and re-align to midnight.
+	if (j >= dj2) {
+	  d1 = j
+	  d2 = dj2
+	} else {
+	  d1 = dj2
+	  d2 = j
 	}
-
-	jb = ja + 1524
-	jc = int (6680. + ((jb - 2439870) - 122.1) / JYEAR)
-	jd = 365 * jc + int (jc / 4)
-	je = int ((jb - jd) / 30.6001)
-	day = jb - jd - int (30.6001 * je)
-	month = je - 1
-	if (month > 12)
-	    month = month - 12
-	year = jc - 4715
-	if (month > 2)
-	    year = year - 1
-	if (year < 0)
-	    year = year - 1
+	d2 = d2 - 0.5
+
+	# Separate day and fraction.
+	f1 = mod(d1, 1.0)
+	f2 = mod(d2, 1.0)
+	f = mod(f1 + f2, 1.0)
+	if (f < 0.0) f = f + 1.0
+	d = nint(d1-f1) + nint(d2-f2) + nint(f1+f2-f)
+	jd = long(nint(d)) + 1
+
+	# Express day in Gregorian calendar.
+	l = jd + 68569
+	n = (4 * l) / 146097
+	l = l - (146097 * n + 3) / 4
+	i = (4000 * (l + 1)) / 1461001
+	l = l - ((1461 * i) / 4 - 31)
+	k = (80 * l) / 2447
+	day = int(l - (2447 * k) / 80)
+	l = k / 11
+	month = int(k + 2 - 12 * l)
+	year = int(100 * (n - 49) + i + l)
+	t = f * 24
 end