Check gdalwarp using mercantile shapes --mercator (no CRS) as cutline

[sgillies]

Follow up from #14 (comment).

[sgillies]

Using the Rasterio testing data.

Cutline straight from mercantile shapes

mercantile shapes "[72, 109, 8]" --mercator \
| gdalwarp ~/code/rasterio/rasterio/tests/data/RGB.byte.tif /tmp/no_crs_72_109_8.tif -crop_to_cutline -cutline /vsistdin/ -t_srs EPSG:3857

produces a properly georeferenced but empty dataset.

Gdalwarp got the right extents from the cutline but didn't transfer any pixels.

Cutline translated by ogr2gr

mercantile shapes "[72, 109, 8]" --mercator \
| ogr2ogr /vsistdout/ /vsistdin/ -f GeoJSON -a_srs EPSG:3857 \
| gdalwarp ~/code/rasterio/rasterio/tests/data/RGB.byte.tif /tmp/72_109_8.tif -crop_to_cutline -cutline /vsistdin/ -t_srs EPSG:3857

produces a dataset with the same dimensions and extents as the one above, but with the data I expect:

[sgillies]

Conclusion

gdalwarp doesn't need to find an SRS in a cutline datasource to crop the output to the cutline's extents, but it does need to find an SRS in the cutline in order for it to be a useful mask.

Unless you really want cutline masking, it's better to use gdalwarp --te west south east north.

I'll make sure mercantile can produce such output.

Also: support for stdin and stdout in GDAL and OGR utilities is super handy, though I wish they would support /dev/stdin and /dev/stdout on Linux.

[sgillies]

The new --extents option fits into gdalwarp like this:

extent=$(mercantile shapes "[72, 109, 8]" --mercator --extents)
gdalwarp -t_srs EPSG:3857 -te $extent $source $dest