# IOOS System Test: Extreme Events Theme: Inundation

### Can we estimate the return period of a water level by comparing modeled and/or observed water levels with NOAA Annual Exceedance Probability Curves?

Methodology:

• Define temporal and spatial bounds of interest, as well as parameters of interest
• Search for available service endpoints in the NGDC CSW catalog, then inform the user of the DAP (model) and SOS (observation) services endpoints available
• Obtain the stations in the spatial boundaries, and processed to obtain observation data for temporal constraints, identifying the yearly max
• Plot observation stations on a map and indicate to the user if the minimum number of years has been met for extreme value analysis (red marker if condition is false)
• Using DAP (model) endpoints find all available models data sets that fall in the area of interest, for the specified time range, and extract a model grid cell closest all the given station locations
• Plot the annual max for each station as a timeseries plot
• Plot the annual exceedance probability curve for the Nantucket Island, Ma station and compare to NOAA Tides and Currents plot
• The annual exceedance probability curve can be used to estimate the return period of a particular water level, whether it be modeled or observed data.

#### import required libraries

``````

In [1]:

from datetime import datetime, timedelta
import requests
from warnings import warn

import folium
from IPython.display import Image
import iris
from iris.exceptions import CoordinateNotFoundError, ConstraintMismatchError
import matplotlib.pyplot as plt
import numpy as np
from owslib import fes
# from shapely.geometry import Point  # For lat lon points.
from owslib.csw import CatalogueServiceWeb
import pandas as pd
import prettyplotlib as ppl
from pyoos.collectors.coops.coops_sos import CoopsSos
import scipy.stats as stats

from utilities import (get_Coops_longName, fes_date_filter, get_coordinates, service_urls, inline_map, coops2data, find_timevar,
find_ij, nearxy, mod_df)

``````

#### Specify Temporal and Spatial conditions.

``````

In [2]:

# Bounding box of interest,[bottom right[lat,lon], top left[lat,lon]].
bounding_box_type = "box"
# Build a bounding box around Long Island Sound
bounding_box = [[-75.94, 39.67],
[-66.94, 41.5]]

# Temporal range. Get 30+ yrs obs data for extremal analysis and recent modeled data to compare
start_date = datetime(1980, 5, 1).strftime('%Y-%m-%d %H:00')
end_date = datetime(2014, 5, 1).strftime('%Y-%m-%d %H:00')
time_date_range = [start_date, end_date]

jd_now = datetime.utcnow()
jd_start,  jd_stop = jd_now - timedelta(days=3), jd_now + timedelta(days=3)

jd_start = datetime.strptime(jd_start.strftime('%Y-%m-%d %H:00'), '%Y-%m-%d %H:%M')
jd_stop = datetime.strptime(jd_stop.strftime('%Y-%m-%d %H:00'), '%Y-%m-%d %H:%M')

print('%s to %s' % (start_date, end_date))

``````
``````

1980-05-01 00:00 to 2014-05-01 00:00

``````

#### Specify a list of standard names to search for

``````

In [3]:

name_list = ['water level',
'sea_surface_height',
'sea_surface_elevation',
'sea_surface_height_above_geoid',
'sea_surface_height_above_sea_level',
'water_surface_height_above_reference_datum',
'sea_surface_height_above_reference_ellipsoid',
'sea_surface_height_above_reference_level']

sos_name = 'water_surface_height_above_reference_datum'

``````

#### Define the CSW geoportal endpoint

``````

In [4]:

endpoint = 'http://www.ngdc.noaa.gov/geoportal/csw'  # NGDC Geoportal.
csw = CatalogueServiceWeb(endpoint, timeout=60)

for oper in csw.operations:
if oper.name == 'GetRecords':
print('\nISO Queryables:\n%s' %
'\n'.join(oper.constraints['SupportedISOQueryables']['values']))

# Put the names in a dict for ease of access.
data_dict = {}
data_dict["water"] = {"names": name_list,
"sos_name": sos_name}

``````
``````

