python4oceanographers

Turning ripples into waves

Playing with pyephem

I'm teaching tides in a Waves and Tides course and, every now and then, I find myself creating new figures to show the students. This post is to show how I created a "moon phase path" using pyephem. The final plot was based on this example.

Just for fun I'll find the city we are in using the machine IP instead of specifying the longitude and latitude. The function below returns the IP.

In [2]:
import re
import urllib

def get_ip(url='http://checkip.dyndns.org'):
    request = urllib.urlopen(url).read().decode('utf-8')
    return re.findall(r"\d{1,3}\.\d{1,3}\.\d{1,3}.\d{1,3}", request)[0]

Now, to find the city, we need the GeoLiteCity database and the GeoIP module.

In [3]:
from pygeoip import GeoIP

def get_location(ip, fname='./data/GeoLiteCity.dat'):
    """Database can be downloaded at:
        http://dev.maxmind.com/geoip/legacy/geolite/"""
    gi = GeoIP(fname)
    location = gi.record_by_addr(ip)
    return location
In [4]:
ip = get_ip()
location = get_location(ip)
location
Out[4]:
{'area_code': 0,
 'city': u'Cama\xe7ari',
 'continent': 'SA',
 'country_code': 'BR',
 'country_code3': 'BRA',
 'country_name': 'Brazil',
 'dma_code': 0,
 'latitude': -12.483300000000014,
 'longitude': -38.2167,
 'metro_code': None,
 'postal_code': None,
 'region_code': u'05',
 'time_zone': 'America/Bahia'}

PyEphem main object is the Observer class. We need to "feed it" with our position and date. Also, since we are interested in Earth's natural satellite, we will start a moon object as well.

In [5]:
import ephem
import numpy as np
from datetime import datetime

date = datetime.utcnow()
moon = ephem.Moon()

obs = ephem.Observer()
obs.date = date
obs.lon = np.deg2rad(location['longitude'])
obs.lat = np.deg2rad(location['latitude'])

We can print some information like rising times, and next full moon date.

In [6]:
print("Next full Moon: %s" % ephem.next_full_moon(date))
print("Previous rising: %s" % obs.previous_rising(moon))
print("Next rising: %s" % obs.next_rising(moon))
Next full Moon: 2016/2/22 18:19:52
Previous rising: 2016/2/2 03:02:17
Next rising: 2016/2/3 03:47:43

Now let's loop over one day saving the azimuth, altitude, hour, and a symbol for the moon phase.

In [7]:
az, alt, symbols, times = [], [], [], []

for k in range(24):
    moon.compute(obs)  
  
    nnm = ephem.next_new_moon(obs.date)  
    pnm = ephem.previous_new_moon(obs.date)  
    lunation = (obs.date - pnm) / (nnm - pnm)  
    symbol = lunation * 26
    if symbol < 0.2 or symbol > 25.8:
        symbol = '1'
    else:  
        symbol = chr(ord('A') + int(symbol + 0.5) - 1)
        times.append(ephem.localtime(obs.date).time().strftime("%H:%M"))
        symbols.append(symbol)
        alt.append(moon.alt)
        az.append(moon.az)
    
    obs.date += ephem.hour
    
az, alt = map(np.rad2deg, (az, alt))

Finally, we can plot the figure with the moon path.

In [8]:
import matplotlib.pyplot as plt
from matplotlib.font_manager import FontProperties

prop = FontProperties(fname='./data/moon_phases.ttf', size=25)

fig, ax = plt.subplots()
offset = 20
ax.set_xlim(0-10, 360+30)
ax.set_ylim(-90-5, 90+5)
[ax.text(x, y, text, fontproperties=prop, alpha=0.3, zorder=1)
 for x, y, text in zip(az, alt, symbols)]
[ax.text(x, y, hour, fontdict=dict(size=11, weight='bold'), zorder=0, color='r')
 for x, y, hour in zip(az, alt, times)]
                            
ax.set_xlabel("Azimuth (compass direction, in degrees)")
ax.set_ylabel("Elevation (degrees above horizon)")
Out[8]:
<matplotlib.text.Text at 0x7f772abae910>

The font moon_font.ttf convert the symbol we crated into a nice moon phase icon.

In [9]:
HTML(html)
Out[9]:

This post was written as an IPython notebook. It is available for download or as a static html.

Creative Commons License
python4oceanographers by Filipe Fernandes is licensed under a Creative Commons Attribution-ShareAlike 4.0 International License.
Based on a work at https://ocefpaf.github.io/.

Comments