Python遥感开发之批量掩膜和裁剪

2024-07-11 14:27| 来源: 网络整理| 查看: 265

Python遥感开发之批量掩膜和裁剪 1.使用arcpy进行批量掩膜1.1 批量掩膜代码1.2 单个掩膜代码 2.使用GDAL进行批量掩膜3.使用rasterio进行批量裁剪

前言：主要介绍了使用arcpy、gdal、rasterio对遥感影像进行批量掩膜和裁剪。

1.使用arcpy进行批量掩膜

注意：arcpy是基于python2.7版本的，也就是arcgis里面自带的，使用的时候添加arcgis的python2.7，即可使用arcpy。环境如图所示：在这里插入图片描述

1.1 批量掩膜代码 filepath：需要掩膜的tif文件夹outfile ：掩膜好的tif保存的文件夹inMaskData ：可以是shp，可以是tifunicode：因为python2.7下直接支持中文路径，使用unicode可以支持中文路径 # coding=UTF-8 import arcpy,os,glob from arcpy import env if __name__ == '__main__': filepath = unicode(r"D:\0000Desktop\AAAAAwork\研究方向\data\河南NDVI", "utf-8") outfile = unicode(r"D:\0000Desktop\AAAAAwork\研究方向\data\鹿邑NDVI", "utf-8") inMaskData = unicode(r"D:\0000Desktop\AAAAAwork\研究方向\矢量图\鹿邑县\luyixian.shp", "utf-8") #shp文件可以换成tif文件 arcpy.CheckOutExtension("Spatial") env.workspace = filepath os.chdir(filepath) names = glob.glob("*.tif") for name in names: arcpy.gp.ExtractByMask_sa(name, inMaskData, outfile + "\\" + name.split(".")[0] + ".tif") print(name + '-----掩膜成功') #删除多余的文件 for file_i in glob.glob(os.path.join(outfile, '*.xml')): os.remove(file_i) for file_i in glob.glob(os.path.join(outfile, '*.tfw')): os.remove(file_i)

在这里插入图片描述

1.2 单个掩膜代码 # coding=UTF-8 import arcpy if __name__ == '__main__': filepath = unicode(r"D:\0000Desktop\AAAAAwork\研究方向\data\河南NDVI\modis_ndvi_2019_10_08.tif", "utf-8") outfile = unicode(r"D:\0000Desktop\AAAAAwork\研究方向\data\鹿邑NDVI\new_modis_ndvi_2019_10_08.tif", "utf-8") inMaskData = unicode(r"D:\0000Desktop\AAAAAwork\研究方向\矢量图\鹿邑县\luyixian.shp", "utf-8") #shp文件可以换成tif文件 arcpy.CheckOutExtension("Spatial") arcpy.gp.ExtractByMask_sa(filepath, inMaskData, outfile)

在这里插入图片描述

2.使用GDAL进行批量掩膜

注意：GDAL是基于Python3.x系列的，需要自己安装GDAL，使用的时候arcpy使用的时候一定要区分开。

from osgeo import gdal import os import glob def mask_gdal(inMaskData,filepath,outfile): """ GDAL掩膜提取 :param inMaskData: 圈选范围的路径 shp文件可以换成tif文件 :param filepath: 要掩膜的tif文件 :param outfile: 掩膜好的tif文件 :return: """ dataset = gdal.Open(filepath) # 打开遥感影像 gdal.Warp(outfile, dataset, cutlineDSName=inMaskData, cropToCutline=True) # 按掩膜提取 print(filepath + '-----掩膜成功') if __name__ == '__main__': filepath = r"D:\0000Desktop\AAAAAwork\研究方向\data\河南NDVI" # 要裁剪的tif文件所在的文件夹 outfile = r"D:\0000Desktop\AAAAAwork\研究方向\data\鹿邑NDVI" inMaskData = r"D:\0000Desktop\AAAAAwork\研究方向\矢量图\鹿邑县\luyixian.shp" # 圈选范围的路径 shp文件可以换成tif文件 os.chdir(filepath) names = glob.glob("*.tif") # 读取文件 for name in names: out = outfile + "\\" + name.split(".")[0] + ".tif" # 按照圈选范围提取出的影像所存放的路径 mask_gdal(inMaskData,name,out)

在这里插入图片描述

3.使用rasterio进行批量裁剪

注意：rasterio也是基于python3.x系列的，需要自己安装rasterio包

import os import glob import rasterio as rio import rasterio.mask import numpy as np from geopandas import GeoDataFrame def clip_raster(inMaskData,filepath,outfile): """ raster裁剪提取 :param inMaskData:圈选范围的路径 shp文件可以换成tif文件 :param filepath:要裁剪的tif文件 :param outfile:裁剪好的tif文件 :return: """ #shp转geojson shpdata = GeoDataFrame.from_file(inMaskData) geo = shpdata.geometry[0] feature = [geo.__geo_interface__] #裁剪 src = rio.open(filepath) out_image, out_transform = rio.mask.mask(src, feature, crop=True, nodata=np.nan) out_meta = src.meta.copy() out_meta.update({"driver": "GTiff","height": out_image.shape[1],"width": out_image.shape[2],"transform": out_transform}) band_mask = rasterio.open(outfile, "w", **out_meta) band_mask.write(out_image) print(filepath + '-----裁剪成功') if __name__ == '__main__': filepath = r"D:\0000Desktop\AAAAAwork\研究方向\data\河南NDVI" # 要裁剪的tif文件所在的文件夹 outfile = r"D:\0000Desktop\AAAAAwork\研究方向\data\鹿邑NDVI" inMaskData = r"D:\0000Desktop\AAAAAwork\研究方向\矢量图\鹿邑县\luyixian.shp" # 圈选范围的路径 shp文件可以换成tif文件 os.chdir(filepath) names = glob.glob("*.tif") # 读取文件 for name in names: out = outfile + "\\" + name.split(".")[0] + ".tif" # 按照圈选范围提取出的影像所存放的路径 clip_raster(inMaskData,name,out)

在这里插入图片描述

【本文地址】

公司简介

联系我们