Degenerate Conic

Algorithms • Modern Fortran Programming • Orbital Mechanics

Feb 01, 2018

IAU SOFA Version 14

iau_wb

The IAU Standards of Fundamental Astronomy (SOFA) library implements standard models used in fundamental astronomy. Version 14 ( 2018-01-30) has just been released. According to the release notes, this update includes the following:

  • Change in the copyright status of the iau_DAT routine. This is to provide for a user-supplied mechanism for updating the number of leap seconds. SOFA doesn't provide any such mechanism, but they are now allowing you to replace or modify this routine to provide your own without having to rename the routine (all other SOFA routines cannot be modified without renaming them).
  • Implementation of two new categories of routines:
    • Three new routines for the horizon/equatorial plane coordinates.
      • iau_AE2HD — (azimuth, altitude) to (hour angle, declination)
      • iau_HD2AE — (hour angle, declination) to (azimuth, altitude)
      • iau_HD2PA — parallactic angle
    • Six new routines dealing with gnomonic (tangent plane) projections.
      • iau_TPORS — solve for tangent point, spherical
      • iau_TPORV — solve for tangent point, vector
      • iau_TPSTS — project tangent plane to celestial, spherical
      • iau_TPSTV — project tangent plane to celestial, vector
      • iau_TPXES — project celestial to tangent plane, spherical
      • iau_TPXEV — project celestial to tangent plane, vector
  • The Astrometry Tools Cookbook, the test program and other supporting files have also been updated.
  • Other minor documentation/typographical corrections to various files.

IAU SOFA has a Fortran 77 and a C version. While this is a great resource, I do wonder why the IAU insists on using a source format for the Fortran version that was rendered obsolete 30 years ago? There's no reason for this. This is new code that is regularly kept up to date, so why not modernize it for those of us living in the 21st century?

Consider this SOFA routine:

      SUBROUTINE iau_TR ( R, RT )

      IMPLICIT NONE

      DOUBLE PRECISION R(3,3), RT(3,3)

      DOUBLE PRECISION WM(3,3)
      INTEGER I, J

      DO 2 I=1,3
        DO 1 J=1,3
          WM(I,J) = R(J,I)
1       CONTINUE
2     CONTINUE
      CALL iau_CR ( WM, RT )

      END SUBROUTINE iau_TR

All this is doing is transposing a 3x3 matrix. Well, since the late 1990s, to do that in Fortran all you need is:

rt = transpose(r)

The original routine is just awful: all caps, line numbers, CONTINUE statements, and 100% unnecessary. Unfortunately, it's code like this that gives Fortran a bad name (reinforcing the perception that Fortran is an arcane and obsolete programming language). Anyone who is compiling this code is probably using a compiler that is at least Fortran 95 compatible. Who on earth is using straight up Fortran 77 compilers at this point? So why restrict such a useful library to a programming style that is so totally out of date? I think it is even probably affecting the efficiency of the code, since the example above required using a temporary array and a function call to do what can now be done with a (presumably much more efficient) intrinsic routine. It wouldn't even be that hard to convert this code to modern Fortran style. It's just low-level math routines, it doesn't have to be object-oriented or anything fancy like that. A while back, I wrote a little Python script to merge all the SOFA files into a single Fortran module (another innovation from the 1990s) so they can take advantage of the automatic interface checking that modules provide. There are also any number of tools out that could be used to automate the fixed to free-form conversion.

See also