Repository navigation
another (much smaller) error in TT to TDB calculation in Time #1938
Description
Hi All (esp @taldcroft and @mhvk),
So @demorest and I found another issue with the way that the conversion to TDB is happening in _get_delta_tdb_tt(). It turns out that there is another documentation bug in the SOFA examples where one of the parameters was being passed off by 0.5 day. That caused an of-order 10 us error in TDB values.
With the following patch (which has the "ut" calculation written in the same style of the SOFA Time guide example) the TDB values we get are now equivalent to those given with TEMPO2 (using this particular TDB model) at the few nanosecond level. So that is excellent. TEMPO2 does use a different TT-TDB model by default (from Irwin & Fukushima 1999) which gives a different TDB from astropy with what appears to be a secular term at the tens of ns level difference and a few ns variation on shorter timescales on top of that.
But it appears that the time code now in astropy (and therefore SOFA/ERFA) is good enough for high-precision MSP work! So that is really cool.
Happy New Year,
Scott
diff --git a/astropy/time/core.py b/astropy/time/core.py
index f1aca8f..b38f188 100644
--- a/astropy/time/core.py
+++ b/astropy/time/core.py
@@ -811,8 +811,10 @@ class Time(object):
# pretty close (few msec?), assume TT.
njd1, njd2 = erfa_time.tt_tai(jd1, jd2)
njd1, njd2 = erfa_time.tai_utc(njd1, njd2)
- # TODO: actually need to go to UT1 which needs DUT.
- ut = njd1 + njd2
+ # ut is the UT1 fraction of a day, referenced where 0.0 is
+ # midnight. Therefore we need to add 0.5 day to the UTC JD.
+ # The difference between using UT1 and UTC is tiny.
+ ut = np.fmod(np.fmod(njd1, 1.0) + np.fmod(njd2, 1.0), 1.0) + 0.5
And here are the 2 new "correct" lines from the SOFA Time example document on page 25 (which is used in the unit tests):
TDB 2006/01/15 21:25:42.684373
TCB 2006/01/15 21:25:56.893952
Activity
demorest commentedon Dec 31, 2013
Also note that
eraDtdbdoesfmodon theutargument interally as well, so all thenp.fmod's in the patch might not really be necessary.scottransom commentedon Dec 31, 2013
Yup. That is true. Although I kind of like using a very similar code to what is in the SOFA example just so people can match what is going on between the two.
mhvk commentedon Dec 31, 2013
@scottransom, @demorest - thanks! I used a slightly different solution to it in #1939.
Correct UT in calculation of TDB; closes astropy#1938
Merge pull request #1939 from mhvk/time-correct-ut-in-tdb
taldcroft commentedon Jan 3, 2014
@demorest @scottransom - just to be absolutely sure, in the original issue report the discrepancy was said to be around 10 us, but the actual change in our test case was more like 1 us.
https://github.com/astropy/astropy/pull/1939/files#diff-1ad7b3f73aab523e9f91b15536384814L575
Is that understood?
As an aside, is doing testing with digits out to the microsec level sufficient, or should we add another digit or two where meaningful? These serve as regression tests and can catch very small changes or possibly platform-dependent issues. Just wondering.
demorest commentedon Jan 3, 2014
@taldcroft - the discrepancy in this report had an annual sinusoidal periodicity with a peak-to-peak amplitude of about 5 us. For a single test point to change by about 1 us sounds reasonable to me. I think adding another couple digits to this test would be a good idea.
Merge pull request #1939 from mhvk/time-correct-ut-in-tdb
Merge pull request #1939 from mhvk/time-correct-ut-in-tdb