ISO Queryables:
apiso:Subject
apiso:Title
apiso:Abstract
apiso:AnyText
apiso:Format
apiso:Identifier
apiso:Modified
apiso:Type
apiso:BoundingBox
apiso:CRS.Authority
apiso:CRS.ID
apiso:CRS.Version
apiso:RevisionDate
apiso:AlternateTitle
apiso:CreationDate
apiso:PublicationDate
apiso:OrganizationName
apiso:HasSecurityConstraints
apiso:Language
apiso:ResourceIdentifier
apiso:ParentIdentifier
apiso:KeywordType
apiso:TopicCategory
apiso:ResourceLanguage
apiso:GeographicDescriptionCode
apiso:Denominator
apiso:DistanceValue
apiso:DistanceUOM
apiso:TempExtent_begin
apiso:TempExtent_end
apiso:ServiceType
apiso:ServiceTypeVersion
apiso:Operation
apiso:OperatesOn
apiso:OperatesOnIdentifier
apiso:OperatesOnName
apiso:CouplingType

``````

#### Convert User Input into FES filters.

``````

In [5]:

start, stop = fes_date_filter(start_date, end_date)
box = []
box.append(bounding_box[0][0])
box.append(bounding_box[0][1])
box.append(bounding_box[1][0])
box.append(bounding_box[1][1])
bbox = fes.BBox(box)

or_filt = fes.Or([fes.PropertyIsLike(propertyname='apiso:AnyText',
literal=('*%s*' % val),
escapeChar='\\',
wildCard='*',
singleChar='?') for val in name_list])
val = 'Averages'
not_filt = fes.Not([fes.PropertyIsLike(propertyname='apiso:AnyText',
literal=('*%s*' % val),
escapeChar='\\',
wildCard='*',
singleChar='?')])

``````
``````

In [6]:

filter_list = [fes.And([bbox, start, stop, or_filt, not_filt])]
# connect to CSW, explore it's properties
# try request using multiple filters "and" syntax: [[filter1,filter2]]
csw.getrecords2(constraints=filter_list,
maxrecords=1000,
esn='full')

``````
``````

In [7]:

# Print records that are available:
print("number of datasets available: %s" % len(csw.records.keys()))

``````
``````

number of datasets available: 417

``````

Print all the records (should you want too).

``````

In [8]:

# print("\n".join(csw.records))

``````

#### DAP URLs:

``````

In [9]:

dap_urls = service_urls(csw.records, service='odp:url')

# Remove duplicates and organize.
dap_urls = sorted(set(dap_urls))
print("Total DAP: %s" % len(dap_urls))

# print the dap urls...
print("\n".join(dap_urls))

``````
``````

