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 [ ]: