使用多边形切割线在 numpy 数组中裁剪网格答案

【问题标题】：Crop the grid in a numpy array using a polygom cutline使用多边形切割线在 numpy 数组中裁剪网格
【发布时间】：2018-11-28 09:25:38
【问题描述】：

如何裁剪存储在带有 wkt 多边形的 numpy 数组中的 EPSG:25832 栅格？

【问题讨论】：

你尝试了什么？
事情是；我的一所大学编写了一些 python 脚本，它们生成了一个插值图，现在有一个使用剪切线裁剪图像的请求。我一直在研究这个解决方案：gis.stackexchange.com/questions/262021/… 但是我必须像 tiff 一样临时存储我的数据，当我准备好将我的数据放在一个 numpy 数组中时，这有点尴尬。

标签： python numpy gis gdal

【解决方案1】：

您不需要将它保存在磁盘中，您可以将其用作临时对象。您可以简单地使用“Mem”格式：

from osgeo import gdal, gdal_array as gdarr

input_raster = "path/to/rgb.tif" 
# or if the input is already a gdal raster object you can use that gdal object
input_raster=gdal.Open("path/to/rgb.tif")
input_kml = "path/to/cropline.kml" # or any other format

ds = gdal.Warp('',
              input_raster,
              format = 'Mem',
              cutlineDSName = input_kml, # or any other file format
              cutlineLayer = 'extent', # if cutline is a shapefile no need for this cutlineLayer 
              dstNodata = -9999) # select your no data value
#do stuff with ds object, it is your cropped dataset
#for example convert it to numpy array
npDs=ds.gdarr.DatasetReadAsArray(ds, 0, 0, ds.RasterXSize, ds.RasterYSize)
print(npDs.shape)

如果你有一个 numpy 数组，那就更棘手了，因为你需要知道那个 numpy 数组的 geotransform 参数。我的建议是将这个 numpy 数组发回给 gdal 对象，然后你可以使用 gdal.warp()

【讨论】：