Total DAP: 106
http://comt.sura.org/thredds/dodsC/data/comt_1_archive/inundation_tropical/VIMS_SELFE/Hurricane_Ike_2D_final_run_with_waves
http://comt.sura.org/thredds/dodsC/data/comt_1_archive/inundation_tropical/VIMS_SELFE/Hurricane_Ike_2D_final_run_without_waves
http://comt.sura.org/thredds/dodsC/data/comt_1_archive/inundation_tropical/VIMS_SELFE/Hurricane_Ike_3D_final_run_with_waves
http://comt.sura.org/thredds/dodsC/data/comt_1_archive/inundation_tropical/VIMS_SELFE/Hurricane_Ike_3D_final_run_without_waves
http://comt.sura.org/thredds/dodsC/data/comt_1_archive/inundation_tropical/VIMS_SELFE/Hurricane_Rita_2D_final_run_with_waves
http://comt.sura.org/thredds/dodsC/data/comt_1_archive/inundation_tropical/VIMS_SELFE/Hurricane_Rita_2D_final_run_without_waves
http://comt.sura.org/thredds/dodsC/data/comt_1_archive/inundation_tropical/VIMS_SELFE/Hurricane_Rita_3D_final_run_with_waves
http://comt.sura.org/thredds/dodsC/data/comt_1_archive/inundation_tropical/VIMS_SELFE/Hurricane_Rita_3D_final_run_without_waves
http://data.ioos.us/thredds/dodsC/rutgers/ru23-20131017T1614.nc3.nc.ncml
http://geoport-dev.whoi.edu/thredds/dodsC/estofs/atlantic
http://geoport.whoi.edu/thredds/dodsC/coawst_4/use/fmrc/coawst_4_use_best.ncd
http://oos.soest.hawaii.edu/thredds/dodsC/dist2coast_1deg
http://oos.soest.hawaii.edu/thredds/dodsC/dist2coast_4deg
http://oos.soest.hawaii.edu/thredds/dodsC/hioos/tide_pac
http://oos.soest.hawaii.edu/thredds/dodsC/pacioos/hycom/global
http://sos.maracoos.org/stable/dodsC/stationHRPIER84-agg.ncml
http://sos.maracoos.org/stable/dodsC/stationHRPMNTH-agg.ncml
http://sos.maracoos.org/stable/dodsC/stationHRWSTPTH-agg.ncml
http://tds.marine.rutgers.edu/thredds/dodsC/roms/espresso/2013_da/his_Best/ESPRESSO_Real-Time_v2_History_Best_Available_best.ncd
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/028p1/028p1_d03.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/028p1/028p1_d07.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/028p1/028p1_d08.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/029p1/029p1_d06.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/029p1/029p1_d13.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/029p1/029p1_historic.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/036p1/036p1_d32.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/045p1/045p1_d12.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/067p1/067p1_d16.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/067p1/067p1_historic.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/071p1/071p1_d07.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/071p1/071p1_d08.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/071p1/071p1_historic.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/091p1/091p1_d14.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/092p1/092p1_d12.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/094p1/094p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/094p1/094p1_d04.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/094p1/094p1_d05.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/094p1/094p1_d06.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/094p1/094p1_d07.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/094p1/094p1_d10.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/094p1/094p1_d11.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/094p1/094p1_d12.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/094p1/094p1_historic.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/095p1/095p1_d05.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/095p1/095p1_d06.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/096p1/096p1_d08.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/096p1/096p1_d09.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/096p1/096p1_historic.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/098p1/098p1_d06.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/098p1/098p1_d08.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/098p1/098p1_d12.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/100p1/100p1_d02.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/100p1/100p1_d03.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/100p1/100p1_d05.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/100p1/100p1_historic.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/101p1/101p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/102p1/102p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/111p1/111p1_d02.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/111p1/111p1_d07.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/114p1/114p1_d04.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/118p1/118p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/126p1/126p1_d02.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/128p1/128p1_d07.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/139p1/139p1_d02.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/139p1/139p1_d06.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/139p1/139p1_historic.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/151p1/151p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/151p1/151p1_historic.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/152p1/152p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/152p1/152p1_d02.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/154p1/154p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/154p1/154p1_d02.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/154p1/154p1_d03.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/154p1/154p1_d04.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/154p1/154p1_d05.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/154p1/154p1_historic.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/156p1/156p1_d02.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/162p1/162p1_d02.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/162p1/162p1_d03.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/166p1/166p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/168p1/168p1_d05.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/172p1/172p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/174p1/174p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/175p1/175p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/175p1/175p1_d02.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/175p1/175p1_d03.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/175p1/175p1_d04.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/179p1/179p1_combo_tmp.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/179p1/179p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/179p1/179p1_d02.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/179p1/179p1_d03.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/179p1/179p1_d04.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/195p1/195p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/204p1/204p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/207p1/207p1_d01.nc
http://thredds.cdip.ucsd.edu/thredds/dodsC/cdip/archive/207p1/207p1_historic.nc

``````

#### SOS URLs:

``````

In [10]:

sos_urls = service_urls(csw.records, service='sos:url')

# Remove duplicates and organize.
sos_urls = sorted(set(sos_urls))
print("Total SOS: %s" % len(sos_urls))
print("\n".join(sos_urls))

``````
``````

