ncep ww3 vs ucar ncep ww3

ncepww3.py

import matplotlib
matplotlib.use('AGG')

import os
import numpy


import xarray as xar
from matplotlib import pyplot as plt
from matplotlib import colors
from mpl_toolkits.basemap import Basemap


''' Determining if the two 'NCEP' models below are infact the same models.

Model 1: Unidata Ncep (UCAR) - http://thredds.ucar.edu/thredds/ncss/grib/NCEP/WW3/Regional_US_West_Coast/Best/dataset.html

Model 2: Ncep's Ncep - http://nomads.ncep.noaa.gov:9090/dods/wave/mww3


1. Access Datasets at the same time
2. Pull the lat, lon, and appropriate variables
3. Save 

'''
time = 1

def get_color_map():
    modified_jet_cmap_dict = {
        'red': ((0., .0, .0),
                (0.3, .5, .5),
                (0.4, .7, .7),
                (0.45, .8, .8),
                (0.5, 1, 1),
                (0.55, 1, 1),
                (0.6, 1, 1),
                (0.65, 1, 1),
                (0.85, 1, 1),
                (1, 0.4, 0.4)),
        'green': ((0., .4, .4),
                  (0.2, 1, 1),
                  (0.5, 1, 1),
                  (0.65, .7, .7),
                  (0.8, .45, .45),
                  (0.92, 0.1, 0.1),
                  (0.99, .0, .0),
                  (1, 0, 0)),
        'blue': ((0., .4, .4),
                 (0.2, 1, 1),
                 (0.4, .3, .3),
                 (0.5, .7, .7),
                 (0.6, .2, .2),
                 (0.75, 0, 0),
                 (1, 0, 0))
    }
    return colors.LinearSegmentedColormap('modified_jet', modified_jet_cmap_dict, 256)

def convert_to_degrees_west(x):
    y = 180 - x; return -(y + 180)

# Open Datasets and pull out latitude and longitude

ucar_url = 'http://thredds.ucar.edu/thredds/dodsC/grib/NCEP/WW3/Regional_US_West_Coast/WW3_Regional_US_West_Coast_20180419_0000.grib2'
ncep_url = 'http://nomads.ncep.noaa.gov:9090/dods/wave/mww3/20180419/multi_1.wc_10m20180419_00z'

ucar = xar.open_dataset(ucar_url, engine='netcdf4')
ncep = xar.open_dataset(ncep_url, engine='netcdf4')

ucar_time = ucar['time'][time]
ncep_time = ncep['time'][time]

ucar_lats = ucar['lat']
ucar_longs = ucar['lon']
ncep_lats =  ncep['lat']
ncep_longs = ncep['lon']

# Grab combined heights
ucar_heights = ucar['Significant_height_of_combined_wind_waves_and_swell_surface'][time]
ncep_heights = ncep['htsgwsfc'][time]

# Set up Basemaps and output files
fig = plt.figure()

min_period = -10  # Meters
max_period = 10 # Meters

num_color_levels = 200

contour_range = max_period - min_period
contour_range_inc = float(contour_range) / num_color_levels

color_levels = []
for i in xrange(num_color_levels + 1):
    color_levels.append(min_period + 1 + i * contour_range_inc)


# UCAR
ax = plt.subplot(221)  # one subplot in the figure

bmap = Basemap(projection='merc',  # A cylindrical, conformal projection.
               resolution='h', area_thresh=1.0,
               llcrnrlat=ucar_lats[0], urcrnrlat=ucar_lats[-1],
               llcrnrlon=ucar_longs[0], urcrnrlon=ucar_longs[-1],
               epsg=4326)

x, y = numpy.meshgrid(ucar_longs, ucar_lats)
print "UCAR Shapes: ", x.shape, y.shape, ucar_heights.shape
ucar_heights = numpy.array(ucar_heights)
bmap.drawcoastlines()
bmap.contourf(x, y, ucar_heights[::-1], color_levels, extend='both', cmap=get_color_map())
bmap.colorbar()

# NCEP
ax = plt.subplot(222)

bmap = Basemap(projection='merc',  # A cylindrical, conformal projection.
               resolution='h', area_thresh=1.0,
               llcrnrlat=ncep_lats[0], urcrnrlat=ncep_lats[-1],
               llcrnrlon=ncep_longs[0], urcrnrlon=ncep_longs[-1],
               epsg=4326)

x, y = numpy.meshgrid(ncep_longs, ncep_lats)
ncep_heights = numpy.array(ncep_heights)

print "NCEP SHAPES: ", x.shape, y.shape, ncep_heights.shape
bmap.drawcoastlines()
bmap.contourf(x, y, ncep_heights, color_levels, extend='both', cmap=get_color_map())
bmap.colorbar()

combo = plt.subplot(223)

bmap = Basemap(projection='merc',  # A cylindrical, conformal projection.
                    resolution='h', area_thresh=1.0,
                    llcrnrlat=ncep_lats[0], urcrnrlat=ncep_lats[-1],
                    llcrnrlon=ncep_longs[0], urcrnrlon=ncep_longs[-1],
                    epsg=4326)

combo_heights = ncep_heights - ucar_heights[::-1]
bmap.drawcoastlines()
bmap.contour(x, y, combo_heights, color_levels, extend='both', cmap=get_color_map())
bmap.colorbar()

combo = plt.subplot(224)

bmap = Basemap(projection='merc',  # A cylindrical, conformal projection.
               resolution='h', area_thresh=1.0,
               llcrnrlat=ncep_lats[0], urcrnrlat=ncep_lats[-1],
               llcrnrlon=ncep_longs[0], urcrnrlon=ncep_longs[-1],
               epsg=4326)

combo_heights =  ucar_heights[::-1] - ncep_heights
bmap.drawcoastlines()
bmap.contour(x, y, combo_heights, color_levels, extend='both', cmap=get_color_map())

print "UCAR TIME: ", ucar_time
print "NCEP TIME: ", ncep_time

bmap.colorbar()

fig.savefig(
    os.path.join("./out.png"),
    dpi=2000, bbox_inches='tight', pad_inches=0,
    transparent=True, frameon=False)
plt.close(fig)