Mathematics for Conversion between ICRS and AltAz

Question

I am writing some software, for my own learning experience, to convert astronomical coordinates to different frames of reference. An example of this would be conversion between ICRS which a lot of NASA data is supplied in automatically now, to Alt Azimuth, which is required for observation, based on the current time and geographic location.

From my prior studies of coordinate frames while I attended university, I know that all coordinate frames can be represented as a series of Matrix transformations from a "base" frame of reference. For instance, rotations around the X,Y,Z axis and a translation. For simplicity, these transformations are often all combined into a single affine transformation matrix like seen below where the blue section is the rotation matrix, the red is the translation, and the green is the scale.

This leads me to believe that I should be able to represent conversions from ICRS to Alt Azimuth in much the same way, though the transformation matrix will change with time/location. Based on this, the transformation matrix should look something like below where each rotation needed to create the overall transformation is a function of time, latitude, and longitude.

M = Rx(t, lat, lon) * Ry(t, lat, lon) * Rz(t, lat, lon)

However, in my googling of this topic all I can find is how to use AstroPy's built in functions to convert between different frames of reference. This is not all that useful to me as I actually want to learn about these transformations, not just be able to perform them. While I could study the AstroPy source code, it's such a massive library with a lot of interconnecting parts that understanding this one topic would take a huge time investment.

Does anyone have any experience with the mathematics involved in performing such a conversion, references to articles, or recommended textbooks that discuss the topic of conversion between different coordinate systems specifically those using a matrix format for the context of astronomy?

Answer 1

Sean E. Urban and P. Kenneth Seidelmann, Explanatory Supplemement to the Astronomical Almanac 3rd ed. (Mill Valley, CA: University Science Books, 2013).

This book provides a detailed treatment, and describes the methods used to create the annual Astronomical Almanac issued by HM Nautical Almanac Office and the US Naval Observatory. Be sure to download the errata.

Answer 2

The number of details involved depends on how accurate you want to be. There are more than just rotations involved though, you also need to account for parallax, aberration, gravitational deflection, light time distance, non-uniform rotation of the Earth, etc. There are some which are pure rotations, like the Earth Rotation angle, Precession, Nutation, and polar motion. The Explanatory Supplement mentioned by Gerard goes through all of the gory details if you really need that level of accuracy.

But most people won't need that level of accuracy. Below is a link to code I wrote which gets to sub-arc second accuracy for recent dates (compared to the results from JPL Horizons). Accounting for the non-uniform rotation of the Earth more accurately, and the ephemeris used (VSOP87) are the current limits to its accuracy.

https://www.celestialprogramming.com/vsop87-multilang/example_alt_az.html

Going through every step would be very tedious, but I'll provide overview of the "reduce" function which applies the rotations and other corrections you are asking about.

It first converts a time in UTC to Terrestrial Time (TT), which accounts for leap seconds, and somewhat for the Earth's non-uniform rotation. Next it gets the position of the Earth relative to the Sun, then gets the position of the requested planet relative to the sun. The Earth's position is subtracted from the planet's position to provide a geocentric position. Then the observer's position relative to the center of the Earth is computed. Rotation matrixes are created for the precession and nutation corrections. The inverse of these matrixes are applied to the observer's location because the apparent effects of precession and nutation are actually caused by the observer moving in the opposite direction. The result, is the XYZ position of the observer in the J2000/ICRF frame. Subtracting this position from the planet's yeilds the astrometric J2000 RA/Dec of the planet. Applying the precession and nutation matrix to the astrometric coordinates yields an "of date" position, and finally applying the affect of aberration yields an apparent position in RA/Dec. The last step is to convert those RA/Dec coordinates into alt/az coordinates, which is a rotation about the Earth's axis based on the current sidereal time.

[I am currently working on an article which walks through all of the steps outlined in the Explanatory Supplement in much less terse terms, I will try to remember to update this when it's done.]