Total SOS: 79
http://data.ioos.us/thredds/sos/rutgers/ru23-20131017T1614.nc3.nc.ncml?service=SOS&version=1.0.0&request=GetCapabilities
http://opendap.co-ops.nos.noaa.gov/ioos-dif-sos/SOS?service=SOS&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/028p1/028p1_d03.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/028p1/028p1_d07.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/028p1/028p1_d08.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/029p1/029p1_d06.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/029p1/029p1_d13.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/029p1/029p1_historic.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/036p1/036p1_d32.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/045p1/045p1_d12.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/067p1/067p1_d16.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/067p1/067p1_historic.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/071p1/071p1_d07.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/071p1/071p1_d08.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/071p1/071p1_historic.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/091p1/091p1_d14.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/092p1/092p1_d12.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/094p1/094p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/094p1/094p1_d04.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/094p1/094p1_d05.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/094p1/094p1_d06.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/094p1/094p1_d07.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/094p1/094p1_d10.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/094p1/094p1_d11.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/094p1/094p1_d12.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/094p1/094p1_historic.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/095p1/095p1_d05.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/095p1/095p1_d06.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/096p1/096p1_d08.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/096p1/096p1_d09.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/096p1/096p1_historic.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/098p1/098p1_d06.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/098p1/098p1_d08.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/098p1/098p1_d12.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/100p1/100p1_d02.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/100p1/100p1_d03.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/100p1/100p1_d05.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/100p1/100p1_historic.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/101p1/101p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/102p1/102p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/111p1/111p1_d02.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/111p1/111p1_d07.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/114p1/114p1_d04.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/118p1/118p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/126p1/126p1_d02.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/128p1/128p1_d07.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/139p1/139p1_d02.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/139p1/139p1_d06.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/139p1/139p1_historic.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/151p1/151p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/151p1/151p1_historic.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/152p1/152p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/152p1/152p1_d02.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/154p1/154p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/154p1/154p1_d02.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/154p1/154p1_d03.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/154p1/154p1_d04.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/154p1/154p1_d05.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/154p1/154p1_historic.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/156p1/156p1_d02.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/162p1/162p1_d02.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/162p1/162p1_d03.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/166p1/166p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/168p1/168p1_d05.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/172p1/172p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/174p1/174p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/175p1/175p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/175p1/175p1_d02.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/175p1/175p1_d03.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/175p1/175p1_d04.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/179p1/179p1_combo_tmp.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/179p1/179p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/179p1/179p1_d02.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/179p1/179p1_d03.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/179p1/179p1_d04.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/195p1/195p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/204p1/204p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/207p1/207p1_d01.nc?service=SOS&version=1.0.0&request=GetCapabilities
http://thredds.cdip.ucsd.edu/thredds/sos/cdip/archive/207p1/207p1_historic.nc?service=SOS&version=1.0.0&request=GetCapabilities

``````

### SOS Requirements

#### Use Pyoos SOS collector to obtain Observation data from COOPS.

``````

In [11]:

start_time = datetime.strptime(start_date, '%Y-%m-%d %H:%M')
end_time = datetime.strptime(end_date, '%Y-%m-%d %H:%M')

``````
``````

In [12]:

iso_start = start_time.strftime('%Y-%m-%dT%H:%M:%SZ')
iso_end = end_time.strftime('%Y-%m-%dT%H:%M:%SZ')

collector = CoopsSos()
collector.set_datum('NAVD')
collector.server.identification.title
collector.start_time = start_time
collector.end_time = end_time
collector.variables = [data_dict["water"]["sos_name"]]

``````
``````

In [13]:

print("Date: %s to %s" % (iso_start, iso_end))
box_str = ','.join(str(e) for e in box)
print("Lat/Lon Box: %s" % box_str)

``````
``````

Date: 1980-05-01T00:00:00Z to 2014-05-01T00:00:00Z
Lat/Lon Box: -75.94,39.67,-66.94,41.5

``````
``````

In [14]:

# Grab the opendap sos url and use it for the service.
for sos_url in sos_urls:
if 'opendap' in sos_url:
break

url = (sos_url.split("?")[0] + '?'
'service=SOS&request=GetObservation&version=1.0.0&'
'observedProperty=%s&'
'offering=urn:ioos:network:NOAA.NOS.CO-OPS:WaterLevelActive&'
'featureOfInterest=BBOX:%s&responseFormat=text/tab-separated-values&'
'eventTime=%s') % (sos_name, box_str, iso_end)

r = requests.get(url)
data = r.text
# Get the headers for the col.
data = data.split("\n")
station_list_dict = dict()
# Parse the headers so I can create a dict.
c = 0
field = h.split(":")[0].split(" ")[0]
station_list_dict[field] = {"id": c}
c += 1

``````
``````

In [15]:

# Create dict of stations.
station_list = []
for i in range(1, len(data)):
station_info = data[i].split("\t")
station = dict()
for field in station_list_dict.keys():
col = station_list_dict[field]["id"]
if col < len(station_info):
station[field] = station_info[col]
station["type"] = "obs"
if "latitude" not in station:
continue
station_list.append(station)

``````

#### Go through stations and parse data

``````

In [16]:

station_yearly_max = []
for s in station_list:
if s["type"] is "obs" and s['station_id']:  # If its an obs station.
# Get the long name.
sta = s['station_id'].split(':')[-1]
s["long_name"] = get_Coops_longName(sta)
s["station_num"] = str(s['station_id']).split(':')[-1]
# This is different than sos name, hourly height is hourly water level.
s["data"] = coops2data(collector, s["station_num"], "high_low")

