Sun and moon's orbit in ecliptic coordinates using skyfield

Question

I am new to using skyfield, is there any doc or help file that can show me on how to get the orbit of Sun and moon in ecliptic coordinate for a particular date and time. This is a follow up question of this question

I don't really see how this question is different from the one you link to. — usernumber, Nov 28 '19 at 10:21
@usernumber that was how can be done, this is for how can be done using skyfield — Rasika, Nov 28 '19 at 10:26

score 3 · Accepted Answer · edited Nov 25 '22 at 10:21

3

I cant help you with skyfield, but I usually use JPL Horizons Web interface. No installation required, you can also print it in a text file if you want:

https://ssd.jpl.nasa.gov/horizons.cgi

Otherwise I found the documentation for skyfield: https://rhodesmill.org/skyfield/toc.html

And if nothing of that works, I made a little astropy script for you

from astropy import units as u
from astropy.coordinates import SkyCoord, EarthLocation, AltAz, get_body
from astropy.time import Time
import numpy as np
Create 1000 Timepoints between Time 1 and Time 2 (one year later)
t = np.linspace(2451545, 2451545+365, 1000)
pointlist = []
#Loop through this times
for tn in t:
    # For every timepoint, create an astropy_time object
    astropy_time = Time(tn, format="jd")
    # Get Planet (as string, "earth", "moon", "mercury" etc. in aequatorial coordinates
    planet_aequatorial = get_body("moon", time = astropy_time)
    #Transform to Barycentric True Ecliptic (relative to the center of mass of the solar system).
    planet_ecliptic =  planet_aequatorial.transform_to("barycentrictrueecliptic")
# Add a point to the orbit. Every point is described as (longitude [deg], latitude (ecliptic coords), distance (km))
pointlist.append([planet_ecliptic.lon.deg, planet_ecliptic.lat.deg, planet_ecliptic.distance.km])
print(planet_ecliptic.distance.km)
# So pointslist is a 2D array. The rows are all the 1000 points of the orbit


In every point there is 3 columns for [Long, Lat, Distance]
print(pointlist)
You can also save the result with
pointlist = np.array(pointlist)

np.save("results.npy", pointlist)

edited Nov 25 '22 at 10:21

Prince Pugs

596
4
23

answered Nov 29 '19 at 17:32

Andres Huster

164
5

thanks or the reply, it gives me a single point, however i want complete orbit, is that possible, my question may look silly, apologise – Rasika Nov 29 '19 at 17:44
create a list of times (in julian days) with

import numpy as np t = np.linspace(2451545, 2451545+365, 2000)

then loop over the above code but create Time object with

t = Time(t1, format="jd")

instead
– Andres Huster Dec 09 '19 at 22:47
did that help or should i do it for you? – Andres Huster Jan 03 '20 at 13:01
yes, please. that will be really helpful – Rasika Jan 18 '20 at 19:01
can you give me your e-mail or something? i cant post all the code here i think... – Andres Huster Jan 23 '20 at 11:11
I edited the first answer up here – Andres Huster Jan 23 '20 at 18:25
Sorry for the late reply, i used it, but the output latitude does not go beyond 2 degree – Rasika Mar 22 '20 at 12:16

Sun and moon's orbit in ecliptic coordinates using skyfield

1 Answers1

Create 1000 Timepoints between Time 1 and Time 2 (one year later)

In every point there is 3 columns for [Long, Lat, Distance]

You can also save the result with