Accessing ncSOS 1.2 on TDS 4.5 with OWSLib


In [1]:
from owslib.sos import SensorObservationService
from owslib.etree import etree
import pandas as pd
import datetime as dt
import matplotlib.pyplot as plt
import numpy as np
%matplotlib inline

In [2]:
# usgs woods hole
# regular time series data 
url='http://geoport-dev.whoi.edu/thredds/sos/usgs/data2/emontgomery/stellwagen/CF-1.6/BUZZ_BAY/2651-A.cdf'
usgs = SensorObservationService(url)
contents = usgs.contents

In [3]:
usgs.contents


Out[3]:
{'network-all': <owslib.swe.observation.sos100.SosObservationOffering at 0x7f1818e96590>,
 'urn_ioos_station_gov.usgs.cmgp_2651-A': <owslib.swe.observation.sos100.SosObservationOffering at 0x7f1818e96950>}

In [4]:
off = usgs.offerings[1]
off.name


Out[4]:
'urn:ioos:station:gov.usgs.cmgp:2651-A'

In [5]:
off.response_formats


Out[5]:
['text/xml;subtype="om/1.0.0"',
 'text/xml;subtype="om/1.0.0/profiles/ioos_sos/1.0"']

In [6]:
off.observed_properties


Out[6]:
['http://mmisw.org/ont/fake/parameter/NEP_56',
 'http://mmisw.org/ont/cf/parameter/sea_water_pressure',
 'http://mmisw.org/ont/cf/parameter/sea_water_pressure',
 'http://mmisw.org/ont/cf/parameter/air_temperature',
 'http://mmisw.org/ont/fake/parameter/V00_1900',
 'http://mmisw.org/ont/fake/parameter/V01_1901',
 'http://mmisw.org/ont/fake/parameter/V02_1902',
 'http://mmisw.org/ont/fake/parameter/rtrn_4012']

In [7]:
#pdb.set_trace()
response = usgs.get_observation(offerings=['2651-A'],
                                 responseFormat='text/xml;subtype="om/1.0.0"',
                                 observedProperties=['http://mmisw.org/ont/cf/parameter/air_temperature'])

In [8]:
print response[0:2000]


<?xml version="1.0" encoding="UTF-8"?>
<om:ObservationCollection xmlns:om="http://www.opengis.net/om/1.0" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xmlns:xlink="http://www.w3.org/1999/xlink" xmlns:gml="http://www.opengis.net/gml" xmlns:swe="http://www.opengis.net/swe/1.0" xsi:schemaLocation="http://www.opengis.net/om/1.0 http://schemas.opengis.net/om/1.0.0/observation.xsd">
  <gml:description>2651-A         start time = 20 Aug 1982  11:31:52  cycles = 27617
               stop  time = 31 Oct 1982  09:31:52  # days = 71
Expt. = 'USGS'  sampling interval = 3.75 minutes
Lat =  41.650890
Lon = -70.687553  File created:  10:01 SEP  1,'87
depth = 14m  Mag.var = -15.000000
#  Variable     Units        Codes Depth   Inst.   Minimum   Maximum
-- ------------ ------------ ----- ------- ------ --------- ---------
 1 TIME         SECONDS      T     
 2 TEMPERATURE  DEG.CELC     R       13.00           12.436    21.954
 3 PRESSURE-I   MILLIBARS    R       13.00         2018.077  2204.906
 4 PRESSURE-B   MILLIBARS    R       13.00         2017.956  2204.963
 5 PSDEV                     R       13.00            0.000    10.370
 6 TRANSON      VOLTS        R       13.00           -0.001     4.996
 7 EXT_COEF                  R       13.00            0.003    24.858
 8 TRANSOFF     VOLTS        R       13.00           -0.203     4.996
 9 NEPHON       VOLTS        R       13.00            0.277     4.996
10 NEPHOFF      VOLTS        R       13.00            0.241     4.996
---------------------------------------------------------------------
 Comments:                                                      MATHOP
 BUZZARDS BAY VARIABLES LISTED IN ORDER SHOWN ABOVE</gml:description>
  <gml:name>BUZZ_BAY (265)</gml:name>
  <gml:boundedBy>
    <gml:Envelope srsName="http://www.opengis.net/def/crs/EPSG/0/4326">
      <!-- overwrite these with your actual offering ROI -->
      <gml:lowerCorner>41.6508903503418 -70.68755340576172</gml:lowerCorner>
   

In [9]:
def parse_om_xml(response):
    # extract data and time from OM-XML response
    root = etree.fromstring(response)
    # root.findall(".//{%(om)s}Observation" % root.nsmap )
    values = root.find(".//{%(swe)s}values" % root.nsmap )
    date_value = np.array( [ (dt.datetime.strptime(d,"%Y-%m-%dT%H:%M:%SZ"),float(v))
            for d,v in [l.split(',') for l in values.text.split()]] )
    time = date_value[:,0]
    data = date_value[:,1]
    return data,time

In [10]:
data, time = parse_om_xml(response)

In [11]:
ts = pd.Series(data,index=time)

In [12]:
ts.plot(figsize=(12,4));


Try adding a subset on time:


In [13]:
start = '1982-10-01T00:00:00Z'
stop = '1982-10-04T00:00:00Z'
response = usgs.get_observation(offerings=['2651-A'],
                                 responseFormat='text/xml;subtype="om/1.0.0"',
                                 observedProperties=['http://mmisw.org/ont/cf/parameter/air_temperature'],
                                 eventTime='%s/%s' % (start,stop))

In [14]:
data, time = parse_om_xml(response)
ts = pd.Series(data,index=time)
ts.plot(figsize=(12,4));