Solar system and lunar ephemerides via jplephem take 2 by StuartLittlefair · Pull Request #4890 · astropy/astropy

StuartLittlefair · 2016-05-13T10:44:53Z

This PR supersedes #4746. @bmorris3 's description of that PR is included here, slightly modified.

This PR implements an interface for grabbing solar system ephemerides from jplehem (MIT licensed), after caching a copy of the DE430 table (see the jplephem docs for info). There are three new public functions in this PR. The two expected to see most use are get_moon(time, [location]) and get_body(time, body_key, [location]). Both return SkyCoord objects with GCRS frames, the former for Earth's moon and the latter for the Sun and solar system "planets" Mercury through Pluto inclusive. The final new functions is get_barycentric_body_position(time, body_key) which returns a CartesianRepresentation of the ICRS coordinate of the solar system body in question.

Here are the new functions in action:

In [1]: from astropy.coordinates import (get_moon, get_body, EarthLocation,
                                                              get_barycentric_body_position)
In [2]: from astropy.time import Time
In [6]: loc = EarthLocation.of_site('Apache Point')
In [7]: moon = get_moon(Time.now(), loc)
In [8]: moon
Out[8]:
<SkyCoord (GCRS: obstime=2016-03-29 06:35:36.927857, obsgeoloc=[-1463969.30185172 -5166673.34223433  3434985.71204565] m, obsgeovel=[ 0.  0.  0.] m / s): (ra, dec, distance) in (deg, deg, km)
    (250.94788165, -17.04585998, 400244.30166804)>

In [9]: jupiter = get_body(Time.now(), 'jupiter', loc)
In [10]: jupiter
Out[10]:
<SkyCoord (GCRS: obstime=2016-03-29 06:36:28.310851, obsgeoloc=[-1463969.30185172 -5166673.34223433  3434985.71204565] m, obsgeovel=[ 0.  0.  0.] m / s): (ra, dec, distance) in (deg, deg, km)
    (167.22360841, 7.07220939, 673151288.66240811)>
In [11]: icrs_moon = get_barycentric_body_position(Time.now(), 'moon')
In [12]: icrs_moon
Out[12]: 
<CartesianRepresentation (x, y, z) in km
    (-90980973.95347136, -110026376.34185831, -47742071.3962695)>

I am comparing the outputs to JPL Horizons and Skyfield. Comparisons against JPL Horizons are awkward to get right because the ra and dec returned are not in the GCRS coordinate frame but are w.r.t to the Earth's true equator and poles. You can get the same results from Skyfield py passing epoch='date' to the radec function. I made a brief wrapper function to convert to ra and dec w.r.t to Earth's true equator and poles for comparison.

For geocentric observers I get compatible answers (to <1 arcsec) to both Horizons and Skyfield. The situation is worse when we have a non-geocentric observer. Errors can be as large as 1 degree for the Moon, but are much smaller for planets. I believe this is a bug in the astropy transformation chain for GCRS frames. I am shortly about to submit another pull request to fix that.

(This PR addresses #3219, #2244, #4269 and overlaps in functionality with PRs #4724, #4736.)

…ess jplephem installed

…ote_data tags

…CRS coordinates uses full astropy machinery. This also includes a change to planet index numbering.

axes align instead with the true equator and poles of the Earth. This should eventually find a home in a new frame in the coordinates package

…x ra and dec before comparison

PEP8 fixes

astrofrog · 2016-05-13T11:19:57Z

+                        unicode_literals)
+
+import numpy as np
+import six


This is causing failures - should be from ...extern import six

astrofrog · 2016-05-13T11:28:33Z

