之前应公司需求需要一个能批量裁剪输出的功能,上网研究了一下各路大神的方法,然后自己根据需要写了个程序,这里记录一下。
前言
基于ArcGIS的arcpy二次开发
一、数据准备
待裁剪影像、标准图框shp文件、ArcGIS10.2、Pycharm
二、代码
代码如下(示例):
import arcpy
import os, time
arcpy.CheckOutExtension("spatial")
typeall = {"tif": "tif", "img": "img", "bmp": "bmp", "jpg": "jpg", "png": "png"}
# arcpy.env.workspace = 'D:/yx'
# clip_Shp = 'D:/1/TFH.shp'
# fieldName = 'TFH'
# output_Path = 'D:/TF/ClipRasterResult'
# type1 = 'tif'
# workSpeed = 'HighSpeed'
# mosaic_first = 'NO'
arcpy.env.workspace = arcpy.GetParameterAsText(0)
type1 = arcpy.GetParameterAsText(1)
mosaic_first = arcpy.GetParameterAsText(2)
clip_Shp = arcpy.GetParameterAsText(3)
fieldName = arcpy.GetParameterAsText(4)
output_Path = arcpy.GetParameterAsText(5) + "\\ClipRasterResult"
workSpeed = arcpy.GetParameterAsText(6)
expression = arcpy.AddFieldDelimiters(clip_Shp, fieldName) + " <> ''"
type1 = typeall[type1]
rasters = arcpy.ListRasters("*", type1)
def create_dir(path):
if os.path.exists(path):
clear_folder(path)
else:
os.mkdir(path)
def clear_folder(root_path):
files_list = os.listdir(root_path)
for f in files_list:
file_path = os.path.join(root_path, f)
if os.path.isdir(file_path):
clear_folder(file_path)
else:
os.remove(file_path)
arcpy.AddMessage('clear ' + root_path)
def bool_point(X, Y, XMin, YMin, XMax, YMax):
if XMin < X < XMax and YMin < Y < YMax:
return True
else:
return False
def bool_extents(extent1, extent2):
extent1_XMin = extent1.XMin
extent1_XMax = extent1.XMax
extent1_YMin = extent1.YMin
extent1_YMax = extent1.YMax
extent2_XMin = extent2.XMin
extent2_XMax = extent2.XMax
extent2_YMin = extent2.YMin
extent2_YMax = extent2.YMax
if (bool_point(extent2_XMin, extent2_YMax, extent1_XMin, extent1_YMin, extent1_XMax, extent1_YMax) == True) or (
bool_point(extent2_XMax, extent2_YMax, extent1_XMin, extent1_YMin, extent1_XMax, extent1_YMax) == True) or (
bool_point(extent2_XMax, extent2_YMin, extent1_XMin, extent1_YMin, extent1_XMax, extent1_YMax) == True) or (
bool_point(extent2_XMin, extent2_YMin, extent1_XMin, extent1_YMin, extent1_XMax, extent1_YMax) == True):
return True
else:
return False
def clip_method2(in_raster, in_extent, clip_shp, in_field_name, out_path, speed):
for row in arcpy.SearchCursor(clip_shp, where_clause=expression):
t1 = time.time()
mask = row.getValue("Shape")
if bool_extents(in_extent, mask.extent):
FID_name = str(row.getValue(in_field_name)) + '.tif'
output_path = os.path.join(out_path + '/', FID_name)
if os.path.exists(output_path):
FID = str(row.getValue("FID"))
(file_path, file_name) = os.path.split(in_raster)
FID_name = str(file_name).replace('.tif', '') + '-' + FID + '-' + FID_name
output_path = os.path.join(out_path + '\\', FID_name)
if speed == "HighQuality":
mask_raster = arcpy.sa.ExtractByMask(in_raster, mask)
arcpy.CopyRaster_management(in_raster=mask_raster,
out_rasterdataset=output_path,
config_keyword="DEFAULTS",
background_value=255,
nodata_value=None,
onebit_to_eightbit="",
colormap_to_RGB="",
pixel_type="8_BIT_UNSIGNED",
scale_pixel_value=None,
RGB_to_Colormap=None)
if speed == "HighSpeed":
arcpy.Clip_management(in_raster, '#', output_path, mask, None, 'ClippingGeometry', "MAINTAIN_EXTENT")
t2 = time.time() - t1
arcpy.AddMessage(FID_name + ' is generated successfully, total time:' + str(round(t2, 2)))
if rasters:
create_dir(output_Path)
if mosaic_first == 'NO':
for raster in rasters:
r = arcpy.sa.Raster(raster)
arcpy.AddMessage(raster + "--is being processed")
clip_method2(raster, r.extent, clip_Shp, fieldName, output_Path, workSpeed)
arcpy.AddMessage(raster + "--has done")
else:
rasterList = ";".join(rasters)
r = arcpy.sa.Raster(rasters[0])
desc = arcpy.Describe(rasters[0])
outCoorSystem = desc.spatialReference
dataType = desc.DataType
cellWidth = r.meanCellWidth
bandCount = r.bandCount
arcpy.MosaicToNewRaster_management(rasterList, output_Path, "MosaicRaster.tif", outCoorSystem,"8_BIT_UNSIGNED", None, bandCount, "MAXIMUM", "FIRST")
arcpy.AddMessage('raster is mosaiced successfully!')
m = arcpy.sa.Raster(output_Path + "/MosaicRaster.tif")
clip_method2(output_Path + "/MosaicRaster.tif", m.extent, clip_Shp, fieldName, output_Path, workSpeed)
else:
arcpy.AddMessage('folder does not exist raster!')
其中注释部分为绝对路径,供测试时使用。
在ArcGIS中导入脚本并设置参数,如图所示: