WMS for Geoscience Australia ARG25

In my last post I showed how to use Python and OWSLib to access the Geoscience Australia ARG25 web services. WMS is just as easy and here's a code snippet for that.

Accessing the Geoscience Australia ARG25 OGC Services from Python

Geoscience Australia (GA) recently released a beta of their Australian Reflectance Grid 25 (ARG25) product using Open Geospatial Consortium (OGC) web services.
The services are comprised of the following:

• a Catalogue Service for the web (CSW) allowing users to easily search ISO 19115 compliant metadata, which ARG25 data adheres to;
• Web Map Services (WMS) allowing users to obtain images of the data as quick looks to determine if the full resolution image will meet requirements; and
• Web Coverage Services (WCS) allowing users to obtain the full resolution data.

 GA also provide a web mapping portal you can use to search and download scences, but the best way is "my way" - and the web services allow you to do exactly that, search and connect directly to the archive using your own tools.

My tool of choice for today is Python using OWSLib 0.7.

OWSLib provides support for CSW, WMS and WCS. I won't be using WMS today as I want actual data from the WCS, not a pretty picture which is what the WMS provides.

Step 1: Searching the catalogue and extracting the WCS URLs

OWSLib makes this trivial. The GA CSW url is used to create a CatalogueServiceWeb object, pick your bounding box (in this case the Australian Capital Territory), then get the 'full' metadata records. 'full' records are needed to ensure we get the various web services URIs and not just a 'summary' metadata record. After that there is a bit of searching through the records URIs to find the WCS endpoints for each scene in the bounding box.

Note: In this example I've limited the number of return records to 5. There are something like 130,000 of them in the full catalogue. The metadatarecords also contain timebounds for the acquisition time of the imagery. OWSLib 0.7-dev has preliminary support for Filter construction including time spans. I've used the latest stable release of OWSLib for my example and haven't constructed a time filter. I may try it with 0.7-dev later - if you should happen to try it let me know and I'll add the extra example.

Step 2: Grab the data from the WCS

Each Coverage service contains multiple bands (e.g. 'Band1', 'Band7' ). You can also obtain related information about the support file formats the WCS can return, and the bounding box for the scene.
The code snippet below loops through the WCS urls obtained in Step 1, prints out some of this information and then gets the scene data for Band7 in GeoTiff format, saving them to a sequence of files as they go.



The GeoTIFF files (foo0.tif, foo1.tif...) can be opened using your favourite viewer. Here's the result for foo2.tif:

If you're not familiar with Landsat you may think foo0 and foo3 are failing as there is a stripping effect. That is the real data, not a problem with the code or services - I'll leave it to the reader to go find out what happen to cause "Landsat 7 stripping".

That's all there is to it, except error handling, unit tests, tidier python code (maybe this blog should be called "rusty middle aged coder").