@StuartLittlefair - thanks for doing this! I have a few comments above, and my only other main comment is that this definitely needs to be mentioned in the main documentation (I think we need a new section in http://docs.astropy.org/en/stable/coordinates/index.html#using-astropy-coordinates, even though it doesn't need to be long). Is this something you can do today to make sure we include this in 1.2?

StuartLittlefair · 2016-05-13T11:47:41Z

I can make these changes in the next few hours. The tests will fail because of the outstanding issue with ICRS->GCRS transformations for non-geocentric observers. Should I remove non-geocentric tests for now?

astrofrog · 2016-05-13T12:47:24Z

@StuartLittlefair - you mentioned you could make a pull request for the ICRS -> GCRS issue, so if that's something you can also do today, no need to skip the tests, as once the fix is merged we can just rebase this. Or do you think the issue will be more intricate and harder to fix?

StuartLittlefair · 2016-05-13T21:29:14Z

As @eteq has pointed out there, PR #4894 was a mistake on my part and it not needed. What is needed as part of this PR is a way of getting the GCRS position and velocity of an EarthLocation at a given time. I have included code for this here.

@eteq - thanks for the review. On your more substantive points:

Yes. I too was worried about the light deflection part of the ICRS->GCRS transformation. However, it's correct to have those light deflections in I think (suppose Jupiter is along the L.O.S to the body you request). As you say the tests seem to work.
Since it's getting late for me and I've made way too many mistakes along the way already I think I will defer allowing the user to specify the kernel until a later date. Although the way you suggest doing it is probably the way forward.

eteq · 2016-05-13T21:50:40Z

+
+    >>> from astropy.coordinates import SOLAR_SYSTEM_BODIES
+    >>> SOLAR_SYSTEM_BODIES
+    ['sun', 'mercury', 'venus', 'earth-moon-barycenter', 'earth', 'moon', 'mars', 'jupiter', 'saturn', 'uranus', 'neptune', 'pluto']


these now need to be (/), instead of [/]

eteq · 2016-05-13T21:51:29Z

Ok, that all sounds good, @StuartLittlefair. I noticed one tiny thing you'll need to fix to get the tests to pass but I think this is good assuming that all works. I might manage to squeeze out a PR to do the 2) I suggested myself depending on how my other last few things go.

eteq · 2016-05-13T22:53:10Z

Oh, but if it wasn't obvious, I'm 👍 on this now once we see the tests pass. Great teamwork, @bmorris3 and @StuartLittlefair ! This is what open source development is all about!

eteq · 2016-05-13T22:54:39Z

+
+def get_body_barycentric(time, body):
+    """
+    Calculate the barycentric position of thesolar system body ``body``


oops, just noticed a missing space here - feel free to fix it in a [ci skip] commit. Or if the tests pass before you do so, can always just fix small things like this as part of the merge.

astrofrog · 2016-05-13T23:58:14Z

(just FYI, I cancelled the Travis build on any commits that weren't the last one, so it should run the last commit soon)

mhvk · 2016-05-14T00:19:24Z

                                     for the location of this object at the
                                     default ``obstime``.""")

+    def get_gcrs_posvel(self, obstime):


Tiny, but can this be just gcrs_posvel -- I think the standard is to assume a calculation gets done (and your not "getting" the information from somewhere).

I just realised there is a get_itrs but there the obstime is optional and we also have the itrs property. Here, a time is obviously needed, so the situation is not quite the same. Still, given that precedent, probably fine just to keep as is.

mhvk · 2016-05-14T00:29:45Z

Yes more comments, though tiny. Overall, this is a wonderful addition!!

Solar system and lunar ephemerides via jplephem take 2

eteq · 2016-05-14T04:39:42Z

I went ahead and implemented @mhvk's final comments as a final commit in the merge process. So this is now in!

bmorris3 · 2016-05-14T19:31:47Z

Thanks to @StuartLittlefair and all the others (@astrofrog @mhvk @eteq ) who got this together in my absence! I'm really glad this will be available now.

Now it's time to revisit astropy/astroplan#111, astropy/astroplan#35.

bmorris3 and others added 15 commits May 13, 2016 09:43

Implementing solar system ephemerides with jplephem dependency

4a59bb8

Cleaning up docstrings, adding changelog entry, turning off tests unl…

aa7e7a3

…ess jplephem installed

Adding jplephem to optional dependencies in travis matrix, adding rem…

51dc414

…ote_data tags

