Faster/more accurate ways to get Earth-Sun vectors in different reference frames than transforming get_sun results

Astropy
1

I have a simple function for getting the Earth-Sun position vector in ICRS reference frame at a particular time:

import datetime
import astropy.units as u
import numpy as np
from astropy.coordinates import ICRS, get_sun
from astropy.time import Time

def get_earth_sun_position_vector_icrs_at_epoch(epoch: datetime.datetime) -> np.ndarray:
    earth_sun_gcrs = get_sun(Time(epoch))
    earth_sun_icrs = earth_sun_gcrs.transform_to(ICRS)
    earth_sun_vector_icrs_at_epoch = earth_sun_icrs.cartesian.xyz.to(u.km).value
    return earth_sun_vector_icrs_at_epoch

The conversion from GCRS to ICRS is incredibly expensive, time-wise. 80,000 iterations took over ten minutes.

Is there a faster way to get Earth-Sun vectors at given epochs in different reference frames than this way?

Also - when providing only a Time to get_sun, is the vector from the center of the Earth to the Sun, or from some other point to the Sun?

2

First, there is a fundamental error in your code. earth_sun_gcrs = get_sun(...) is the location of Sun center in the GCRS frame, so earth_sun_icrs = earth_sun_gcrs.transform_to(ICRS) is the location of Sun center in the ICRS frame. That is, the vector you return is merely the vector from the origin of ICRS (the solar-system barycenter) to Sun center, with no connection to the Earth at all.

Second, if what you want is the apparent Earth-Sun vector relative to ICRF axes, including the effects of stellar and planetary aberration, then all you need to do is to skip the transformation to ICRS. The origin of GCRS is Earth center, so earth_sun_gcrs.cartesian is the Earth-Sun vector relative to ICRF axes, including the above observer effects.

Alternatively, if what you want is instead the true/geometric Earth-Sun vector relative to ICRF axes, you shouldn’t use get_sun() at all. Instead, you can call get_body_barycentric() to get the true location of a solar-system body in ICRS, e.g.:

from astropy.coordinates import get_body_barycentric
epoch = Time(epoch)
earth_icrs = get_body_barycentric('earth', epoch)
sun_icrs = get_body_barycentric('sun', epoch)
earth_sun_vector_icrs_at_epoch = (sun_icrs - earth_icrs).xyz.to_value(u.km)

Hope that helps!

1 Like
3

It does make sense, thank you.

It does bring up related questions for my project though.

You say ICRF axes use the solar system barycenter - it was my understanding that ICRS is closely analogous to MJ2000 axes/origin. Is that incorrect?

And how might I do the same sort of thing, but getting the TEME earth-sun vector?

4

You can use ICRF axes with any origin. Setting aside some details:

  • ICRS = ICRF axes with the origin at the solar-system barycenter
  • GCRS = ICRF axes with the origin at Earth center
  • HCRS = ICRF axes with the origin at Sun center
  • etc.

My understanding of MJ2000 is that it’d be very similar to GCRS.

If you have the Sun’s location as a SkyCoord, you can transform to Astropy’s TEME, which has the origin at Earth center, so the vector will be the Earth-Sun. I’ll caution that there are multiple definitions of TEME, so it’d be prudent for you to confirm whether Astropy’s definition of TEME is what you need.

2 Likes
5

This has helped clear up a lot.

The code you have for getting the true/geometric Earth-Sun vector relative to ICRF axes - how might I modify it to use GCRS, to be Earth-centered? Like this?:

from astropy.coordinates import get_sun
epoch = Time(epoch)
earth_sun_gcrs = get_sun(epoch)
earth_sun_vector_gcrs = earth_sun_gcrs.cartesian
6

Again, if what you want is the true/geometric Earth-Sun vector relative to ICRF axes, you shouldn’t use get_sun() at all. Note that as the difference of two positions, this vector is the same irrespective of the notion of the coordinate-frame origin, e.g., SSB for ICRS or Earth center for GCRS.

However, if what you want is the Earth-Sun vector relative to ICRF axes, but including the effects of stellar and planetary aberration for an observer at Earth center, then yes, your new code block gives you that.