Notebook

In [1]:

from IPython.core.display import HTML

with open('creative_commons.txt', 'r') as f:
    html = f.read()
    
name = '2014-08-11-supermoon'

html = '''
<small>
<p> This post was written as an IPython notebook.
 It is available for <a href='https://ocefpaf.github.com/python4oceanographers/downloads/notebooks/%s.ipynb'>download</a>
 or as a static <a href='https://nbviewer.ipython.org/url/ocefpaf.github.com/python4oceanographers/downloads/notebooks/%s.ipynb'>html</a>.</p>
<p></p>
%s''' % (name, name, html)

%matplotlib inline
from matplotlib import style
style.use('ggplot')

Yesterday we had the super "moon event". Unfortunately Salvador was cloudy and the show was not as impressive as we expected. All the same, watching the sunset while the moon was rising by the Barra lighthouse was quite impressive. The super moon event inspired me to go back there next full moon to take a fell pictures. Sure I could just Google for the date and time. But let's do this the hard way.

The next two functions are part of a previous post, so I will not get into the details. Just to remember: get_ip returns the local machine ip address and get_location gets, among other things, the geographical position.

In [2]:

import re
import urllib
from pygeoip import GeoIP

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

def get_location(ip, fname):
    """Database can be downloaded at:
        http://dev.maxmind.com/geoip/legacy/geolite/"""
    gi = GeoIP(fname)
    location = gi.record_by_addr(ip)
    return location

We get a nice dictionary back with all the information we need.

In [3]:

ip = get_ip()
fname = './data/GeoLiteCity.dat'
location = get_location(ip, fname)
location

Out[3]:

{'country_name': 'Brazil',
 'continent': 'SA',
 'longitude': -38.516699999999986,
 'country_code3': 'BRA',
 'country_code': 'BR',
 'city': 'Salvador',
 'time_zone': 'America/Bahia',
 'area_code': 0,
 'region_code': '05',
 'latitude': -12.983300000000014,
 'metro_code': None,
 'dma_code': 0,
 'postal_code': None}

Now we can create a pyephem moon object,

In [4]:

import ephem
import numpy as np
from datetime import date

date = date(2014, 8, 10)
moon = ephem.Moon()

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

Out[4]:

<ephem.Observer date='2014/8/10 00:00:00' epoch='2000/1/1 12:00:00' lon=-38:31:00.1 lat=-12:58:59.9 elevation=0.0m horizon=0:00:00.0 temp=15.0C pressure=1010.0mBar>

and query for the local time super moon rising.

In [5]:

print("Next full Moon: %s" % ephem.next_full_moon(date))
print("Next rising: %s" % ephem.localtime(obs.next_rising(moon)))

Next full Moon: 2014/8/10 18:09:21
Next rising: 2014-08-10 17:29:43.000006

Since we could not see the moon rising let's just think positive and wait for the next one.

In [6]:

az, alt, symbols, times = [], [], [], []
for k in range(30):
    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).strftime("%d-%m-%Y"))
        symbols.append(symbol)
        alt.append(moon.alt)
        az.append(moon.az)
    obs.date += ephem.hour * 24  # Increment in days.
az, alt = map(np.rad2deg, (az, alt))

In [7]:

import matplotlib.pyplot as plt
from matplotlib.font_manager import FontProperties

prop = FontProperties(fname='./data/moon_phases.ttf', size=25)
sub = slice(0, -1, 2)
fig, ax = plt.subplots(figsize=(6, 3.75))
[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='#0B3B0B', alpha=0.75)
 for x, y, hour in zip(az[sub], alt[sub], times[sub])]

ax.set_xlim(0-10, 360+30)
ax.set_ylim(-90-5, 90+15)
ax.set_ylabel("Elevation (degrees above horizon)")
ax.set_xlabel("Azimuth (compass direction, in degrees)")
ax.grid(False)
seaborn.despine(trim=True)

And to display the result let use the moon-phases graph. This is probably a very confusing way to plot this. I should stack the plots horizontally every 360 degrees to pass the idea of time. Anyways, we can that the next full moon will be around 2014-09-08, let's hope it does not rain!

In [8]:

HTML(html)

Out[8]:

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

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/.