``````
``````

In [17]:

# Store data in a pandas DataFrame
col = 'water_surface_height_above_reference_datum (m)'

obs_df = []
for s in station_list:
years = s["data"].keys()
xx = []
yx = []
for y in years:
xx.append(int(y))
val = s["data"][y]["max"]
yx.append(val)
col = 'Observed Data'
df = pd.DataFrame(np.array(yx),columns=[col],index=years)
df.name = s["long_name"]
obs_df.append(df)
obs_df[-1].name = df.name

``````

### Plot stations on an interactive map

``````

In [18]:

m.circle_marker(location=[s["latitude"], s["longitude"]],
popup=popup_string, fill_color='#ff0000',

``````
``````

In [19]:

# Number of years required for analysis
num_years_required = 30

#find center of bounding box
lat_center = abs(bounding_box[1][1] - bounding_box[0][1])/2 + bounding_box[0][1]
lon_center = abs(bounding_box[0][0]-bounding_box[1][0])/2 + bounding_box[0][0]
m = folium.Map(location=[lat_center, lon_center], zoom_start=6)

for s in station_list:
if "latitude" in s:
if len(s["data"].keys()) >= num_years_required:
popup_string = ('<b>Station:</b><br>' + str(s['station_id']) +
"<br><b>Long Name:</b><br>" +
str(s["long_name"]))
m.simple_marker(location=[s["latitude"], s["longitude"]],
popup=popup_string)
else:
popup_string = ('<b>Not Enough Station Data for Num of years'
'requested</b><br><br>Num requested:' +
str(num_years_required) +
'<br>Num Available:' +
str(len(s["data"].keys())) +
'<br><b>Station:</b><br>' +
str(s['station_id']) +
"<br><b>Long Name:</b><br>" +
str(s["long_name"]))

m.line(get_coordinates(bounding_box,bounding_box_type), line_color='#FF0000', line_weight=5)

inline_map(m)

``````
``````

Out[19]:

``````

### Creates a time series plot with stations that have enough data

``````

In [20]:

# Set the random seed for consistency
np.random.seed(12)

fig, ax = plt.subplots()

# Show the whole color range
for s in station_list:
if "data" in s:
years = s["data"].keys()
# Only show the stations with enough data.
if len(s["data"].keys()) >= num_years_required:
xx = []
yx = []
for y in years:
xx.append(int(y))
val = s["data"][y]["max"]
yx.append(val)

ax.scatter(xx, yx, marker='o')
ppl.scatter(ax, xx, yx, alpha=0.8, edgecolor='black',
linewidth=0.15, label=str(s["station_num"]))
ppl.legend(ax, loc='right', ncol=1)
ax.set_xlabel('Year')
ax.set_ylabel('water level (m)')

ax.set_title("Stations exceeding " +
str(num_years_required) +
" years worth of water level data (MHHW)")
fig.set_size_inches(14, 8)

``````
``````

``````

### Number of stations available by number of years

``````

In [21]:

fig, ax = plt.subplots(1)
year_list_map = []
for s in station_list:
if "data" in s:
years = s["data"].keys()
year_list_map.append(len(years))

ppl.hist(ax, np.array(year_list_map), grid='y')
plt.plot([num_years_required, num_years_required], [0, 8], 'r-', lw=2)
ax.set_ylabel("Number of Stations")
ax.set_xlabel("Number of Years Available")
ax.set_title('Number of available stations vs available years\n'
'(for bounding box) - red is minimum requested years')

``````
``````

Out[21]:

<matplotlib.text.Text at 0x11f900e90>

``````

## Get model output from OPeNDAP URLS

### Try to open all the OPeNDAP URLS using Iris from the British Met Office. If we can open in Iris, we know it's a model result.

#### Construct an Iris contraint to load only cubes that match the standard names:

``````

In [22]:

print data_dict['water']['names']
name_in_list = lambda cube: cube.standard_name in data_dict['water']['names']
constraint = iris.Constraint(cube_func=name_in_list)

# station_list[0]

``````
``````

['water level', 'sea_surface_height', 'sea_surface_elevation', 'sea_surface_height_above_geoid', 'sea_surface_height_above_sea_level', 'water_surface_height_above_reference_datum', 'sea_surface_height_above_reference_ellipsoid', 'sea_surface_height_above_reference_level']

``````
``````

In [23]:

# Create time index for model DataFrame
ts_rng = pd.date_range(start=jd_start, end=jd_stop, freq='6Min')
ts = pd.DataFrame(index=ts_rng)

#Get the station lat/lon into lists and create list of model DataFrames for each station
obs_lon = []
obs_lat = []
model_df = []
for sta in station_list:
obs_lon.append(float(sta['longitude']))
obs_lat.append(float(sta['latitude']))
model_df.append(pd.DataFrame(index=ts.index))
model_df[-1].name = sta['long_name']

# Use only data within 0.04 degrees (about 4 km).
max_dist = 0.04
# Use only data where the standard deviation of the time series exceeds 0.01 m (1 cm).
# This eliminates flat line model time series that come from land points that should have had missing values.
min_var = 0.01
for url in dap_urls:
try:
# convert to units of meters
# a.convert_units('m')     # this isn't working for unstructured data
# take first 20 chars for model name
mod_name = a.attributes['title'][0:20]
r = a.shape
timevar = find_timevar(a)
lat = a.coord(axis='Y').points
lon = a.coord(axis='X').points
jd = timevar.units.num2date(timevar.points)
start = timevar.units.date2num(jd_start)
istart = timevar.nearest_neighbour_index(start)
stop = timevar.units.date2num(jd_stop)
istop = timevar.nearest_neighbour_index(stop)

# Only proceed if we have data in the range requested.
if istart != istop:
nsta = len(station_list)
if len(r) == 3:
print('[Structured grid model]:', url)
d = a[0, :, :].data
# Find the closest non-land point from a structured grid model.
if len(lon.shape) == 1:
lon, lat = np.meshgrid(lon, lat)
j, i, dd = find_ij(lon, lat, d, obs_lon, obs_lat)
for n in range(nsta):
# Only use if model cell is within 0.04 degree of requested
# location.
if dd[n] <= max_dist:
arr = a[istart:istop, j[n], i[n]].data
if arr.std() >= min_var:
c = mod_df(arr, timevar, istart, istop,
mod_name, ts)
name = station_list[n]['long_name']
model_df[n] = pd.concat([model_df[n], c], axis=1)
model_df[n].name = name

elif len(r) == 2:
print('[Unstructured grid model]:', url)
# Find the closest point from an unstructured grid model.
index, dd = nearxy(lon.flatten(), lat.flatten(),
obs_lon, obs_lat)
for n in range(nsta):
# Only use if model cell is within 0.04 degree of requested
# location.
if dd[n] <= max_dist:
arr = a[istart:istop, index[n]].data
if arr.std() >= min_var:
c = mod_df(arr, timevar, istart, istop,
mod_name, ts)
name = station_list[n]['long_name']
model_df[n] = pd.concat([model_df[n], c], axis=1)
model_df[n].name = name
elif len(r) == 1:
print('[Data]:', url)
except (ValueError, RuntimeError, CoordinateNotFoundError,
ConstraintMismatchError) as e:
warn("\n%s\n" % e)
pass

``````
``````

/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1196: UserWarning: Ignoring netCDF variable 'ele' invalid units 'non-dimensional'
warnings.warn(msg.format(msg_name, msg_units))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1196: UserWarning: Ignoring netCDF variable 'Cs' invalid units 'non-dimensional'
warnings.warn(msg.format(msg_name, msg_units))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1196: UserWarning: Ignoring netCDF variable 'sigma' invalid units 'non-dimensional'
warnings.warn(msg.format(msg_name, msg_units))
-c:80: UserWarning:
failed to merge into a single cube.
cube.long_name differs: u'free_surface_elevation' != u'maximum_free_surface_elevation'
cube.var_name differs: u'elev' != u'elevmax'
cube.attributes keys differ: 'coverage_content_type'
cube.shape differs: (720, 592761) != (592761,)

-c:80: UserWarning:
Missing data dimensions for multi-valued coordinate u'time'

-c:80: UserWarning:
failed to merge into a single cube.
cube.long_name differs: u'maximum_free_surface_elevation' != u'free_surface_elevation'
cube.var_name differs: u'elevmax' != u'elev'
cube.attributes keys differ: 'coverage_content_type'
cube.shape differs: (592761,) != (168, 592761)

-c:80: UserWarning:
failed to merge into a single cube.
cube.long_name differs: u'maximum_free_surface_elevation' != u'free_surface_elevation'
cube.var_name differs: u'elevmax' != u'elev'
cube.attributes keys differ: 'coverage_content_type'
cube.shape differs: (592761,) != (1008, 592761)

/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/cf.py:269: UserWarning: Missing CF-netCDF ancillary data variable u'lat_qc', referenced by netCDF variable u'precise_lat'
warnings.warn(message % (name, nc_var_name))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/cf.py:269: UserWarning: Missing CF-netCDF ancillary data variable u'lon_qc', referenced by netCDF variable u'precise_lon'
warnings.warn(message % (name, nc_var_name))
-c:80: UserWarning:
no cubes found

/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1196: UserWarning: Ignoring netCDF variable 'nbdv' invalid units 'nondimensional'
warnings.warn(msg.format(msg_name, msg_units))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1196: UserWarning: Ignoring netCDF variable 'neta' invalid units 'nondimensional'
warnings.warn(msg.format(msg_name, msg_units))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1196: UserWarning: Ignoring netCDF variable 'ibtypee' invalid units 'nondimensional'
warnings.warn(msg.format(msg_name, msg_units))

('[Unstructured grid model]:', 'http://geoport-dev.whoi.edu/thredds/dodsC/estofs/atlantic')

/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1196: UserWarning: Ignoring netCDF variable 'nvell' invalid units 'nondimensional'
warnings.warn(msg.format(msg_name, msg_units))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1196: UserWarning: Ignoring netCDF variable 'ibtype' invalid units 'nondimensional'
warnings.warn(msg.format(msg_name, msg_units))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1196: UserWarning: Ignoring netCDF variable 'nbvv' invalid units 'nondimensional'
warnings.warn(msg.format(msg_name, msg_units))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1196: UserWarning: Ignoring netCDF variable 'nvel' invalid units 'nondimensional'
warnings.warn(msg.format(msg_name, msg_units))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1196: UserWarning: Ignoring netCDF variable 'nvdll' invalid units 'nondimensional'
warnings.warn(msg.format(msg_name, msg_units))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/cf.py:1038: UserWarning: Ignoring formula terms variable u'h' referenced by data variable u'v' via variable u's_rho': Dimensions (u'eta_rho', u'xi_rho') do not span (u'time', u's_rho', u'eta_v', u'xi_v')
warnings.warn(msg)
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/cf.py:1038: UserWarning: Ignoring formula terms variable u'zeta' referenced by data variable u'v' via variable u's_rho': Dimensions (u'time', u'eta_rho', u'xi_rho') do not span (u'time', u's_rho', u'eta_v', u'xi_v')
warnings.warn(msg)
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/cf.py:1038: UserWarning: Ignoring formula terms variable u'h' referenced by data variable u'u' via variable u's_rho': Dimensions (u'eta_rho', u'xi_rho') do not span (u'time', u's_rho', u'eta_u', u'xi_u')
warnings.warn(msg)
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/cf.py:1038: UserWarning: Ignoring formula terms variable u'zeta' referenced by data variable u'u' via variable u's_rho': Dimensions (u'time', u'eta_rho', u'xi_rho') do not span (u'time', u's_rho', u'eta_u', u'xi_u')
warnings.warn(msg)
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1291: UserWarning: Gracefully filling 'time' dimension coordinate masked points
warnings.warn(msg.format(str(cf_coord_var.cf_name)))

('[Structured grid model]:', 'http://oos.soest.hawaii.edu/thredds/dodsC/hioos/tide_pac')
('[Structured grid model]:', 'http://oos.soest.hawaii.edu/thredds/dodsC/pacioos/hycom/global')

-c:80: UserWarning:
NetCDF: Malformed or inaccessible DAP DDS

Exception

Exception

Exception

Exception

Exception

Exception

Exception

/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1337: UserWarning: Failed to create 'waveTime' dimension coordinate: The points array must be strictly monotonic.
Gracefully creating 'waveTime' auxiliary coordinate instead.
error=e_msg))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1337: UserWarning: Failed to create 'sstTime' dimension coordinate: The points array must be strictly monotonic.
Gracefully creating 'sstTime' auxiliary coordinate instead.
error=e_msg))
/Users/bobfratantonio/anaconda/envs/waves_notebook/lib/python2.7/site-packages/Iris-1.8.0_DEV-py2.7.egg/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1337: UserWarning: Failed to create 'dwrTime' dimension coordinate: The points array must be strictly monotonic.
Gracefully creating 'dwrTime' auxiliary coordinate instead.
error=e_msg))

