Constants used in equation_of_the_equinoxes_complimentary_terms()

[glangford]

In the calculation of equation_of_the_equinoxes_complimentary_terms()

python-skyfield/skyfield/nutationlib.py

Line 111 in eda5fea

(715923.2178 +

some of the values used seem to vary from other published sources. The value 715923.2178 appears in the code, but it appears as 1717915923.2178 in SOFA and in other places eg.
http://www.iausofa.org/2010_1201_C/sofa/fal03.c

There are similar variances in the calculation of the other fa values. Assuming this is intended, should there be code comments added to explain the source of the values and the discrepancies?

[brandon-rhodes]

Good question! Skyfield cites NOVAS, not SOFA, in its bibliography:

https://rhodesmill.org/skyfield/bibliography.html

And NOVAS in its C code uses the constant that you've cited in particular:

https://github.com/brandon-rhodes/python-novas/blob/a40c662493da86da51e9b09361b5415a4fbbbf9e/Cdist/novas.c#L4662

Let me know if there's a way the documentation can be cleared about the astrometry on which Skyfield is based.

In any case, those two constants look so very nearly similar that you might check whether the two libraries are using a different scale for the time?

[glangford]

The time parameter provided in the SOFA calculation of mean anomaly of the moon is 'Julian centuries since J2000.0' according to the code comments. That code cites IERS Technical Note 32, section 5.8.2 of that gives the constants used.

The Skyfield code looks similar, but this 1325 * t is unique and I don't understand that

python-skyfield/skyfield/nutationlib.py

Line 116 in eda5fea

+ (1325.0*t % 1.0) * tau)

Let me know if there's a way the documentation can be cleared about the astrometry on which Skyfield is based.

Is there a reference other than another source code base which describes the calculation and the constants used (similar to the way SOFA points to IERS)?

[glangford]

I came across this in my travels in case this is helpful for someone looking at this thread in future

A comparison of SOFA and NOVAS astrometric software libraries
https://ui.adsabs.harvard.edu/abs/2018SPIE10707E..1ZP/abstract

[brandon-rhodes]

To clarify: are you getting different values out of the routine and are trying to track down why their behavior is different? How big a difference, in that case?

[brandon-rhodes]

Oh, and: the NOVAS guide has an Appendix D that compares SOFA and NOVAS.

https://github.com/brandon-rhodes/python-novas/blob/a40c662493da86da51e9b09361b5415a4fbbbf9e/Cdist/NOVAS_C3.1_Guide.pdf

[brandon-rhodes]

(The Appendix D bibliography might also be worth a look)

[glangford]

No, I haven't actually run the two and compared behaviors - I was just trying to understand where the magic numbers came from, and saw that the two implementations used many of the same constants but not identical ones.

[glangford]

I think I see what might be going on...from Appendix D that you recommended

"Function ee_ct, which evaluates the “complementary terms” in the equation of the equinoxes, is based on IERS function EECT2000"

Here is the Fortran code for EECT2000, which computes the fa (fundamental arguments) per IERS 2000, and it looks the same as Skyfield

ftp://tai.bipm.org/iers/conv2003/chapter5/EECT2000.f

*  Fundamental Arguments (from IERS Conventions 2000)

*  Mean Anomaly of the Moon.
      FA(1) = iau_ANPM ( ( 485868.249036D0 +
     :                   ( 715923.2178D0 +
     :                   (     31.8792D0 +
     :                   (      0.051635D0 +
     :                   (     -0.00024470D0 )
     :                   * T ) * T ) * T ) * T ) * DAS2R
     :                   + MOD ( 1325D0*T, 1D0 ) * D2PI )

...whereas the code comments in SOFA suggests they are aligned with IERS 2003. So both SOFA and NOVAS code bases have similar heritage but something changed in the IERS recommendations between 2000 and 2003.

[brandon-rhodes]

Thanks for the update on your findings! Skyfield has so far been largely feature-driven — for example, at the moment two different historians are asking about discrepancies in results, so I undertaking work to see if Skyfield can be taught to use the Julian calendar for dates before the Gregorian calendar reform, to match textbooks and papers in history (none of which, it seems, ever uses the proleptic Gregorian calendar!).

All of which leaves Skyfield without a planned roadmap for upgrades to its fundamental astrometry. I have been aware that the IAU continues apace with ever more precise formulae, but have not had the time to undertake a several-week review of which routines now have updated formulae available from the authorities.

I suppose I could just wait until people see discrepancies and open issues. But that would be rather haphazard, and in some cases it might not be clear which number precisely is different between, say, Skyfield and HORIZONS.

Another issue is that Skyfield currently uses the old precession-based approach under the hood, while, I believe, modern practice now uses the CIO as the pivot when computing the rotation from ICRS to equator-and-equinox-of-date. So there are three reasons, I suppose, that I have not been pursuing new IAU updates: (a) limited time; (b) the issues folks open are almost always about practical features that could be added, not about underlying fundamental routines, so I'm pushed towards those activities when I do have time to code; and (c) not wanting to burn a couple of weeks figuring out the maneuvers that will be necessary to switch from a nutation / precession to a CIO-based approach.

[glangford]

That's understandable, and thanks as always for your quick responses! I was coming at this from a learning perspective, and trying to better understand the foundations. I have found Skyfield very easy to use for several applications and I ended up wanting to learn the underlying code, and hoping to contribute code and help where I can.

I don't think I am quite at the level to contribute meaningfully to upgrading the fundamental astrometry just yet, but hopefully soon. Thanks again for your guidance!