Moving jplephem to pip install

5e1c476

Correcting dependencies with advice from @bsipocz

301b942

Refactor code so that transformation from JPL (ICRS) coordinates to G…

067e85a

…CRS coordinates uses full astropy machinery. This also includes a change to planet index numbering.

PEP8 consistency

d123fc2

changes to allow indexing by body name, e.g 'moon'

cd408c2

changes so _get_earth_body_vector works with non-scalar times

2099197

A helper routine to rotate a GCRS frame into one in which the

0dedc66

axes align instead with the true equator and poles of the Earth. This should eventually find a home in a new frame in the coordinates package

tests w.r.t skyfield added and coordinates converted into true equino…

cd6a6d7

…x ra and dec before comparison

PEP8 fixes and leading underscore on apparent_position routine

7d0701b

PEP8 fixes

updated name of routines

de7ac8c

use different times for body and Earth in light travel time calculations

a7fd6e1

publically expose get_barycentric_body_position function

48f6765

StuartLittlefair mentioned this pull request May 13, 2016

Solar system and lunar ephemerides via jplephem #4746

Closed

astrofrog reviewed May 13, 2016
View reviewed changes

pllim added coordinates Enhancement labels May 13, 2016

Changes requested by @astrofrog

10e15f5

StuartLittlefair added 2 commits May 13, 2016 22:17

typo in docstring

32bb4bb

changes to function name

6e5f48a

eteq reviewed May 13, 2016
View reviewed changes

minor changes to tests

3f8876b

eteq assigned eteq and unassigned mhvk May 13, 2016

eteq reviewed May 13, 2016
View reviewed changes

StuartLittlefair added 2 commits May 14, 2016 00:03

typo [ci skip]

79ec245

document and PEP8 changes

1e23780

mhvk reviewed May 14, 2016
View reviewed changes

eteq merged commit 1e23780 into astropy:master May 14, 2016

eteq added a commit that referenced this pull request May 14, 2016

Merge pull request #4890 from StuartLittlefair/solarsystem

9d569a6

Solar system and lunar ephemerides via jplephem take 2

eteq mentioned this pull request May 14, 2016

Solar system epherides doc update and kernel URL setting option #4908

Merged

This was referenced May 14, 2016

Wishlist: support JPL ephemerides #3219

Closed

Lunar calculations not yet supported astropy/astroplan#35

Closed

Add get_moon function to astropy.coordinates using the erfa.plan94 #4736

Closed

Feat/get moon #4724

Closed

mhvk mentioned this pull request May 14, 2016

Approximate planet positions using erfa.plan94 #4911

Closed

StuartLittlefair mentioned this pull request May 16, 2016

Cannot create GCRS frames with obsgeoloc shapes other than (3,) #4926

Closed

mhvk mentioned this pull request May 18, 2016

Use solar_system_ephemeris for barycentric corrections and coordinate transforms? #4942

Closed

Conversation

StuartLittlefair commented May 13, 2016

Uh oh!

astrofrog May 13, 2016 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

astrofrog commented May 13, 2016

Uh oh!

StuartLittlefair commented May 13, 2016

Uh oh!

astrofrog commented May 13, 2016 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

StuartLittlefair commented May 13, 2016 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

eteq May 13, 2016

Choose a reason for hiding this comment

Uh oh!

eteq commented May 13, 2016

Uh oh!

eteq commented May 13, 2016

Uh oh!

eteq May 13, 2016

Choose a reason for hiding this comment

Uh oh!

astrofrog commented May 13, 2016

Uh oh!

mhvk May 14, 2016

Choose a reason for hiding this comment

Uh oh!

mhvk May 14, 2016

Choose a reason for hiding this comment

Uh oh!

mhvk commented May 14, 2016

Uh oh!

eteq commented May 14, 2016

Uh oh!

bmorris3 commented May 14, 2016

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

6 participants

astrofrog May 13, 2016 •

edited

Loading

astrofrog commented May 13, 2016 •

edited

Loading

StuartLittlefair commented May 13, 2016 •

edited

Loading