WCS GetCoverage for Grid returning CSV Example

From Earth Science Information Partners (ESIP)

Simple WCS Python Example

This example uses popular python scripting language to call Datafed WCS via SOAP. See WCS page SOAPifying_WCS

Getting python

Windows users: go to http://www.activestate.com and download Python 2.4 or later. You can also go straight to download

The main site for python is http://www.python.org but I prefer to install via the activestate package. Especially Mark Hammonds PythonWin is a nice IDE.

Getting Demo Script

Download http://datafed.net/demo/soapdemo.txt and save it with extension .py

type from command line python soapdemo.py and watch it run.

There are no parameters, the script is hard coded for demo purposes.

How it works

After some declarations there is the WCS query wrapped in SOAP envelope.

The triple quotes allow multiline strig constants in python.

We have now variable example_soap_in

#   this WCS query queries USA temperatures 2005-12-06T12:00:00 at 1350 m height.
#   The Envelope and Body nodes are standard SOAP wrapping.
#   The actual query is the http-post query defined by WCS schema.

example_soap_in = """
<soap:Envelope xmlns:soap="http://schemas.xmlsoap.org/soap/envelope/">
    <soap:Body>
        <wcs:GetCoverage version="1.0.0" service="WCS"
                xmlns:gml="http://www.opengis.net/gml" xmlns:wcs="http://www.opengis.net/wcs">
            <!-- datafed dataset_abbr.param_abbr defines the coverage name -->
            <wcs:sourceCoverage>THREDDS_GFS.T</wcs:sourceCoverage> 
            <wcs:domainSubset>
                <wcs:spatialSubset>
                    <!--
                        This element queries the USA at 1350 m elevation.
                        If the dataset has no elevation, only lat and lon limits
                        are needed. 
                    -->
                    <gml:Envelope srsName="WGS84(DD)">
                        <!-- lon_min lat_min elev_min -->
                        <gml:pos>-130 24 1350</gml:pos>
                        <!-- lon_max lat_max elev_max -->
                        <gml:pos>-65 50 1350</gml:pos>
                    </gml:Envelope>
                    <gml:Grid dimension="3">
                        <gml:limits>
                            <!--
                                grid size. Datafed does not support resizing, so
                                these numbers have no meaning for grid datasets.
                            -->
                            <gml:GridEnvelope>
                                <gml:low>0 0 0</gml:low>
                                <gml:high>99 99 0</gml:high>
                            </gml:GridEnvelope>
                        </gml:limits>
                        <gml:axisName>lat</gml:axisName>
                        <gml:axisName>lon</gml:axisName>
                        <gml:axisName>elev</gml:axisName>
                    </gml:Grid>
                </wcs:spatialSubset>
                <wcs:temporalSubset>
                    <!--
                        query data for one time only.
                        Datafed currently supports also time ranges,
                        but only for stored time interval.
                        Enumerated times are not yet supported.                        
                    -->
                    <gml:timePosition>2005-12-06T12:00:00</gml:timePosition>
                </wcs:temporalSubset>
            </wcs:domainSubset>
            <wcs:output>
                <!--
                    For test purposes, use Comma Separated Values.
                    NetCDF and GML are also supported.
                -->
                <wcs:format>CSV</wcs:format>
            </wcs:output>
        </wcs:GetCoverage>
    </soap:Body>
</soap:Envelope>
"""

Function query_datafed_wcs makes a WCS query to Datafed. Remove the tracing print statements if you use this for real.

It does it by doing all the low level http stuff: open connection, send headers, send data, inspect return, parse return xml, get the real data with another http-get request.

# This function performs WCS SOAP query to datafed
# INPUT: the soap query, similar to example_soap_in
# output: the data stream

def query_datafed_wcs(soap_in_text):
    # open connection to std port 80
    c = httplib.HTTPConnection("webapps.datafed.net", 80)
    c.connect()
    try:
        # the Web Service Definition Language definition is at
        # http://webapps.datafed.net/WCS.asmx?WSDL
        c.putrequest("POST", "/WCS.asmx")
        c.putheader("soapAction", "http://datafed.net/WCS/GetCoverage")
        c.putheader("content-type", "text/xml")
        c.putheader("content-length", repr(len(soap_in_text)))
        c.endheaders()
        c.send(soap_in_text)
        r = c.getresponse()
        print 'response status ' + repr(r.status) + ' ' + repr(r.reason)

        # 200 means OK. Anything other is a failure.
        if r.status <> 200:
            raise 'http-error', repr(r.status) + ' ' + repr(r.reason)

        x = xml.dom.minidom.parse(r)
        print x.toxml()

        # We have the output. It does not contain data, just
        # metadata and an uri pointer to the real data.
        try:
            # The uri is written into 'data' node, so let's look for it.
            for u in x.getElementsByTagNameNS("http://datafed.net/xs/wcs", "data"):
                print 'fetching data ' + u.firstChild.nodeValue
                return urllib2.urlopen(urllib2.Request(u.firstChild.nodeValue))
            # We did not find the uri, so the result is an error            
            raise 'soap-error', x
        finally:
            # due to circular references, you have to have the next line.
            # Otherwise memory will leak. Other python implementations
            # like Jython or IronPython rely on platform provided
            # garbage collection, and don't need it. 
            x.unlink()
    finally:
        #closing expensive resources ASAP is a good practice
        c.close()

Python has no main program, but if invoked from command line the special variable __name__ is set to __main__. So if you issue command python soapdemo.py then the following section will be executed.

Alternative way to run this is to

  • start python
  • type import soapdemo without the .py extension
  • type csv = query_datafed_wcs(example_soap_in) You now have http stream named csv.
  • print it by print csv.read()
#
#   If you run this script from command line, the next will be executed
#
if __name__ == "__main__":
    print 'querying datafed'
    csv = query_datafed_wcs(example_soap_in)
    print "here's the data"
    print ''
    # we have the csv stream, so let's look at it
    print csv.read()