Adding Tides to Wave Chart
Will I ever be able to rest on my laurels? I was barely finished with my wave chart that I decided it would be even more awesome if it showed me the tides. Of course, the tides for the past are not all too exciting, but if you want to know when to go, today’s tides are mighty fine.
So, to adjust the chart, I needed to find a library that handles tide. There are probably more than just one, but this one seemed to do the job just fine, so I used it. Figuring it out was not hard, thanks to examples provided, and adding the functionality to the chart was mostly trivial. Here is the 411.
The library I used is called pytides. You can simply install it by running
sudo pip install pytides
It will run off after you provide it with your password and gather all dependencies and then install them onto your system. If you already installed matplotlib, then you have numpy already on your system. If not, pip will install numpy along with pytides, which can take forever. If you are on a Debian/Ubuntu system, you are better off installing numpy first, and then using the line above.
The changes to the code are in several places:
1. We need an import statement for pytides. We need a few additional libraries, so they all go into a long change.
from pytides.tide import Tide
import pytides.constituent as cons
import numpy as np
from pyquery import PyQuery
import urllib
I added the second one because the author’s example uses the shorthand cons, and I didn’t want to change anything there, but you can remove it. Same deal for numpy.
I decided to add an HTML parsing library, in this case pyquery, because it simplifies reading the HTML. You don’t strictly need it, but it’s convenient and virtually painless.
2. We need to figure out the format of the tide station we are interested in. To do that, we first select the tide station from the NOAA site, which sadly doesn’t offer a “nearest station to lat/lon” function. Instead, you look on the map and pick the one you like (La Jolla is station 9410230). Next we retrieve the station tide metadata and store it into a file (so we don’t always go back for that purpose).
If the file is there, we read it. If not, we fetch it. Then we take the contents of the file and look up the constituents. You may not know what those are, so go and read the NOAA site on the matter (PDF file – read only if on Linux 😉 ).
tidestation='9410230'
urlfile = os.path.expanduser('~/.waveplot-%s.html' % tidestation)
if not os.path.exists(urlfile):
data = urllib.urlretrieve('<a href="http://tidesandcurrents.noaa.gov/harcon.html?unit=0&timezone=0&id=%s'">http://tidesandcurrents.noaa.gov/harcon.html?unit=0&timezone=0&id=%s'</a> % tidestation, urlfile)
txt = open(urlfile,'r').read()
html=PyQuery(txt)
tb=html('table.table-striped')
The table we just loaded has the format: Index, Name, Amplitude, Phase, Speed, Description. Sadly, the Name component doesn’t always match what pytides thinks it should be. Not sadly, though, the indexes are not sequential but are tied uniquely to the Names. So we will use those to map to pytides constituents.
indexs = []
phases = []
amplis = []
if tb:
trs=tb('tr')
for tri,trx in enumerate(trs):
tr=trs.eq(tri)
tds=tr('td')
if not tds:
continue
indexs.append(tds.eq(0).text())
amplis.append(tds.eq(2).text())
phases.append(tds.eq(3).text())
Here i create three lists: one for the indexes, one for the phases, and one for the amplitudes. The indexes into the lists are the same, so that the constituent with the index given at position 4 in the first list will match the phase and amplitude at positions 4.
Notice how the variable tb contains a pyquery object that can be used much like jQuery, the HMTL/JavaScript library. The request html(‘table.table-striped’) gives us the table element with class “table-striped”. The request tb(‘tr’) gives us a collection of tr elements (table rows) inside the table.
Notice how I used enumerate(trs) instead of just looping over them. That’s because the objects returned by the iterator on a collection are (mysteriously) not the trs, but something else. If you want the actual tr object, you need to request it with the operator .eq(). That seems stupid, and it’s one of the things that Python people call, “not pythonic.” I used to balk at the term, but I now realize that “pythonic” is simply short-hand for “predictable and easy to use.” I can get behind that!
Once we have the table row, we get the table cells (td). Here we use the same trick, requesting the cell with tr.eq(index). We know the indices (0 for index, 2 for phase, 3 for amplitude) from the discussion above.
At this point, we have the constituents listed on the page. Turns out not all stations publish all consituents, so we have to build a model that leave the other ones empty. To do so, we first look up which of the constituents that NOAA defines (conveniently stored in pytides.constituents.noaa) is present:
constituents = [c for i,c in enumerate(cons.noaa) if str(i+1) in indexs]
Now, that looks complicated, but it’s actually fairly simple. We enumerate the constituents from the collection. Then we check if the index is present in our list of indices. If so, we add it to the list of constituents. Notice how this is done by adding one to the index, since NOAA indices start at 1. Also, in the list indexs we stored strings, not numbers, so we have to convert using the standard str() function.
3. Now we can generate the Tide object. We will feed it out model:
model = np.zeros(len(constituents), dtype = Tide.dtype)
model['constituent'] = constituents
model['amplitude'] = amplis
model['phase'] = phases
tide = Tide(model = model, radians = False)
What this does is create a numpy array of constituents, filled with default Tide types. Then we fill in the model with the constituent and the corresponding amplitude and phases. Finally, we generate the Tide object and specify that the phases are not in radians (which presumably means they are in degrees).
4. Now we generate the tide heights. To do so, we need a little more than just a mapping, because we are interested in future tides, too. So we take the last 24 hours of wave data and project it ahead of us, adding a series of data points that are not present in the wave heights.
tidetimes=times[:i]
last = tidetimes[0]
for t in tidetimes:
if (last - t).days < 1:
tidetimes.append(t + datetime.timedelta(1))
else:
break
tidetimes.sort()
tides = tide.at(tidetimes)
See? First we get the set of times we used for the other charts. Then we go through the times in our list, and if they are less than a day old, we add that same time but a day later to our chart. That will add a whole day’s worth of times to our chart and will use the same time resolution as the chart overall.
At the end, a stupid bug in pytides prevents it from dealing with unsorted times, so we have to sort them before we feed them to the function. Once that is done, the call to tide.at() handles a list of datetimes perfectly fine, and returns a list of tide values at that time.
5. Finally, we have to add the tides to the chart. That’s a little involved, because we are not just adding new values to the chart, we are also expanding the X axis.
Remember how we created a twin axes object to display the wave heights in feet? Well, I stole that axes and repurposed it for my evil intents. In essence, all I had to do was to plot the tide heights on this second axes object, then use the min_x and max_x properties on it to reset the corresponding properties on the primary axes. That way, the wave data moves to the left and tide and waves are on the same time scale.
I also modified the y min/max on this tide axes because the tides swing wild and it’s hard to follow them if they go all the way. You can modify how small you want them to be, I let them roam in the center half of the plot (you could put them at the bottom, at the top, or however you like, by just shifting the ylim).
That’s pretty much it. My surfing buddy wanted to have the times of the extremes plotted, but that would be too easy and not particularly exciting, so I am leaving it out for now.