[系统运维]python gdal使用shp裁剪遥感影像

最近有个需求，在云服务器上运行shp裁剪遥感影像的功能。服务器内存2G。

影像十分巨大，例如二十万行，三十万列。

目前的shp裁剪影像方法有：

? ? （1）arcgis裁剪：服务器上肯定不能用。

? ? （2）基于rasterio的mask.mask()建立掩膜裁剪。这个方法比较快，以前一直用这个，但是有个缺点就是把数据一次性读入内存，在内存不足时就报错，不能正常运行。基于rasterio的裁剪代码如下。

import fiona
import rasterio as rio
import rasterio.mask
import pyproj
import sys
from osgeo import gdal
import os

def clipRasterByShapefile(src, shpdatafile, dst,nodata=0):

    # 读取shp
    with fiona.open(shpdatafile, "r") as shapefile:
        features = [feature["geometry"] for feature in shapefile]
    # 读取原始影像
    src = rio.open(src)

    # 调用函数执行裁剪
    out_image, out_transform = rio.mask.mask(src, features,
                                             all_touched=True,
                                             crop=True,
                                             nodata=nodata)
    # 元数据信息复制
    out_meta = src.meta.copy()
    out_meta.update({"driver": "GTiff",
                     "height": out_image.shape[1],
                     "width": out_image.shape[2],
                     "transform": out_transform})
    # 输出文件
    output_file = rasterio.open(dst, "w", compress="LZW",**out_meta)
    output_file.write(out_image)
    output_file.close()

? ? （3）gdal.Warp()。这个方法不会有内存不足问题，但是速度会很慢。

from osgeo import  gdal
import time
t0 = time.time()


dst = "LUCC30巴.tif"
shp = r"utm.shp"
src = "LUCC30.tif"

ds = gdal.Warp(dst, # 裁剪图像保存完整路径（包括文件名）
                   src, # 待裁剪的影像
                    # warpMemoryLimit=500 内存大小M
                   format='GTiff', # 保存图像的格式
                   cutlineDSName = shp, # 矢量文件的完整路径
                   cropToCutline = True,
                   copyMetadata=True,
                   creationOptions=['COMPRESS=LZW',"TILED=True"])

t1 = time.time()
print("gdal time is %s"%(t1-t0))

以上方法不能满足小内存服务器上对大影像裁剪的需求。因此本文写了一个速度相对于arcgis和gdal.Warp更快所需内存更小的分块裁剪方法。思路比较简单，首先是使用gdal.Rasterize()对矢量栅格二值化，其次分块对原始影像应用二值化掩膜，写出本地。由于采用了分块处理的思想，不会有内存溢出的问题。

代码链接：clip_gdal_fast.py-Web服务器文档类资源-CSDN下载pythongdal使用shp遥感影像裁剪-基于分块裁剪策略，速度快，所需内存小，速度快于arcg更多下载资源、学习资料请访问CSDN下载频道.https://download.csdn.net/download/weixin_40450867/39710243

使用巴音郭楞蒙古自治州的矢量裁剪30米土地利用的实验表明：gdal.Warp耗时五分钟，而本代码仅需23秒。rasterio更快点，仅需要20秒，但在小内存电脑rasterio不能正常裁剪。

使用如下的矢量进行测试均达到了预期效果。

?