[GIS] Trying to clip raster in QGIS, “Cutline is not Valid”

clipgdalgdalwarpqgis

I am using QGIS 2.14.5-Essen and I am trying to simply clip an image using a square polygon. I am encountering an error that I have never seen before. This is what I am attempting to do:

enter image description here

"clipper" is a shapefile containing a rectangular polygon feature over my desired study area.

The image I am trying to clip is a LANDSAT-8 Surface Reflectance image. When I load it into QGIS, it is already in spatial reference system

+proj=utm +zone=35 +datum=WGS84 +units=m +no_defs

and I verified that it is showing up in the correct area by loading in the OSM basemap and confirming that it is oriented properly.

My "clipper" shapefile is also in

+proj=utm +zone=35 +datum=WGS84 +units=m +no_defs

However, when I try to run the operation,which corresponds to the gdal command

gdalwarp -q -cutline D:/african_towns_identifying/data/original_image/clipper.shp -tr 30.0 30.0 -of GTiff D:/african_towns_identifying/data/original_image/LC81720712016212LGN00_sr_band2.tif D:/african_towns_identifying/data/clipped_image/band2.tif

I get this error:

enter image description here

I don't understand why it just won't clip normally, because I have done this before without any problems. I don't understand the error at all and cannot find anything helpful with a Google search

Trying to clip the image manually using the "Extent" in the dialog box and dragging a box will execute the tool successfully, but I get an image ranging from NaN to NaN

Best Answer

(Option 1) You could try the option in gdalwarp to ignore invalid polygons. This is invoked using "--config GDALWARP_IGNORE_BAD_CUTLINE YES":

gdalwarp -of GTiff -co "COMPRESS=DEFLATE" -tr 30 30 -cutline VectorName -crop_to_cutline InputRaster OutputRaster --config GDALWARP_IGNORE_BAD_CUTLINE YES

(Option 2) Since your polygon is a square, you could read the shapefile in Python/OGR to get the bounding area of the layer/square, and then use the layer coordinates to clip your raster:

from osgeo import ogr, gdal
import subprocess

# Define I/O files
InputImage = "Input.tif"
OutImage = "Out.tif"
Shapefile = "SomeName.shp"`

# Open the shapefile in read-only mode and retrieve the spatial reference
Driver = ogr.GetDriverByName("ESRI Shapefile")
VectorDataset = Driver.Open(Shapefile, 0)
layer = VectorDataset.GetLayer()
feature = layer.GetNextFeature()
geom = feature.GetGeometryRef()
spatialRef = geom.GetSpatialReference()`

# Get bounding area of the layer
minX, maxX, minY, maxY = geom.GetEnvelope()`

# Clip the InputImage to the bounding area of the input feature (using Nearest Neighbour interpolation).
OutTile = gdal.Warp(OutImage, InputImage, format="GTiff", outputBounds=[minX, minY, maxX, maxY], xRes=30, yRes=30, dstSRS=spatialRef, resampleAlg=gdal.GRA_NearestNeighbour)`

# Build overviews for the OutImage
subprocess.call("gdaladdo --config COMPRESS_OVERVIEW DEFLATE "+OutImage+" 2 4 8 16 32 64", shell=True)
print("Done.")
Related Question