``````

### Plot all the model Time Series data near each obs station

``````

In [24]:

for df in model_df:
if 'Nantucket Island, MA' in df.name:
ax = df.plot(figsize=(14, 6), title=df.name, legend=False)
plt.setp(ax.lines[0], linewidth=4.0, color='0.7', zorder=1)
ax.legend()
ax.set_ylabel('m')
break
model_max = list(set(df.max().values))

``````
``````

``````

## Extreme Value Analysis:

#### Plot return periods for the observed data at the Nantucket Island, Ma station

``````

In [25]:

# Get the annual maximums from the DataFrame
for df in obs_df:
if 'Nantucket Island, MA' in df.name:
annual_max_levels = df["Observed Data"].values
name = df.name

``````

#### Fit data to GEV distribution

``````

In [26]:

def sea_levels_gev_pdf(x):
return stats.genextreme.pdf(x, xi, loc=mu, scale=sigma)

``````
``````

In [27]:

mle = stats.genextreme.fit(sorted(annual_max_levels), 0)
mu = mle[1]
sigma = mle[2]
xi = mle[0]
print("The mean, sigma, and shape parameters are %s, %s, and %s, resp." %
(mu, sigma, xi))

``````
``````

The mean, sigma, and shape parameters are 0.536747818328, 0.150074678247, and 0.00121088811397, resp.

