In [18]:

    

%matplotlib inline
from IPython.display import HTML

import matplotlib.animation
import netCDF4
import numpy as np
import datetime
import matplotlib.pyplot as plt


    

In [9]:

    

url = 'http://opendap.nationaleregenradar.nl/thredds/dodsC/radar/TF0005_A/2013/10/01/RAD_TF0005_A_20131001000000.h5'
#url = 'http://opendap.nationaleregenradar.nl/thredds/dodsC/radar/TF0100_A/2013/01/01/RAD_TF0100_A_20130101000000.h5'
import logging
logger = logging.getLogger(__name__)
logger.setLevel(logging.DEBUG)
logging.basicConfig()
import os


    

In [10]:

    

VIDEO_TAG = """<video controls>
 <source src="{0}" type="video/ogg">
 Your browser does not support the video tag.
</video>"""
def display_animation(video):
    return HTML(VIDEO_TAG.format("files/" + video))


    

In [11]:

    

ds = netCDF4.Dataset(url)


    

In [12]:

    

!ncdump -h $url


    

    




netcdf RAD_TF0005_A_20131001000000 {
dimensions:
	east = 500 ;
	north = 490 ;
	time = 8928 ;
variables:
	byte available(time) ;
		available:_lastModified = "1970-01-25T08:16:08Z" ;
		available:_FillValue = 0b ;
		available:_Unsigned = "true" ;
		available:_ChunkSize = 8928 ;
	double east(east) ;
		east:_lastModified = "1970-01-25T08:16:08Z" ;
		east:_ChunkSize = 500 ;
	double north(north) ;
		north:_lastModified = "1970-01-25T08:16:08Z" ;
		north:_ChunkSize = 490 ;
	int time(time) ;
		time:long_name = "time" ;
		time:calendar = "gregorian" ;
		time:units = "seconds since 2013-10-01" ;
		time:standard_name = "time" ;
		time:_lastModified = "1970-01-06T15:16:55Z" ;
		time:_Unsigned = "true" ;
		time:_ChunkSize = 2232 ;
	float precipitation(north, east, time) ;
		precipitation:_lastModified = "1970-01-25T08:16:08Z" ;
		precipitation:_FillValue = -9999.f ;
		precipitation:_ChunkSize = 20, 20, 24 ;
}

In [14]:

    

time = netCDF4.num2date(ds.variables['time'][:], ds.variables['time'].units)
x = ds.variables['east'][:]
y = ds.variables['north'][:]
P = ds.variables['precipitation']
molen = {}
molen['idx'] = (np.abs((x - 110245)).argmin(), np.abs((y - 548771)).argmin())
pseries = P[molen['idx'][0], molen['idx'][1],:]


    

In [34]:

    

import pandas

t0 = datetime.datetime(2013,10,10,0, 0)
t1 = datetime.datetime(2013,10,15, 2, 0)
series = pandas.Series(data=pseries.filled(0), index=time)
np.diff(time)
fig, axes = plt.subplots(1,2, figsize=(20,4))
ax = axes[0]
ax.plot(series.index, series, 'k.-')
pandas.rolling_sum(series, window=24*60/5.0, center=True).plot(title="rolling day sum [mm/day]", ax=ax)
ax.fill_betweenx(ax.get_ylim(), [t0]*2, [t1]*2, alpha=0.3)
ax.set_xlim(t0, t1)

ax = axes[1]
ax.plot(series.index, series*60/5.0, 'k.-')
pandas.rolling_sum(series, window=24*60/5.0, center=True).plot(title="rolling day sum [mm/day]", ax=ax)
ax.set_xlim(t0, t1)

time_idx = np.where(np.logical_and(time >= t0, time < t1))[0]
t0_idx = time_idx.min()
t0_idx = time_idx.max()


    

In [32]:

    

rgb = matplotlib.colors.hex2color('#87CEEB')
segmentdata = {'red': [], 'green': [], 'blue': []}
for value, color in zip((0, 0.2, 0.9, 1), ('#87CEEB', '#FFFFFF', '#888888', '#444444')):
    rgb = matplotlib.colors.hex2color(color)
    segmentdata['red'].append((value, rgb[0], rgb[0]))
    segmentdata['green'].append((value, rgb[1], rgb[1]))
    segmentdata['blue'].append((value, rgb[2], rgb[2]))
Rain = matplotlib.colors.LinearSegmentedColormap('rain', segmentdata)    
segmentdata = {'red': [], 'green': [], 'blue': [], 'alpha': []}
for value, color, alpha in zip((0, 0.2, 0.9, 1), ('#FFFFFF', '#FFFFFF', '#888888', '#444444'), (0,0.3,0.7,0.9)):
    rgb = matplotlib.colors.hex2color(color)
    segmentdata['red'].append((value, rgb[0], rgb[0]))
    segmentdata['green'].append((value, rgb[1], rgb[1]))
    segmentdata['blue'].append((value, rgb[2], rgb[2]))
    segmentdata['alpha'].append((value, alpha, alpha))
Rain_alpha = matplotlib.colors.LinearSegmentedColormap('rain', segmentdata)


    

In [33]:

    

metadata = dict(title='Movie Test', artist='Matplotlib',
        comment='Movie support!')
writer = matplotlib.animation.FFMpegWriter(fps=10, metadata=metadata, codec='libtheora')

fig, ax = plt.subplots(figsize=(10,6))
P = ds.variables["precipitation"]
im = ax.pcolorfast(x, y, P[:,:,time_idx[0]].filled(0), vmin=0, vmax=2, cmap=Rain)
ax.set_xlabel('projected x coordinate [m] epsg:28992')
ax.set_ylabel('projected y coordinate [m] epsg:28992')
cb = plt.colorbar(im, ax=ax)
cb.set_label('precipitation [mm/5min]')

with writer.saving(fig, "movie.ogv", 100):
    for i in time_idx[:]:
        data = P[:,:,i].filled(0)
        im.set_data(data)
        ax.set_title(str(time[i]))
        writer.grab_frame()


    

In [21]:

    

display_animation("movie.ogv")


    

    
Out[21]:

In [ ]:

    

rgb = matplotlib.colors.hex2color('#87CEEB')

plt.imshow(P[:,:,800].filled(0), cmap=Rain, vmax=0.3)
plt.colorbar()


    

In [ ]: