¿Mejor manera de obtener la pista terrestre aproximada para un satélite usando Skyfield?

Quería hacer un gráfico rápido de seguimiento en tierra, así que usé el paquete Skyfield de Python para propagar un TLE y devolver las coordenadas inerciales centradas en la Tierra. Luego creé un Topos()objeto fijado a la superficie de la Tierra en lat/lon = 0, 0 y lo usé para generar la rotación de la tierra para "desenrollar" la posición de la ISS en la superficie.

nota al margen: esta no es una buena manera de hacer esto, pero da resultados lo suficientemente buenos como para hacer un pequeño mapa para una ilustración. Los problemas incluyen asumir que el eje de la tierra todavía está en la dirección z y, por supuesto, tratar la superficie de la tierra como esférica.

¿Hay una mejor manera de hacer esto dentro de Skyfield sin usar un método que comience con un guión bajo, en otras palabras, usando métodos destinados a ser completamente públicos para el usuario? Además, ¿hay alguna forma de centrar la Tierra, fijar las coordenadas de la Tierra (es decir, rotar la Tierra) directamente en Skyfield sin desenrollarse así?

EDITAR: he ajustado el guión desde que ha pasado más de un año y se lanzó Skyfield v 1.0.

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"""
L1, L2 = ISS_TLE.splitlines()

import numpy as np
import matplotlib.pyplot as plt
from skyfield.api import Loader, EarthSatellite, Topos

degs     = 180./np.pi

r_earth  = 6371.  # for ROUGH approx. ground track, just use a spherical Earth

load    = Loader('~/Documents/YourNameHere/SkyData')
data    = load('de421.bsp')
earth   = data['earth']
topoZZ  = Topos(latitude_degrees=0.0, longitude_degrees=0.0)

location = earth + topoZZ

ISS      = earth + EarthSatellite(L1, L2)

ts       = load.timescale()

minutes  = np.arange(0, 140, 0.5)
time     = ts.utc(2016, 12, 7, 12, minutes)

Epos     = earth.at(time).position.km
ZZpos    = topoZZ.at(time).position.km    ## Position of (0.0N, 0.0E) to get rotation
ISSpos   = ISS.at(time).position.km - Epos

theta_ZZ = np.arctan2(ZZpos[1], ZZpos[0])   # calculate Earth's rotaion

sth, cth         = np.sin(-theta_ZZ), np.cos(-theta_ZZ) # unwind
xISS, yISS, zISS = ISSpos
xISSnew, yISSnew = xISS*cth - yISS*sth, xISS*sth + yISS*cth # rotate ISS data to match Earth
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()

ingrese la descripción de la imagen aquí

Esto también se ha preguntado aquí: github.com/skyfielders/python-skyfield/issues/121

Respuestas (1)

Sé que esta es una vieja pregunta, pero para los divertidos, aquí hay un guión rápido.

Esta pregunta se hizo cuando la .subpoint()característica no era compatible con Skyfield para captar la longitud/latitud de una órbita de satélite proyectada sobre la Tierra.

Aquí hay una secuencia de comandos rápida que muestra cómo usar las funciones integradas de Skyfield con el trazado usando cartopy para construir el mapa y las proyecciones.

from skyfield.api import load, EarthSatellite
import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs

ts = load.timescale(builtin=True)


TLE = """ISS (ZARYA)             
1 25544U 98067A   19203.81086311  .00000606  00000-0  18099-4 0  9996
2 25544  51.6423 184.5274 0006740 168.1171 264.4057 15.50995519180787"""

name, L1, L2 = TLE.splitlines()

sat = EarthSatellite(L1, L2)

minutes = np.arange(0, 200, 0.1) # about two orbits
times   = ts.utc(2019, 7, 23, 0, minutes)

geocentric = sat.at(times)
subsat = geocentric.subpoint()

fig = plt.figure(figsize=(20, 10))
ax = fig.add_subplot(1, 1, 1, projection=ccrs.PlateCarree())

ax.stock_img()

plt.scatter(subsat.longitude.degrees, subsat.latitude.degrees, transform=ccrs.PlateCarree(),
            color='red')
plt.show()
¡Excelente! Para diversión adicional, ¿qué tal una forma de construir un subpunto terrestre para otro objeto? Los comentarios debajo de esta pregunta también en este chat sugieren que al OP eventualmente le gustaría trazar un rastro de los puntos del subplaneta a través de la superficie de la Tierra. Puedo pensar en dos maneras; 1) usar reverse_terra()en esta respuesta (cont.)
...o 2) construya un objeto geocéntrico en Skyfield basado en las posiciones del planeta y luego aplíquelo .subpoint(). Lamentablemente, el texto de esa pregunta aún no se ha actualizado para reflejar la intención del OP. (nota: la pregunta vinculada está en Astronomy SE, no aquí en Space SE], es común confundir los dos)
Para tu información, ahora que tienes más de 50 puntos de reputación, ¡puedes dejar comentarios en cualquier publicación que te guste!
Descubrí que la línea "ax.stock_img() da como resultado el siguiente error de Python: FileNotFoundError: [Errno 2] No such file or directory: '/usr/lib/python3/dist-packages/cartopy/data/raster/ natural_earth/50-natural-earth-1-downsampled.png': parece que la estructura de carpetas ha cambiado y la imagen se ha movido. ¿Dónde más está disponible el mismo archivo?
¿Mejor manera de obtener la pista terrestre aproximada para un satélite usando Skyfield?
RespuestasAquíPolítica de privacidad1y, 0v