``````

#### Probability Density Plot

``````

In [28]:

min_x = min(annual_max_levels) - 0.5
max_x = max(annual_max_levels) + 0.5
x = np.linspace(min_x, max_x, num=100)
y = [sea_levels_gev_pdf(z) for z in x]

fig = plt.figure(figsize=(12, 6))
axes = fig.add_axes([0.1, 0.1, 0.8, 0.8])
xlabel = (name + " - Annual max water level (m)")
axes.set_title("Probability Density & Normalized Histogram")
axes.set_xlabel(xlabel)
axes.plot(x, y, color='Red')
axes.hist(annual_max_levels,
bins=np.arange(min_x, max_x, abs((max_x-min_x)/10)),
normed=1, color='Yellow')

``````
``````

Out[28]:

(array([ 0.        ,  0.        ,  0.86931324,  1.88351203,  1.59374094,
0.57954216,  0.        ,  0.14488554,  0.        ]),
array([-0.163 ,  0.0342,  0.2314,  0.4286,  0.6258,  0.823 ,  1.0202,
1.2174,  1.4146,  1.6118]),
<a list of 9 Patch objects>)

``````

#### Return Value Plot

``````

In [29]:

noaa_station_id = 8449130
Image(url='http://tidesandcurrents.noaa.gov/est/curves/high/' +
str(noaa_station_id)+'.png')

