qgis – How to Save WMS Layer to Local Raster in PyQGIS

pyqgispythonqgisrasterwms

How would one save a WMS view to a local raster from the python command line?

I've seen the answer to Download WMS to harddrive and that works for me but I'd like to know how to replicate that process exactly from the console.

ie How to progamatically set the layer, the extent to save, and then save that extent for that layer to a tif?

All of the searches I did threw up lots of results of people suggesting that downloading via wfs was the best solution. Or the linked solution via the gui. But none pointed me to the OWSLib solution.

Here is the working code I made on the back of the answer from @nmtoken I got.

from owslib.wms import WebMapService
from osgeo import gdal, osr
import processing

url='http://<mapserverurl>/MapServer/WMSServer?service=WMS'
folder = '~/output'
wms = WebMapService(url, version='1.1.1')
wms_layers= list(wms.contents)


for layerno in range(0,len(wms_layers)):
    raw_tiff=os.path.join(folder,'_raw'+str(layerno)+'.tif')
    georeferenced_tiff=os.path.join(folder,'_georeferenced'+str(layerno)+'.tif')
    defaultepsg=4326 
    #Note that available epsgs can be determined with:
    print wms[wms_layers[0]].crsOptions
    srs_string='EPSG:'+str(defaultepsg)
    lonmin=1
    lonmax=2
    latmin=4
    latmax=5
    Xpix=400
    Ypix=300
    #Note output options can be determined with:
    print wms.getOperationByName('GetMap').formatOptions
    img = wms.getmap( layers=[wms_layers[layerno]], styles=['default'], srs=srs_string, bbox=( lonmin,latmin , lonmax,latmax), size=(Xpix, Ypix), format='image/tiff' )
    out = open(raw_tiff, 'wb')
    out.write(img.read())
    out.close()
    out = None

    #Georeferencing
    src_ds = gdal.Open(raw_tiff)
    format = "GTiff"
    driver = gdal.GetDriverByName(format)
    dst_ds = driver.CreateCopy(georeferenced_tiff, src_ds, 0)
    gt = [lonmin, (lonmax-lonmin)/Xpix, 0, latmax, 0, -(latmax-latmin)/Ypix]
    dst_ds.SetGeoTransform(gt)
    srs = osr.SpatialReference()
    srs.ImportFromEPSG(defaultepsg)
    dest_wkt = srs.ExportToWkt()
    dst_ds.SetProjection(dest_wkt)
    dst_ds = None
    src_ds = None

Best Answer

You can use OWSLib to work with WMS and other OGC services.

OWSLib is a Python package for client programming with Open Geospatial Consortium (OGC) web service (hence OWS) interface standards, and their related content models.

ref: OsGeo projects pages

There are a number of posts on this site related to working with OWSlib so hopefully enough to get you started.

and then save that extent for that layer to a tif?

You will only be able to do this directly if the WMS service itself supports TIFF as an output format (as advertised in the GetCapabilities response); though many WMS do, so probably won't be an issue.

Related Question