In [1]:
import ephem
import numpy as np
import datetime
import matplotlib.pyplot as plt
%matplotlib inline
plt.style.use('ggplot')
In [5]:
iss = "ISS (ZARYA)", "1 25544U 98067A 16145.69686492 .00003239 00000-0 55059-4 0 9992", "2 25544 51.6424 175.2139 0001726 113.6054 353.3802 15.54685316 1333"
start = datetime.datetime.now()
end = start + datetime.timedelta(30) #end time: +30 days
step = datetime.timedelta(0,5) #time step in s
sat = ephem.readtle(*iss)
pdx = ephem.Observer()
pdx.lat = '45.5231'
pdx.lon = '122.6765'
pdx.elevation = 30
critical_az_angle = ephem.degrees(30)
In [6]:
pass_happening = False
now = start
pass_list = []
elevations = []
while now < end:
pdx.date = now
sat.compute(pdx)
alt = np.degrees(sat.alt)
if pass_happening:
if alt < critical_az_angle:
pass_happening = False
pass_end = now
pass_list.append(((pass_end-pass_start).total_seconds(), pass_start, pass_end))
else:
if alt > critical_az_angle:
pass_happening = True
pass_start = now
elevations.append(alt)
now += step
In [7]:
pass_lengths, pass_starts, pass_ends = zip(*pass_list)
pass_lengths = [x/60.0 for x in pass_lengths]
plt.figure(figsize=(15,10))
n, bins, patches = plt.hist(pass_lengths, 10)
plt.xlabel("Length of pass (minutes)")
plt.ylabel("Number of passes")
plt.title(u"ISS pass lengths over PDX for %g days above critical elevation angle = %g°"%((end-start).days, critical_az_angle))
plt.show()
In [ ]: