Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Better document that newer Skyfield versions carry more recent ∆T data, which affects observatory position and orientation #962

Open
Dave-II opened this issue May 3, 2024 · 2 comments

Comments

@Dave-II
Copy link

Dave-II commented May 3, 2024

I have run across a discrepancy between 1.4.6 and 1.4.8 in satellite RA/Dec reported on a satellite-observer difference. Here is a example piece of code run and the result under the two versions.

import datetime
from skyfield.api import wgs84, load, EarthSatellite
 
def test() :
    name = 'COSMOS 2506'
    tle1 = '1 40699U 15029A   24111.15280151  .00000755  00000-0  19062-3 0  9995\n'
    tle2 = '2 40699  98.1576 165.8491 0008633  79.0387 281.1781 14.53727601468365\n'
    sat = EarthSatellite(tle1, tle2, name, None)    
    obs = wgs84.latlon(43.192909, -81.315655)
    dt = datetime.datetime(2024, 4, 20, 2, 13, 27, 320000) - datetime.timedelta(seconds=0.846)
    t = load.timescale().utc(dt.year, dt.month, dt.day, dt.hour, dt.minute, dt.second +dt.microsecond / 1000000)
    difference = sat - obs
    topocentric = [difference.at](http://difference.at/)(t)
    ra, dec, distance = topocentric.radec()
    print(t)
    print(ra._degrees)
    print(dec._degrees)
    print([distance.km](http://distance.km/))

V1.4.6:

121.58953952934377
56.35798447940125
796.8656666304348

V1.4.8:

121.59166803867389
56.35877469948038
796.8583740425258

@brandon-rhodes
Copy link
Member

Hey, Dave! If you'll check the version 1.47 changelog entry:

https://rhodesmill.org/skyfield/installation.html#v1-48-2024-february-7

— you'll see:

Skyfield’s internal table for the ∆T Earth orientation parameter has been updated, so that its predictions now extend to 2025-01-18.

Because the more recent versions of Skyfield have better data for where the Earth was actually pointing on April 20, the answers it gives for positions relative to a ground station are slightly different. You can cancel out this effect by setting ∆T artificially to a static value, like zero:

import datetime
from skyfield.api import wgs84, load, EarthSatellite

def test() :
    name = 'COSMOS 2506'
    tle1 = '1 40699U 15029A   24111.15280151  .00000755  00000-0  19062-3 0  9995\n'
    tle2 = '2 40699  98.1576 165.8491 0008633  79.0387 281.1781 14.53727601468365\n'
    sat = EarthSatellite(tle1, tle2, name, None)
    obs = wgs84.latlon(43.192909, -81.315655)
    dt = datetime.datetime(2024, 4, 20, 2, 13, 27, 320000) - datetime.timedelta(seconds=0.846)
    ts = load.timescale(delta_t=0.0)
    t = ts.utc(dt.year, dt.month, dt.day, dt.hour, dt.minute, dt.second +dt.microsecond / 1000000)
    difference = sat - obs
    topocentric = difference.at(t)
    ra, dec, distance = topocentric.radec()
    print(t)
    print(ra._degrees)
    print(dec._degrees)
    print(distance.km)

test()

This script prints the same thing for me under both 1.46 and 1.48; see if it's the same for you too.

You can see the difference directly by asking a time for its .delta_t value, which is measured in seconds. For version 1.46:

from skyfield.api import load

ts = load.timescale()
t = ts.utc(2024, 4, 20, 2, 13)
print(t.delta_t)
69.13558934754973

Whereas version 1.48 prints:

69.197558494097

The difference of about 62 milliseconds is the extra time the Earth needed to rotate to the position that Skyfield had thought it would have reached by 2:13:00 on April 20. (Or is it the other way around, and it's ahead of schedule? I'd need to draw a diagram to figure out which.)

To address this issue, I'll plan to write up a new little section of the "Accuracy and Efficiency" page to explain what sized difference this often makes between Skyfield versions. If you ever want to check what version of the ∆T table your Skyfield is using, you can run this command:

$ python -m skyfield
Skyfield version: 1.48
jplephem version: 2.22
sgp4 version: 2.20
Built-in leap seconds table ends with leap second at: 2016-12-31 23:59:60 UTC
Built-in #T table from finals2000A.all covers: 1973-01-01 to 2025-01-18

@brandon-rhodes brandon-rhodes changed the title Topocentric.radec discrepancy between 1.4.6 and 1.4.8 for satellite-observatory difference Better document that newer Skyfield versions carry more recent ∆T data, which affects observatory position and orientation May 3, 2024
@Dave-II
Copy link
Author

Dave-II commented May 4, 2024

Thanks for the quick response sir. Yes I am getting the same results in 1.46 and 1.48 if using dt=0. Also, as privately suggested, I have used a data 2 years in the past with no dt override, and the the results for 1.46 and 1.48 are very similar, varying by only the 8th or 9th decimal in RA and Dec degrees.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

2 participants