Jpl horizons algorithm delta t#4673
Conversation
…lder with ERFA sources and CMakeLists.txt- Included license files from ERFA project- Updated top-level and module CMakeLists.txt to build ERFA as static library
|
Great PR! Please pay attention to the following items before merging: Files matching
This is an automatically generated QA checklist based on modified files. |
|
Thanks, definitely interesting, although I have no time to fully testdrive this in these weeks. Hopefully just around turn-of-the-year. I think, when we integrate ERFA, we should set this up via CPM like many other tools we changed from earlier "exernal"izing. But I am no expert in this setup, @alex-w, @10110111 ? You wrote you used cmake to build it, indicating its own installation instructions around yet another tool, meson, is not required? And then we might eventually utilize more of ERFA as time comes. |
|
I analyzed the failure on Linux and macOS and I believe that it is caused by an invalid value of TTmTDB from test case JD 2437684.50000, with DeltaT expected at 34.051667 seconds. This test case happens to be the very first one in the time range 1962-01-20, which mean that we now need TTmTDB extracted from linux_p1550p2650.440t to compute DeltaT to +/- 10 microseconds of JPL Horizons results. If the request to get TTmTDB is not available when we reach this first test case, then GetDe440Coor() will return deOK false (see line 2017 of StelUtils.cpp) and TTm TDB will fall back to 0.0. On my Windows-based system, I debugged the code of testDeltaT.cpp (with breakpoints) and saw TTmTDB = -0.000467609 second at JD = 2437684.50000, with lead seconds = 1.8671988 seconds. DeltaT is computed as follow : deltaT = 32.184 - TTmTDB + leapSeconds. If TTmTDB = 0, we would get DeltaT = 34.05120, which is not within the +/- 10 microseconds criteria, but if TTmTDB = -0.000467609, we would get DeltaT = 34.05167, which is within the +/- 10 microseconds criteria. So bottom line, for some reason, Linux and macOS are not ready to extract TT-TDB for this very first test case in the range of 1962-01-20 and up, but works for all the others test cases within that rnage... So the question is why Linux and macOS fail to be ready accessing the linux_p1550p2650.440t file on this very first request while Windows seems to be ready ? Note: while numbers reported by this automated review match my analysis, it remains to be proved... Someone should investigate on a Linux and macOS computer if this is really the case. I don't have a Linux and macOS computer, but I am more than willing to help out if needed. |
|
@gzotti ... You wrote you used cmake to build it, indicating its own installation instructions around yet another tool, meson, is not required? I did not need to use meson tool. I personally do not know this tool. I simply copied all the *.c and *.h files of ERFA in a directory, and first used Visual Studio 2026 to build a static library out of them. I realized within this development environment that I had to set some macros manually (which probably meson was taking care...), and fixed the code in some files (those were not related to SOFA, but very specific to ERFA, and these changes are very minor and have no impact on astrometry calculation). Then I figured out how this building process could be automated by CMake within Stellarium. I more or less followed the same approach that I used to integrate the original SOFA library as I mentioned in issue #4670. The first step of the Pull Request is specific to the integration of ERFA into a static library, with all the changes in the CMakeLists.txt files needed. I understand that your team would prefer to integrate via a different approach (like CPM) if ERFA is an acceptable source of code. I have no issue with that. I do hope that you will eventually integrate this new functionality in Stellarium. If not, well, with this Pull Request you have all the details on how JPL Horizons compute Delta_T, so Enjoy! |
…izons ΔT algorithm
|
I run testDeltaT.cpp in Debug, Release and RelWithDebInfo build configurations on my Windows PC. The following three files are the respective output of test results that I get from these runs: TestDeltaT_Debug_2025-11-27_2136.txt What I am noting is that the failure report from Linux and macOS build by CI Workflow run do not show the QINFO that I normally get during the initTestCase() execution. Most notably QINFO : TestDeltaT::initTestCase() Use DE440 ephemeris file "C:\Users\simar\AppData\Roaming\Stellarium/ephem/linux_p1550p2650.440t" Are these automated test cases, aware that they need to use the t version of the DE440 files ? The QINFO also shows if the DE440 initialization was successful, but it does not appear on the CI Workflow run... I just update testDeltaT.cpp and pushed it to the same branch. I only added more comments to identify the calendar dates (using mixed calendars) to each tested Julian Date. The code has not changed so CI Workflow run should report the same failures. We'll see... |
ERFA doesn't use meson, it uses autotools (which is based on Bourne shell, and I'm not sure if our Windows build environment works well with it). This makes it less suitable for use with CPM. But at the same time, we don't want to always use the bundled version, since ERFA is packaged in at least some Linux distributions (e.g. Ubuntu since 22.04). Given all this, I suppose we should check for a local Another option would be to host a version of it converted to CMake, and then use it with CPM. |
|
I tried to push the ephem folder which contains all my DE files needed for properly testing the new algorithm . This folder includes very big ones (more than 2GB each). The git push fails pathetically... Is there any way someone at your end can ensure these files are accessible in stellarium/ephem before the automated CI test is executed ? Maybe this is the reason why it fails on some platform like Linux and macOS, but not on Windows..
|
And rightly so. Who will ever want to download useless extra 4GB of data when fetching the source repo? Normally all these files are downloadable during runtime (Configuration → Extras). I suppose the unit tests that require these files should be optional. Or, if it's possible, they could make use of some small substitutes (or even mock functions) that contain all the information necessary for the test. |
|
@ssimard2504 where we should get these files? |
|
The location and download of the optional DE4* files is described in the User Guide. No, these certainly won't go into the repo, and are always optional. Most users won't ever need them, and I expect most users also won't bother. (Of course I do like the prospects of getting so much closer...) I also know nothing about *.440t files. When the data files are missing, the test should still succeed with just a warning that it could not be performed on the original data files. Testing locally on a dev PC is then more helpful.
The building instructions say so, but I did not go further. The fewer extra tools the better. And if erfa is an available installable library on Debian-based systems, then yes, either use that or "USE_BUNDLED_ERFA" with sources either fetched ad-hoc (CPM?) and prepared for building with addition of e.g. one extra CMakefile, or we can still add a snapshot into external (and update this every 2 years or so?) and use that on Windows or whereever it is not pre-installed. |
|
The DE4* files are indeed described in the user guide (page 65) and a link is even provided at the bottom of the page to find them: https://ssd.jpl.nasa.gov/ftp/eph/planets/Linux/ . As per the 'current' user guide instructions, they have to be downloaded separately and are optional. They should be saved either in the user directory or installation directory, in a separate folder named /ephem. Or they could be saved in another location specified in the config.ini file as described in the user guide on page 65. However, this 3rd option does not work for some test cases (e.g. testEphemeris and testDeltaT that I recently modified to support testing JPL Horizons algorithm). In their cases, they have to be either in the user directory or installation directory. The new DE4* 't' files required for JPL Horizons are also available from the same JPL FTP website. The 't' version of DE430 and DE440 contains the data to extract TT - TDB from the ephemerides, which the JPL Horizons algorithm needs from date 1962-01-20 and Present. @10110111 ..And rightly so. Who will ever want to download useless extra 4GB of data when fetching the source repo? I am fairly new to all this business of git, CI, repo and whatnot. I had no idea that pushing large file was a no-no in this business. I was just attempting to figure out why testDeltaT was failing on Linux and macOS while it works fine on my Windows system. Since the test will not discover the files unless they are in the user or installation directory, I thought, why not upload them in the installation directory for the CI testing ? From all these discussions, I am reading between the lines that there are two philosophies at play here: The 'good enough' and the 'accurate'. I personally falls in the 'accurate' group and would much prefer that Stellarium aims for this goal in future release, though not if there is a significant penalty on performance. @gzotti ... When the data files are missing, the test should still succeed with just a warning that it could not be performed on the original data files. Testing locally on a dev PC is then more helpful... I guess I could modify the portion of testDeltaT testing JPL Horizons algorithm to act like testEphemeris if the DE4* files cannot be found. I'll do that later today. |
|
I guess 80-95% of users want today's sky and nothing else, and are happy to find and observe objects in their modest instruments. This requires "good enough" and worked many years ago. Here VSOP87 should still be sufficient, and is small enough to run on SBCs with limited storage/power. Here probably even devs don't want to run the tests with huge files, although building customized installs should always work. And then there are the scientific users who want more, be it higher accuracy or simulation in historical context. This is why we did a lot to improve accuracy, including DE4* use, but I stopped with aberration and don't care for gravitative light deflection or similar sub-arcsecond things. (This should not prevent others from doing it...) Some of those users may have the files already in their systems which is why you can configure other locations via config.ini to avoid duplicate gigabytes. Else, the user data directory's ephem subdirectory is the expected place. The program install directory happens to work by the way the directory searching works but should actually not be recommended (We should fix that sentence. Users should not mess around copying stuff there...). To test the Linux issues on a Windows PC, you can build everything under WSL2. If I install the .430t/.440t files, do I still need the non-t as well or do they include all data? In this case we can then change instructions... |
|
@gzotti ... If I install the .430t/.440t files, do I still need the non-t as well or do they include all data? In this case we can then change instructions... My understanding is that the .430t/.440t files are a superset of .430/.440. That is they should give the same position/velocity information of sun & planets but in addition, the 't' version includes Chebyshev polynomial This is what I read from the README file for de430t that you can find here: https://ssd.jpl.nasa.gov/ftp/eph/planets/Linux/de430t/ There is no README file for the de440t, but I would imagine this is also applicable to the de440t version. I kept both type of files for this development on my PC and so far I did not see any differences in RA/Dec when switching between them. However this is very limited check of mine.... On the other hand test case "testEphemeris" give a PASS with the 't' version, so it is pointing in the right direction. If you want to use these 't' version files, you must download them in the same location as the other de files. But you also have to modify two MACROs in EphemWrapper.hpp and rebuild Stellarium. #define DE430_FILENAME "linux_p1550p2650.430t" Bottom line, to answer your question, changing instructions in the user guide to use the 't' version should do it for now. IF, however, the Stellarium team go ahead to incorporate the new JPL Horizons Delta-T algorithm as I propose, the user guide will need a bit more changes and/or the Stellarium GUI may need to hook the activation of DE44x automatically when the JPL Horizons Delta-T algorithm is selected...otherwise all Delta-T after 1962-01-20 will be off by a huge.... 2 milliseconds ;-) |
|
hmm... what about JPL Delta-T when VSOP87 selected? |
|
IIRC we activate JPL DE* ΔT when DE* is used, not the other way. @ssimard2504 indicates when setting JPL ΔT results will be off by 2ms which I don't notice given likely other noise. What I am still unclear after all these years (or have forgotten) is the actual effect of the different ndots. Where do they impact? Lunar position? |
|
@alex-w, @gzotti when the JPL Horizons algo is selected and VSOP87 ephemerides are used (instead of DE44x) Stellarium will still compute Delta-T using the JPL algorithm, however, the correction for Moon Secular Acceleration for all years before 1955.5 will give a correction different from the JPL Horizons app. JPL 44x ephemerides n_dot (i.e. variable ephND) = -25.936 "/cy/cy, while VSOP87 ephemerides n_dot (ephND) =-23.8946 "/cy/cy. Since JPL Horizons algo is based on n_dot (nd) = -25.82 "/cy/cy (the same as for DE43x), the correction between VSOP87 and DE44x will be different - being computed as a function of (ephND - nd). Note that you could have selected DE440 (only this one...) from the Ephemeris settings GUI and still fallback to VSOP87 if you end-up selecting a date outside the range of DE440. DE440 covers 1549-12-21 to 2650-01-25 (assuming mixed calendar). So for a date like 1500-01-01 the Ephemeris used will be those of VSOP87, not DE440. For years after 1955.5 the Moon Secular acceleration correction should not be applied according to my email with Jon Giorgini (JPL). This is how JPL Horizons app work. I therefore added a line of code in getMoonSecularAcceleration() to return 0.0 for this correction for all years after 1955.5. I HAVEN'T VERIFIED IF THIS IS ALSO APPLICABLE TO THE OTHER DELTA-T ALGORITHMS, which is why I flagged it in the PR. It will be important to validate that all algorithms expect no Moon Secular Acceleration correction after 1955.5 Finally, for dates after 1962-01-20, the JPL Horizons algo compute TDB - UTC and use the DE440 file to get TT-TDB and eraDat() to get leap seconds. If DE440 is not selected in the Ephemeris settings, then we loose a small correction of around +/-2 milli-seconds (i.e. TT - TDB) which is why we originally fail the CI test since my criteria was +/- 10 microseconds of matching JPL Horizon app. I have update testDeltaT.cpp and testDeltaT.hpp to bypass the test if file linux_p1550p2650.440t is not found. So far Linux, macOS and Windows are successful. So it looks like the automated CI test just didn't have access to the file for Linux and macOS... but Windows had (very strange...) |
No test has. Windows CI simply doesn't run any tests. |
|
|
||
| #define DE430_FILENAME "linux_p1550p2650.430" | ||
| #define DE430_FILENAME "linux_p1550p2650.430t" // (SS) 2025-11-27 use "t" version to access TT - TDB corrections (if needed) | ||
| #define DE431_FILENAME "lnxm13000p17000.431" |
There was a problem hiding this comment.
Are there also updated files for the long-time solutions where this IMO actually matters more?
There was a problem hiding this comment.
To @gzotti :
...Are there also updated files for the long-time solutions where this IMO actually matters more?
Based on the JPL Public FTP File Browser site https://ssd.jpl.nasa.gov/ftp/eph/planets/Linux/, I could only find the de430t and de440t files with the TT - TDB information. The de431 and de441 'long-time' ephemerides solutions are 'non-t ' version, which means that these long-time solutions do not include the TT - TDB data. Note that the TT - TDB data has a maximum impact of +/- 2 milliseconds on DeltaT.
| // JD_UT >= 1458077.5 or 721BC-01-01 00:00 | ||
| // JD_UT < 2437684.5 or 1962-01-20 00:00 | ||
| // deltaT values correspond to DeltaT = TDB - UT1 in this time range. | ||
| if (1458077.5 <= jDay && jDay < 2437684.5) |
There was a problem hiding this comment.
should be "else if... " (also in all cases below) for clarity.
|
After re-reading our exchange, I wonder (again) about the global benefit of this. Avoid a 2 milliseconds error when our time granularity is around 1 second? (We are using a double-precision JD date that loses accuracy far from JD0, and users cannot easily set sub-seconds (or run slow-motion time for satellites or such).) Are there any long-time errors/effects that are improved by this? I am willing to run the table results for F.9 again. What is very apparent is a slight but noticeable difference in positions in b2500. Will this be (noticeably) improved by this PR? Whenever you wrote "someone will have to change the user guide entry about XY", please be that someone. Apparently it is part of your current activity. The whole ΔT explanation deserves a much clearer discussion in the User Guide, including this "no more changes after 1955.5" which was new to me. Also discuss the relevant time frames. It is a difficult subject with all the definition changes, GMT, UT, UTC, UT1, UT0? Why do we have so many models? (We try to be compatible with many historical studies.) My old classic source book, Meeus' "Astronomical Algorithms", writes that for all practical purposes we can set TT==TDT=TDB, with a difference of 0.0017s, is that your 2ms? I would seriously prefer to fetch the ERFA via CPM if that was possible, and even keep it optional with some sane fallback, instead of ingesting a full copy into this repo. Maybe the ERFA creators can provide a version that can more easily be fetched&built by CPM? Else, do we need everything or just 2-3 relevant source files that must be hosted, together with a README that provides the info where the rest of ERFA can be found? (Uh, is it OK, license-wise, to just grab a part of it?) |
|
To @gzotti : _...After re-reading our exchange, I wonder (again) about the global benefit of this.... Are there any long-time errors/effects that are improved by this? It is not clear to me what you mean by 'this'... Are you just referring to the 2 milliseconds correction ? You known that this correction is only used for dates after 20 January 1962, and therefore it has nothing to do with historical dates such as b2500. For these historical dates, the revised JPL Horizons DeltaT algorithms was developed to match what JPL HORIZONS app generates with a difference below +/-0.1 second. As I pointed out in issue #4670, The current JPL Horizons DeltaT algorithm do not match the JPL app, and is limited in date range, since it does not even cover post telescopic times. Maybe, the current algorithm has been developed only for estimating Delta T for historical dates...(i.e. pre-telescopic time), and therefore there is no need to extend the date range up to present time....When I tackled this effort, my objective was to reproduce the JPL HORIZONS app as closely as possible, not focusing on specific date ranges. ...Whenever you wrote "someone will have to change the user guide entry about XY", please be that someone. Please be aware that I did work on the Stellarium User Guide to update the section relevant to the JPL Horizons model. Refer to commits 43d25e2 and d763df0. The issue with "no more changes after 1955.5" that I raised on 29 Nov 2025 (see above) is a tricky one... It does make sense to me, but I am not an expert in this field, so I was turning to your (expert) opinion on it. If this is also new to you... oops... The 2 msec is indeed the difference between TT and TDB. JPL HORIZONS app take this into account for dates after 20 January 1962, which is why I included it as part of the whole algorithm revision. It could be done easily, since JPL provides the necessary files (DE430T and DE440T) and Stellarium code has the logic to extract the TT-TDB data in these files. But if we want to ignore this aspect to Delta T, we can certainly do it, since it is only a contribution of +/-2 msec. For the issue of the ERFA library, I prefer not to get involved as it is license-related. However, if you resolve this license aspect, why do it for only 1-3 files.. go for the whole library. There is quite a lot in this library, that could become very useful for future development of Stellarium... |
|
"this" = "this huge pull request that promises to improve present-time computations by 2 milliseconds". Now, aren't there long-time effects from using different n_dot values? |
|
To @gzotti : ..."this" = "this huge pull request that promises to improve present-time computations by 2 milliseconds". The original JPL Horizons algorithm for computing Delta T has its year bounded between b3000 and 1620. It is therefore giving 0 seconds as a return value for Delta T for present-time. The revised algorithm is not just improving present-time computations by 2 milliseconds. It computes a Delta-T and if it has access to the 't' version of the JPL DE files can include the difference between TT and TDB which amount to a small correction of +/- 2 millseconds max. As I wrote before, this can be dropped, but otherwise the code was already there in Stellarium, so it was not a big deal to add this correction to my point of view. ... Now, aren't there long-time effects from using different n_dot values? My understanding on the various n_dot values is that those are dependent on the ephemerides used. Each ephemerides (VSOP87, ELP200-82B (for the moon), DE430, DE440, etc) has its own n_dot value. Now a set of equations (or model) is provided to compute a Delta T. These equations are based on a certain a moon secular acceleration parameter (i.e. n_dot). If the assumed n_dot for the equations (or model) is different from the n_dot of the ephemerides, then you have to apply a correction. This is done by getMoonSecularAcceleration(). Having said that, I learned, while working on the JPL Horizons algo, that we should not apply this correction for years after 1955.5, but otherwise it must always be applied before year 1955.5. So yes, using different n_dot values matter and will have long-time effects the further back you go in time. Georg, I understand that this is a huge PR. I worked on it because I was curious to learn on this subject and thought I could contribute to Stellarium by the same token. BUT it remains your call to move on and incorporate it in an official release. I am not forcing anybody... Best Regards, Sylvain Simard |
|
I just gave it a try. Comparing the results from Appendix F.9 for Venus, b2500-01-01. Our current solution is perfect in RA and deviates by half an arcsecond in δ, the new results are far worse. On the other hand, at least this single reading is better for the Moon. So, some good elements may be included, but not everything became better. This needs more testing. I care much more for a bad minute (that shifts a solar eclipse zone) in antiquity than 2ms in the present era. |
|
Just to confirm we see the same thing...and few more questions on this subject. Here is what I get from the JPL HORIZONS app and from Stellarium with JPL_horizons_algorithm_deltaT mods of this PR for Venus @ Earth geocentric at JD(UT) = 808298.5, i.e. b2500-Jan-01, 00:00 (Julian Calendar). From JPL HORIZONS app: From Stellarium + JPL_horizons_algorithm_deltaT: Note : DE441 is selected for the ephemerides for both JPL app and for Stellarium According to Stellarium User Guide, the entry in Table F.2 of section F.9 report the following results from Stellarium 0.21.2b From Stellarium 0.21.0b: There is indeed now a difference 1.52s in RA and 9" in DE, while from the previous version of Stellarium there was no difference in RA and only 0.5" for DE for Venus at this specific date. It should be noted, however that Delta T has significantly been improved, since the mods give a difference of only 0.025s for this parameter, while version 0.21.2b had a difference of 954.5s Do you see the same numbers as well ? Here are the support file of my observations, including three screenshots to show my setting of Stellarium, I have the Config.ini and Log.txt file if needed. horizons_results_Venus@Earth_b2500-Jan-01_00-00.txt 
I conceed that for this particular data point, the predictions are worse, but we are not talking of minutes or arc-minutes, but rather few seconds and arc-seconds. Now that Delta-T is similar between Stellarium and JPL HORIZONS app, we could rule out this parameter as an explanation of the differences. I would personally suspect that models of Precession and Nutation may be different between JPL HORIZONS app and Stellarium. According the the user manual of JPL HORIZONS (https://ssd.jpl.nasa.gov/horizons/manual.html#refs) JPL uses the Owen Precession model before 1799-Jan-1 and after 2202-Jan-1, and the Nutation (IAU 1980) model, as per the following references:
Does Stellarium uses similar models ? If not, could the different models of Precession and Nutation explain the variation in RA and DE seen between Stellarium and JPL HORIZONS app ? I do not know enough on the inner working of Stellarium to judge on this point. I found a paper from Owen in 1991 about his Precession model (a summary of his thesis) that I could share if you don't have it. NOTE: I find strange that JPL uses Precession and Nutation models developed in year 1990 and 1980 for their app. I would have thought that new models have been developed since then.... Please feel free to share your thought on this ... if you are still interested on this subject ;-) |
|
We are using Vondrák 2011 (with typofix correction from 2012) for precession over 200.000 years. And yes, I am interested in improving things, but while I mourned the lack of DeltaT 15 years ago, I was not involved in the original implementation, and admittedly some details may have evaded my attention. If our n_dot values for some models are off, at least those can be taken from here and fixed. But sorry, I am afraid it's not to be rushed into 26.1 at this stage, as we still have other PRs and bugfixes to work on. |
|
I checked the first 4 entries of Table F.2 (Sun @ b2500 and b2499) with Stellarium 26.0+ JPL_horizons_algorithm_DeltaT mods and in all cases the new predictions of RA and DE are worst than those from the previous version of Stellarium (21.2b). Not by much, but still worst. On the other hand, the Delta T is spot on with the JPL HORIZONS app (to +/- 0.1s) , while Stellarium 21.2b was almost 1000s smaller at this historical epoch. I wonder if the old JPL Horizons algorithm of Stellarium (i.e. before the mods) was not somehow compensating something else, and this way providing Stellarium predictions of RA and DE closer to those of JPL HORIZONS app. Stellarium is using a different (and more recent) Precession model - Vondrák 2011 (with typofix correction from 2012). JPL HORIZONS app is using an older Precession model (Owen 1990). An effort is needed to evaluate the impact of this difference between the two models on historical RA/DE predictions. I'll see if I can contribute something useful on this. aspect... We should all keep in mind that seeing differences between Stellarium and JPL HORIZONS app is to be expected if different models of Precession (and Nutation and Aberrations) are used in the two implementations. And, I should add more... that it does not mean that the JPL HORIZONS app predictions are necessarily 'more true' than those of Stellarium. JPL uses a Precession model two decades older that Stellarium.... strange isn't it! Anyhow, I fully understand that there are other PRs and bugs that require your attention. As you said, at this stage we should not rush this PR in release 26.1. Stay cool |
|
Yes, I had seen the divergences in the displayed values of DeltaT, esp. when the coordinate fit so well. This is a point that should give me headaches. It has to do with the chicken/egg problem between n_dot and DeltaT. One earlier contributor was important and may be able to solve this mystery. I will try to catch up with him. |
|
To @gzotti : I've done few more checks of differences in RA/Dec and Delta-T between Stellarium 26.1+ and JPL HORIZONS app for the Sun and Venus. While this is still work in progress and quite limited in test coverage, I can report the following: a) Stellarium with the new JPL Horizons algorithm for Delta-T matches JPL Horizons app for RA/Dec for Astrometric coordinates. That is the RA/Dec (J2000) field without Parallax, Aberrations and Nutation checkbox in configuration panel gives exactly the same values between Stellarium and JPL HORIZONS app. I could only check a handful of dates between 9999BC-Mar-20 (UT) to 2026-Jan-01 (UT) so far but they match perfectly. The value of Delta-T matches to the precision that I expected (i.e. better than +/- 0.1 second between 9999BC up to 1960 and then better than +/- 10 micro-second from 1960 to today. I am documenting these results in an Excel file. This is an improvement over Stellarium 25.4 which does not match RA/Dec (J2000) for Astrometric coordinates. Moreover, 25.4 has a valid range limited between 3000BC and 1620, while 26.1+ can cover 9999BC up to today. b) For Apparent coordinates, which correspond to the RA/Dec (on date) field I am seeing a trend over time. Differences between Stellarium and JPL HORIZONS app are small (less than 0.01 sec) near the epoch J2000.0 but increasing up to -10 seconds at 9999BC. It does not seem to be linear over time and we could get positive and negative values of differences. More data points are needed to have a better idea of this behavior. The Apparent declination has a similar trend, being small (0.1 arcsec) near J2000.0 but increasing up to +/- 50 arcsec at 9999BC. Interestingly, it varies over a year in some kind of oscillating pattern, possibly being depending on the sine or cosine of the declination. I think, having Astrometric coordinates matching perfectly is a good sign. Moreover, it is a clue that differences between Stellarium and JPL Horizons app for Apparent coordinates may be simply due to differences in models used for Precession and Nutation. I personally suspect that the two Precession models may simply drift over time, and what we see as differences in RA/Dec is the result of this drift. For sure there ought to be a contribution from the different Nutation models as well (IAU 1980 for JPL HORIZONS app and IAU2000B for Stellarium), but I suspect this will be a second order effect. More work is needed on this hunch before I become fully convinced. Moreover since I was able only to check few data points for the Sun and Venus a new approach using scripts will be needed to automate data capture and analysis. Note: I also contacted Jon Giorgini (JPL) by email for support on the Owen(1990) Precession model used by JPL HORIZONS app. Although he was not able to help much on this side, he gave me some comments and advices. See transcript of the email in attachment. Re JPL Horizons app and Delta-T algorithms - Jon.D.Giorgini@jpl.nasa.gov - 2026-03-11 1408.md Out of curiosity, are differences of +/- 10 sec in RA and +/- 50" in Dec between Stellarium and JPL Horizons app important for archeoastronomy ? |
|
To @gzotti . Here is another interesting comparison between Stellarium 25.4, Stellarium 26.1+JPL Horizons Algorithm and JPL HORIZONS app. Observation of the Earth from an observer located on the Sun. a) Measurements of RA/Dec (J2000) with and w/o Parallax+Abberation+Nutation These observations match perfectly between JPL Horizons app and Stellarium 26.1+JPL Horizons Algorothm, but not with Stellarium 25.4. Also, these measurements are not affected by Parallax, Aberration and Nutation settings in the Configuration panel. Somehow these effects do not impact these measurements. Possibly because we are on the ecliptic (TBC). I checked these measurements over a time span covering 9999BC and today. We have to keep in mind however that for Stellarium 25.4 time outside 3000BC to 1620 are not valid (in principle) as documented in the User Guide and Stellarium. This is another strong indication that 26.1+ is capable to reproduce perfectly some simulations of the JPL HORIZONS app, while 25.4 can't. In addition, this kind of measurements seem to be independent of the Precession/Nutation correction, so therefore independent of the models used for these effects. Therefore any residuals differences between 26.1+ and JPL HORIZONS app for other observing configurations (e.g. observing the Sun while located on the Earth) is very likely linked to the differences in these models (TBC). Here are some files and screenshots documenting these results: a) JPL HORIZONS app results: b) Log and config.ini files of the two test cases: config.ini_26.1+.txt config.ini-25.4.txt c) Screenshots from test cases with Stellarium 26.1+ JPL Horizons Algorothm,: 
d) Screenshots from test cases with Stellarium 25.4: 
|
|
Thank you, this is at least promising. However, magically, several other PRs suddenly appeared which are easier to handle but also need time, in addition to an important bit I am biting heavily at. As said, we need to better understand the relation between DE4xx and our understanding of how ΔT is computed and applied. Yes, precession/nutation is a likely candidate for slight deviations. Sorry, I am not able to do all this within a week (I do have a daytime occupation...), but please send further observations if you find something relevant. I will continue to respond after 26.1 is out. |
|
To @gzotti. Here is an update of my progress. I worked on software tools (Stellarium script and Excel) to be able to compare various Stellarium and JPL Horizons predictions (RA, Dec, Astrometric and Apparent - with and without Aberration and Nutation, as well as Apparent Ecliptic Longitude and Latitude) over time spanned (mainly) between today and 4000BC. I did look at data at year 9999BC, but since Stellarium stop using its Nutation model after 4000BC, it does not make sense to compare Stellarium and JPL between 9999BC and 4000BC. I processed data only for Venus for now, but could eventually do it for the other planets and moons as listed in Table F.2. Here is the Stellarium script that I used to capture Stellarium data (it is a slight modification of your own script). Remove the .txt file extension to get the script. After capturing the various output.txt files and a similar one from JPL Horizons app covering the same time span at the same # of days per step, I combine the data into an excel spreadsheet and generate figures of differences between them. Here are few figures for a Time Span between 4000BC and 2026, at 100 days per step.
Well that's it for now. My next steps are as follow: A) Develop C/C++ code (using SOFA and JPL DE441 Ephemerides) to reproduce Stellarium predictions for Venus using Stellarium Aberration, Precession and Nutation models. B) Add Gravitational light deflection to this code and see if this account for differences seen at 3) above C) Modify the C/C++ code to use JPL Horizons Precession and Nutation models instead of Stellarium models. Compare results with A (and B) to see if we reproduce the drift seen in 4) above. D) Try to explain results in 5) with this new code. For example, study how Ecliptic Latitude is sensitive to aberration models and gravitational light deflection. This is a lot of work and I could post my progress here if you're OK with this. Let me know otherwise.
|
|
To @gzotti. Here is another update of my progress. As I described in the precedent comment I had come up with a program for investigating why Stellarium predictions of on-date Apparent RA/Dec and Ecliptic Longitude/Latitude of a planet position (Venus for the moment) where at odds with the JPL Horizons application in the order of few 10's of arcsec over the time range between 4000BC and today. While this objective remains, I have slightly changed my program of investigation. A) Develop C/C++ code (using SOFA and JPL DE441 Ephemerides) to reproduce JPL Horizons predictions for Venus using JPL Precession and Nutation models, as well as impact of Light Deflection by the Sun and annual aberration. B) Disable Light Deflection and assess degradation in the predictions from JPL Precession & Nutation models. This step is needed to get closer to Stellarium implementation, since Light Deflection is not handled by Stellarium. C) Replace the JPL Precession and Nutation models by IAU models (from SOFA) which (in principle) should be equivalent to those used by Stellarium. I used IAU Precession model iauLtpb() which should be the equivalent to the long-term Vondrak (2011-12) model, and iauNum00b() for Nutation which should be the same as Nutation 2000B used by Stellarium. D) Discuss results of RA/Dec in the true on-date Equatorial and Equinox reference plane and Longitude/Latitude in the on-date Ecliptic reference plane between JPL Horizons and Stellarium and between JPL Horizons and this C/C++ code. While this program is slightly different from the precedent comment, the spirit remains the same. Step A) After many tests and try-out I more or less decipher how JPL Horizons compute apparent RA/Declination as well as Longitude and Latitude coordinates for planet Venus. JPL Horizons seems to use different Precession & Nutation models over 5 time zones: Zone #1: From 9999BC-Mar-20 to 1799-01-01 they use the Owen(90) long-term theory for the Precession model. This theory is implemented by Chebyshev coefficients to compute angles moving between the mean J2000.0 equatorial plane to/from the invariable plane of the solar system and to the final on-date mean equatorial plane. Then JPL Horizons uses Nutation model 1980 to add Nutation. I personally add a small Rz rotation (about 100 mas) to the final vector to better center results and be closer to JPL Horizon predictions. My code do not fully reproduce the location of the origin in the on-date referential frame and this small rotation correct for this. Zone #2: From 1799-01-01 to 1962-01-20 they use the Precession/Nutation models 1976/80. I also add a small Rz rotation (about 100 mas) to the final vector to better center results and be closer to JPL Horizon predictions. Zone #3: From 1962-01-20 to Today + 78 days. Here they claim to use the Precession/Nutation models of 1976/80 with EOP corrections. I could not find the EOP correction files on JPL site, but JPL Horizons manual claims that IAU models 2006/2000A reproduce their results with only a shift of about 53 mas in the origin of the on-date equatorial plane. So I used those Precession/Nutation models from SOFA with an added Rz correction of about 53 mas to account for the shift in origin. I developed a parabolic fit of the correction for a better match with JPL Horizons within this time zone. Zone #4: From Today + 78 days to 2202-01-01 they simply re-use the Precession/Nutation models 1976/80. I also add a small Rz rotation (about 100 mas) to the final vector to better center results and be closer to JPL Horizon predictions. Zone #5: From 2202-01-01 to 9999-Dec-30, they re-used the Owen(90) long-term theory for precession and model 1980 for Nutation. Again I add a small Rz rotation (about 100 mas) to the final vector to better center results and be closer to JPL Horizon predictions. Here are the final results in graphic format covering the whole period between 9999BC to 9999 for differences between my code and JPL Horizons
As you can see, Right Ascension and Longitude are within +/- 20 mas and Declination and Latitude are +/- 80 mas. This is over the whole time period provided by JPL Horizons app (i.e. 9999BC to 9999AD). Step B) If we disable light deflection by the sun in my C/C++ code we get a degradation by a factor of 10. Here are the same results with the same precession and nutation models, but with only light deflection turn off
Step C) I then proceeded to replace in my custom C/C++ code the Precession & Nutation models of JPL by those from SOFA which should correspond with those of Stellarium, i.e. Vondrak (2011-12) for Precession and IAU 2000B for Nutation. I kept Light Deflection turn-off since this is what Stellarium does. I also limited the time period between 4000BC and 2026 since Stellarium does not apply Nutation before 4000BC. Here are another series of graphs showing predictions from my C/C++ code and from Stellarium. Note that those from Stellarium are (more or less) the same ones shown in the precedent comment.
Step D) The Apparent RA and Longitude graphs from my C/C++ code using SOFA code for Vondrak (2011-12) and IAU 2000B are quite different from those from Stellarium. However, while the form of the curves are very different, the deviation from JPL Horizons is more or less of the same order for the absolute magnitude. We get roughly 2 seconds in RA with JPL Horizons (though Stellarium is +2 sec and my code is -2 seconds....) and +/- 15 arcsec for declination. Longitude covers about 35 arcsec for Stellarium and 30 arcsec for my code, while Latitude get to +/- 25 arcsec at 4000BC starting close to 0 arcsec in 2026, while Stellarium is a constant +/- 10 arcsec over the whole period. From these results I think we can conclude that Precession and Nutation models being different between Stellarium and JPL are the main cause of the discrepancies between them. Although my simple C/C++ code was not able to reproduce the behavior of Stellarium it does give results of the same order of magnitude. I WAS SURPRISED, HOWEVER BY THE RESULTS IN THE SHAPE OF THE CURVES. I WAS EXPECTING TO GET MUCH CLOSER TO STELLARIUM BY USING THE SOFA VONDRAK AND 2000B ROUTINES. THIS REQUIRE MORE INVESTIGATION AND A BETTER UNDERSTANDING ON MY PART ON HOW STELLARIUM COMPUTE APPARENT POSITIONS OF PLANETS IN ITS GUT... The good side of this investigation is that we now have C/C++ code that implement how JPL Horizons app computes apparent position for planet Venus to an accuracy of +/- 20 mas in RA and Longitude and +/- 80 mas in Declination and Latitude, or better than +/- 0.8 arcsec if we do not include Light Deflection by the Sun. It will be your call if you want to incorporate this code into Stellarium. As we mentioned earlier, more accuracy may not be required for Stellarium considered who the users are... if you decide to incorporate it, I remain available to support this effort, though help on your side will definitively be required as I have obviously not a good grasp of the inner working of Stellarium in these matters. Any comments will be welcome. |
|
Oh boy, that is serious stuff! I simply need the 3-6 months well-rested daytime to incorporate all this, my current 2-hours-in-the-evening-plus-weekends self-exploit on top of doing other interesting but unrelated things is simply not sustainable. But honestly, this work has potential, until maybe mid-2027 (IAU GA Rome?), of
All this on top of the other user/community requests, satellite data file format switchover, adding Alpaca into an old plugin, cultural issues, artistic requests, ... |
|
That's quite a program you have in mind! and you're right it will take a 'good' effort to bring Stellarium within tens of milli-arcsec differences with JPL Horizons. Is it worth it ? if you ask me, I'll say 'yes' because I am on the 'more accurate' side, but the 'good enough' side has their say too... In the upcoming weeks I'll study your second bullet. I am curious to understand why I got this kink in the differences in RA between IAU routines and Stellarium. I may not be able to provide feedback until June as I will be on vacation in May (visiting Morocco and hoping to see the milky way in real, this time ;-). I used SOFA for this investigation but ERFA is the same code. The Owen90 Precession model is not very complicated once you know how to compute the needed angles using Chebyshev representation. Similarly adding light deflection to Stellarium is easy using the available SOFA/ERFA routine. Note that, I found by chance, that JPL Horizons simply turn-off light deflection when the elongation of the planet Venus reaches beyond 179.36 degree. Venus is then passing behind the Sun, so this avoid unrealistic modeling of the SOFA iauLd() routine. Anyhow, I'll be glad to support you if you decide to engage in this endeavour. Considering your limited availability, you should consider to act as the program lead, meaning you need to build-up a team of programmer(s) who knows the inner working of Stellarium and supervise their work instead of doing it. My 5 cents. Stellarium is already a very nice product but we can always improve on it. That's the beauty of it. Stay cool. |
|
"need to build up a team... " Sure. I was happy that we became a stable 3 after 10 years... Now maybe 4 and 5? |
|
I read again the conclusion of Vondrak et al (2011) paper: They clearly state that this model is only accurate to few arcsecond during the historical period. If we assume that the historical period meant 'back to 3000BC' then this is more or less what I saw in my investigation... +/- 10 arcseconds back to -50 Julian centuries from J2000.0 (i.e. back to 3000BC). Note that Vondrak et al (2011) are fitting a numerical integration of solar dynamic over +/- 200 kyrs to represent general precession and obliquity by cubic polynomials plus 8 to 14 periodic terms. They might achieve good fit over +/- 200 kyrs, but when compared to the numerical integration they state that they only achieve few arcsecond during the historical period.... Owen (90) is also based on a numerical integration, but over +/- 500 kyrs this time. However in his case he uses Chebyshev polynomials (of order nine) to map the precession angles within time period spanning 80 centuries each. That is, he has 125 tables of 10 Chebyshev coefficients for each angles. Although the numerical integration is a bit outdated mapping it with Chebyshev polynomials over many short time periods (i.e. 80 centuries each) should results in much less errors during the historical period. I cannot say by how much less, but I would guess that using Chebyshev polynomials is the main reason to get better fit with a numerical integration. After all, DE44x ephemeris are based on Chebyshev polynomials and they are definitively a better fit to a numerical integration of solar bodies position/velocity than the VSOP87 model for example. We should not be surprised, then, to see few 10's of arcsecond discrepancies with JPL Horizons back in 4000BC (beyond the historical period) when using Vondrak (2011-12) precession model. It is the claimed accuracy of this model. Does this make sense ? |
|
Totally. When I started, Meeus gave a formulation for ε that is applicable for J2000±10ky. When I found Vondrák's paper I was happy to have a solution that provided precession to a few arcseconds over far more (±200ky) than I found really necessary for (pre)historical purposes (±25.000 years). My overall accuracy goal was "good sub-arcminute", matching or surpassing naked-eye requirements, but with known effects (nutation, aberration) that could sum up to larger errors taken into account. I did not look further, and did not know about Owen90 before you mentioned it here. When you want to enter milliarcseconds territory, this requires another even more thorough review of all computing steps (including the sytematic offset between VSOP87 and "ecliptic frame of J2000.0" that we have somewhere as tiny rotation matrix. Is it applied correctly everywhere? The effect is negligible for me, but important for your goal. Again, finding a reference solution to compare against is the major challenge. ), and further refinement of the components. Presumably Owen90 is available somewhere, and when using it as switchable alternative to the Vondrák formulation, we could compare its processing impact. If it covers an even longer period, some folks might want to study stellar proper motion for even more extended periods (with all problems that come from extrapolating of just ~10 years of observation...). With such extended precession, people will ask about planetary orbits 250.000 years back as well. There are published orbital evolution curves going back tens of millions of years concerning overall stability of the Solar system, but I have found no numerical models of even mean orbital elements that could be used for approximate Kepler orbits for times before DE441. (Not that anybody really has a practical use for knowing about constellations, seasons or planetary conjunctions back then. And sure, "Skies over T.rex" remains a Hollywood requirement...) |
|
I found a short description of Owen90 in a paper published in 1991: owen1991.pdf The full thesis (which includes the long-term theory and all the Chebyshev coefficients) can be downloaded from the Internet Archive website: https://archive.org/details/theoryofearthspr00owen I also reached out to Dr Owen recently via its publisher, the cie Wiley. He answered back to my first enquiry, providing a very quick resumé of his theory that he developed in 1990. He did not seem to be aware that JPL Horizons app uses his long-term theory since he mentioned in his email that he thought the app uses the IAU precession model ... very strange since he worked for JPL... I tried to reach out again with more questions, but he did not answer them yet. Here is the exchanged emails, so far: To Kurzman, Brett of Wiley.com.txt I coded only three intervals of 80 centuries each in my C/C++ code to map time spanned between 9999BC to 9999AD; intervals 62, 63 and 64 of the 125 intervals covering the +/- 500 kyrs period. So yes we can extend to more than +/- 10kyrs the precession model by encoding more intervals. Keep in mind, though, that this will be good only for the precession model. Note: JPL Horizons uses the Nutation model of 1980 along with Owen90 precession model. They do not use the IAU 2000B model. If we integrate Owen90 Precession model into Stellarium as a switchable alternative to Vondrák formulation, we might also have to switch to the 1980 Nutation model to keep thing compatible.... to be further investigated... |
Need a minor modification to getDeltaTByJPLHorizons() to adjust an upper bound in a else if test case
…nsition JD 2437683.99999 to JD 2437684.0
4 participants
This pull request addresses a single issue:
Update the JPL Horizons Delta-T algorithm to be in-line with JPL Horizons application. The revised algorithm is based on initial guidelines generously provided by Jon Giorgini (JPL) and from information available within the JPL Horizons app user manual:
Jon Giorgini
Senior Analyst
Solar System Dynamics Group (392R)
Navigation and Mission Design Section
Jet Propulsion Laboratory
Pasadena, CA USA
Jon.D.Giorgini@jpl.nasa.gov
source : JPL Horizons app user manual: https://ssd.jpl.nasa.gov/horizons/manual.html
Description
Changes required to three (3) .cpp files: StelCore.cpp, StelUtils.cpp, StelGuiItems.cpp
Changes required to two (2) .hpp files: EphemWrapper.hpp, StelUtils.hpp,
Changes required to one test case: testDeltaT.cpp
As documented in issue #4670, the motivation is essentially to resolve unrealistic deltaT values (i.e.. 0.0) outside the time range (3000BC to 1620) of the previous algorithm and to improve accuracy and concordance with JPL Horizons app.
This change requires access to a routine which can compute leap seconds from 1962-10-20 to present. Software from the ERFA library was integrated within Stellarium and build by CMake as a static external library. It still remains to the Stellarium team to decide if this integration is legally valid.
Screenshots (if appropriate):
None
Type of change
NOTE Stellarium User Guide needs to be revised to clarify the need for the use of new JPL DE* files and how to properly use Stellarium and interpret Delta-T results.
How Has This Been Tested?
Section relevant to JPL Horizons algorithm within testDeltaT.cpp has been revised. Reference data were extracted from JPL Horizons app and compared with values obtained by Stellarium. Differences between Stellarium and JPL Horizons are as follow:
from 9999BC to 1962-01-20 : les than +/- 0.1 seconds differences
from 1962-01-20 to 2650-01-25 : less than +/- 10 microseconds differences
from 2650-01-25 to 9999AD : less than +/- 2 milliseconds differences
It is not expected that this change affect other area of the code. HOWEVER, a modification to getMoonSecularAcceleration() limiting this effect to JD < 1955.5 has been added which may impact other DeltaT algorithm if they are not based on this same assumption. This need to be investigated further.
Test Configuration:
• Operating system : Windows 11 Pro, 25H2, OS build 26200.7171
• Processor : 13th Gen Intel(R) Core(TM) i7-13700K (3.40 GHz)
• Graphic Card : NVIDIA GeForce RTX 5050 Ti
• Driver : 581.80
• Display : ASUS PA329CV, 3840x2160 | 60 Hz | Highest (32-bit) | HDCP
Checklist: