Access Oregon State TPX08 Tidal Data via OPeNDAP

To test the efficiency of viability of accessing the Oregon State TPX08 Tidal Data using OPeNDAP, I first downloaded the TPX08 NetCDF files for elevation from the home page below. The original files were NetCDF3 files, so I converted these to NetCDF4 files using nccopy, which comes with the netcdf distribution from Unidata. Using level 6 deflation, the size of the 9 files from 4014MB to 74MB without loss of information. I also used 128x128 chunks to optimize performance for extracting subregions.

Here's the actual script I used:

#!/bin/bash
for CON in m2 s2 n2 k2 k1 o1 p1 q1 m4
do
  echo $CON
  nccopy -d6 -c nx/128,ny/128 hf.${CON}_tpxo8_atlas_30c.nc ${CON}_nc4.nc
done

I then placed these NetCDF4 files in a folder scanned by our THREDDS Data Server so they could be accessed via OPeNDAP. In the python script below, I illustrate extracting the M2 amplitude and phase data using OPeNDAP for a specified bounding box.



In [1]:

    
from pylab import *
import netCDF4
import time



In [2]:

    
# OPeNDAP Data URLs for TPX08 Data on Geoport THREDDS Server
cons=['m2','n2','s2','k2','k1','o1','p1','q1','m4']
url={}
for con in cons:
    url[con]='http://geoport.whoi.edu/thredds/dodsC/usgs/data2/rsignell/data/tpx08/nc4/%s_nc4.nc' % con



In [3]:

    
url['m2']









    Out[3]:





'http://geoport.whoi.edu/thredds/dodsC/usgs/data2/rsignell/data/tpx08/nc4/m2_nc4.nc'



In [4]:

    
# Open the OPeNDAP Dataset using the NetCDF4-Python library
nc = netCDF4.Dataset(url['m2']).variables
nc.keys()









    Out[4]:





[u'hIm', u'hRe', u'lat_z', u'lon_z', u'con']



In [5]:

    
print(nc['hRe'])









    



<type 'netCDF4.Variable'>
int32 hRe(u'nx', u'ny')
    units: millimeter
    option_0: land
    option_1: ~=0 - water
    field: Re(z), scalar, amp=abs(hRe+i*hIm);GMT phase=atan2(-hIm,hRe)/pi*180;
    long_name: Tidal elevation complex amplitude, Real part
    _ChunkSize: [128 128]
unlimited dimensions = ()
current size = (10800, 5401)



In [12]:

    
# Specify region to extract [lon_min, lon_max, lat_min, lat_max]
box = [-71.5+360, -63+360, 39.5, 46.0]  # gulf of maine
box = [-75.+360,    -60.+360., 12.,  22.]



In [13]:

    
lon = nc['lon_z'][:]
lat = nc['lat_z'][:]
bi = (lon>=box[0]) & (lon<=box[1])
bj = (lat>=box[2]) & (lat<=box[3])



In [14]:

    
# Extract the complex amplitude components for this subregion
time0=time.time()
zi = nc['hIm'][bi,bj].T    # index dimensions in NetCDF file are reversed from usual [j,i] ordering, so transpose
zr = nc['hRe'][bi,bj].T
print 'elapsed time = %f seconds' % (time.time()-time0)
lon = nc['lon_z'][bi]
lat = nc['lat_z'][bj]









    



elapsed time = 0.075691 seconds



In [15]:

    
# Compute amplitude and phase from real and imaginary parts
amp = abs(zr+1j*zi)/1000.   # convert millibars to meters              
phase = (arctan2(-zi,zr)/pi*180)  
amp = ma.masked_equal(amp,0.)  # mask where amp = 0
phase = ma.masked_equal(phase,0.)  # mask where phase = 0



In [16]:

    
# Make a plot
figure(figsize=(20,8))
subplot(121,aspect=1.0/cos(lat.mean() * pi / 180.0))
pcolormesh(lon,lat,amp)
axis(box)
colorbar()
grid()
title('TPX08: M2 amplitude (m)');
subplot(122,aspect=1.0/cos(lat.mean() * pi / 180.0))
pcolormesh(lon,lat,phase)
axis(box)
colorbar()
grid()
title('TPX08: M2 phase (deg)');



In [11]: