I am trying to clip a raster file with a multi-polygon shapefile to get separate clipped-raster file for each polygon with its 'id' as part of the name of clipped-raster file. I tried the following code by it throws an error.
code:
import ogr
import subprocess
shapefile = ("/shape.shp")
in_raster = ("/whole.tif")
ds = ogr.Open(shapefile)
lyr = ds.GetLayer(0)
lyr.ResetReading()
ft = lyr.GetNextFeature()
while ft:
shape_id= ft.GetFieldAsString('id')
out_raster = in_raster.replace('.tif', '_%s.tif' % shape_id.replace(' ', '_'))
subprocess.call(['gdalwarp', in_raster, out_raster, '-cutline', shapefile,
'-crop_to_cutline', '-cwhere', "'id'='%s'" % shape_id])
ft = lyr.GetNextFeature()
ds = None
Error:
ERROR 1: Did not get any cutline features
Best Answer
It's often easier to use the GDAL python wrappers. See my attempt below. If you still get errors, check that both raster and shapefile have geographic information (i.e. you can overlap them in QGIS or something like that).