``````
``````

Out[29]:

``````
``````

In [30]:

Image(url='http://tidesandcurrents.noaa.gov/est/images/color_legend.png')

``````
``````

Out[30]:

``````
``````

In [31]:

fig = plt.figure(figsize=(20, 6))
axes = fig.add_axes([0.1, 0.1, 0.8, 0.8])
T = np.r_[1:250]
sT = stats.genextreme.isf(1./T, 0, mu, sigma)
axes.semilogx(T, sT, 'r')
N = np.r_[1:len(annual_max_levels)+1]
Nmax = max(N)
axes.plot(Nmax/N, sorted(annual_max_levels)[::-1], 'bo')
title = name
axes.set_title(title)
axes.set_xlabel('Return Period (yrs)')
axes.set_ylabel('Meters above MHHW')
axes.set_xticklabels([0, 1, 10, 100, 1000])
axes.set_xlim([0, 260])
axes.set_ylim([0, 1.8])
axes.grid(True)

# Now plot the max height from the event onto this plot
event = np.ones(len(T))*model_max
axes.plot(T, event, 'g--')

``````
``````

Out[31]:

[<matplotlib.lines.Line2D at 0x11f5f8810>]

``````

This plot does not match exactly. NOAA's curves were calculated using the Extremes Toolkit software package in R whereas this notebook uses scipy. There is a python package based on the Extremes Toolkit called pywafo but this is experimental and isn't building properly on Mac OS X.