skyfielders / python-skyfield

Elegant astronomy for Python
MIT License
1.41k stars 211 forks source link

How to get a simple (approx) satellite groundtrack without touching (semi)private method? #121

Closed davidmikolas closed 6 years ago

davidmikolas commented 7 years ago

I did the following to make a quick ground-track plot of the ISS.

  1. I used this, ISSpos = ISS._position_and_velocity_TEME_km(time)[0] # semi-private a method with a single underscore which seems to be a very handy method. Is this OK or is there a non-underscored way to get the same information in the Earth-inertial or earth-fixed coordinates?

1.1 Neither earth.geometry_of(satellite) nor satellite.geometry_of(earth) is working for me.

  1. Is there a better way to get ground track? (spherical earth is fine for me now)

note: I've asked here also: http://space.stackexchange.com/q/19339/12102

ISS_TLE = """1 25544U 98067A   16341.96974289  .00003303  00000-0  57769-4 0  9996
2 25544  51.6456 276.4739 0005937 300.1004 104.8148 15.53811586 31866"""

import numpy as np
import matplotlib.pyplot as plt
from skyfield.api import load

pi, twopi, halfpi  = np.pi, 2.*np.pi, 0.5*np.pi
degs, rads         = 180./pi, pi/180.

r_earth = 6371.  # approximately, though it's a bit oblate in reality

data     = load('de421.bsp')
Earth    = data['earth']
ZeroZero = Earth.topos(0.0, 0.0)
ISS      = Earth.satellite(ISS_TLE)

ts       = load.timescale()
minutes  = np.arange(0, 12*60, 1)
time     = ts.utc(2016, 12, 7, 12, minutes)

Epos      = Earth.at(time).position.km
ZZpos     = ZeroZero.at(time).position.km - Epos   ## Position of (0N, 0E) to get rotation
ISSpos    = ISS._position_and_velocity_TEME_km(time)[0]   # semi-private

theta_ZZ  = np.arctan2(ZZpos[1], ZZpos[0])

sth, cth         = np.sin(-theta_ZZ), np.cos(-theta_ZZ) # calculate Earth's rotaion
xISS, yISS, zISS = ISSpos
xISSnew, yISSnew = xISS*cth - yISS*sth, xISS*sth + yISS*cth # rotate ISS to match Earthn
ISSnew           = np.vstack((xISSnew, yISSnew, zISS))

x, y, z = ISSnew
r       = np.sqrt((ISSpos**2).sum(axis=0))
rxy     = np.sqrt(x**2 + y**2)
ISSlat, ISSlon = np.arctan2(z, rxy), np.arctan2(y, x)

plt.figure()
plt.plot(degs*ISSlon, degs*ISSlat, 'ok')
plt.show()
yeus commented 7 years ago

Same problem here. I would like to get satellite coordinates (lat,lon) from skyfield.

I used this, ISSpos = ISS._position_and_velocity_TEME_km(time)[0] # semi-private a method with a single underscore which seems to be a very handy method. Is this OK or is there a non-underscored way to get the same information in the Earth-inertial or earth-fixed coordinates?

--> wouldn't the correct way to get the position be (instead of _position_and_velocity_TEME_km):

 ISSpos = ISS.gcrs(self.time).position.km

as ZeroZero.at(time).position.km is also in ITRF?

And also: is there a better way to get the position from a certain point on earths surface than doing:

 ZZpos     = ZeroZero.at(time).position.km - Epos

? something like: ZeroZero.at(time).gcrf.position.km ...

thank you for all the effort!

digital-idiot commented 6 years ago

Stuck with the same limitation. Using PyEphem for now. Is there any update regarding progress or solution?

astrojuanlu commented 6 years ago

You might be interested in trying out orbit-predictor:

http://nbviewer.jupyter.org/github/Juanlu001/groundtrack-plotting/blob/master/Groundtrack.ipynb

digital-idiot commented 6 years ago

Did not know about this. I will give it a try. Thanks

brandon-rhodes commented 6 years ago

There is now a supported method for getting ground track latitudes and longitudes. I've added it to the documentation here:

http://rhodesmill.org/skyfield/api-position.html#skyfield.positionlib.Geocentric.subpoint

If you have any problems with it, feel free to open a new issue.