skyfielders / python-skyfield

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

Issue with my Programm on propagating Satellite Orbits, plotting them and calculating minimal distance #848

Closed callmeArlette closed 1 year ago

callmeArlette commented 1 year ago

Hi, I have got a problem. My goal is to programm a code, which reads in two TLE files and then recreates an Event with two satellites involved. I have got a Code but can't figure out, how i can plot the timeframe which i entered and then calculate the minimal distance. I somehow always end up in plotting a straight line instead of an orbit.... Also i want to mark the points from the data in the orbit to see in which direction my satellite moves... Does somebody know how i can do that?

import skyfield
from skyfield.api import load, EarthSatellite, wgs84
#from AltAzRange TODO install
from skyfield.api import Topos
import numpy as np
from skyfield.timelib import Time
from skyfield.api import utc
import skyfield.positionlib
from datetime import datetime
import csv
import os
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from mpl_toolkits.mplot3d import proj3d
from matplotlib.text import Annotation

#read in TLE and TXT files (for event Astrocast x Fossasat)
sat1_file =  r'Astrocast_0103_TLE_history.txt'
sat2_file =  r'FOSSASAT-2E6_06-2022.tle'
name = ['Astrocast_0103','FOSSASAT-2E6_06-2022']

ts = load.timescale()

# Load TLE data and create EarthSatellite objects
def load_sat_files1(sat_file):
    satellite_file = []
    with open(sat1_file, 'r') as f:
        lines = f.readlines()
        for i in range(0, len(lines), 2):
            line1 = lines[i]
            line2 = lines[i+1]
            satellite = EarthSatellite(line1, line2, name[0], ts)
            satellite_file.append(satellite)
    return satellite_file

def load_sat_files2(sat_file):
    satellite_file = []
    with open(sat1_file, 'r') as f:
        lines = f.readlines()
        for i in range(0, len(lines), 2):
            line1 = lines[i]
            line2 = lines[i+1]
            satellite = EarthSatellite(line1, line2, name[1], ts)
            satellite_file.append(satellite)
    return satellite_file

satellite1 = load_sat_files1(sat1_file)

satellite2 = load_sat_files2(sat2_file)

# Zeitintervall
start_time = datetime(2022, 6, 14, 0, 0, 0)
start_time = start_time.replace(tzinfo=utc)
end_time = datetime(2022, 6, 20, 0, 0, 0)
end_time = end_time.replace(tzinfo= utc)

#Time instances 
start_t = ts.from_datetime(start_time)
end_t = ts.from_datetime(end_time)

#filter of the times needed (shortly before and after event)
def filter_times(time_start, time_end, satellit):
    satellite_filtered = []
    for sat in satellit:
        epoch_t = sat.epoch.utc_datetime()
        #print(epoch_t)
        if start_t.utc_datetime() <= epoch_t <= end_t.utc_datetime():
            satellite_filtered.append(sat)
    return satellite_filtered

sat1_filtered = filter_times(start_t, end_t, satellite1)
sat2_filtered = filter_times(start_t, end_t, satellite2)

#get positions at times from above from both sats
def positions_sat(sat_filter):
    satellite_positions = []
    for i in range(len(sat_filter)):
        satellite_pos = sat_filter[i].at(sat_filter[i].epoch).position.km
        satellite_positions.append(satellite_pos)
    return satellite_positions

time_step = np.arange(0, 24*60*60, 10)

#print(time_step)
sat1_pos_ep_0 = sat1_filtered[0].at(ts.utc(2022, 6, 14, time_step)).position.km
#print(satellite1_positions)
sat2_pos_ep_0 = sat2_filtered[0].at(ts.utc(2022, 6, 14, time_step)).position.km

#plotting
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot(sat1_pos_ep_0[0,:], sat1_pos_ep_0[1,:], sat1_pos_ep_0[2,:], label=name[0])
ax.plot(sat2_pos_ep_0[0,:], sat2_pos_ep_0[1,:], sat2_pos_ep_0[2,:], label=name[1])
ax.set_xlabel('X (km)')
ax.set_ylabel('Y (km)')
ax.set_zlabel('Z (km)')

plt.legend()
plt.show()

#TODO Calculate minimal distance 
brandon-rhodes commented 1 year ago
  1. Alas, I don't think I'm going to have time this month to try running and debugging third-party scripts.
  2. But, I have at least edited your question to make it easier for other contributors to help you if they are interested: I have added triple-back-ticks around your code, so that it forms successfully and so that folks can cut and paste it.
  3. Your code won't run for people, if they want to try it out, unless you provide the data files too, so you might think about providing them — or else maybe sharing a screenshot of the output plot, so people can see the result without running it?
  4. Since this Issue doesn't report a bug with Skyfield, but rather with your own script, I'm going to change it to a Discussion so that folks find it in the correct category.

We will see if anyone responds. I hope you're able to get something working!