利用shp点文件对影像按指定尺寸进行裁剪
现有一份shp点文件,文件包含要裁剪的中心点经纬度坐标,如下图所示,以及包含坐标点的遥感影像。现要利用已知坐标按规定大小对影像进行裁剪。本文指定的大小为256*256像素。
ArcGIS在安装时已自动安装Python2.7以及arcpy包。在用编辑器编辑的时候需要将编译环境设置为python2.7。本文使用Pycharm进行编辑。设置为python2.7的操作如下:文件->设置->项目->Project Interpreter将项目翻译设置为Python2.7即可。同时在该目录下也可查看当前Python版本安装了哪些包。
#coding=utf-8
import os
import arcpy
import string
import scipy
import os
from arcpy import env
from arcpy.sa import *
#allow overwrite the exists file
arcpy.gp.overwriteOutput = 1
#setting workspace
env.workspace = "D:/pycharm_workspace"
#利用数据集中单个波段查询波段信息,因为影像不同波段可能信息不一样
raster_band = "tif5.tif/Band_1"
desc = arcpy.Describe(raster_band) #获取波段信息
ext = desc.extent #影像范围
x_ext = ext.XMax - ext.XMin
y_ext = ext.YMax - ext.YMin
height = desc.height
width = desc.width
spatialref = desc.spatialReference
#cell size
x_mean = x_ext/width
y_mean = y_ext/height
#desc.meanCellHeight
#desc.meanCellWidth
#为shp文件设置一个搜索游标,获得X,Y坐标
fc = "D:/pycharm_workspace/Export_Output_proj.shp"
cursor = arcpy.da.SearchCursor(fc, ["FID", "X坐标", "Y坐标"])
#map_doc = arcpy.mapping.MapDocument("D:/pycharm_workspace/forclip.mxd")
#利用shp文件里的坐标,以坐标为中心,创建一个要素类,要素类大小为256*256像素,利用要素对影像进行裁剪
for row in cursor:
#print "坐标序号 = {0}, X= {1}, Y= {2}".format(row[0],row[1],row[2])
feature = "feature" + "_" + str(row[0]) + ".shp"
arcpy.CreateFeatureclass_management("D:/pycharm_workspace", feature, "Polygon","","","",spatialref)
point_array = arcpy.Array()
left_top_point = arcpy.Point()
right_top_point = arcpy.Point()
right_bottom_point = arcpy.Point()
left_bottom_point = arcpy.Point()
end_point = arcpy.Point()
#计算裁剪框四个角点坐标
left_top_point.ID = 0
left_top_point.X = row[1] - 128 * x_mean
left_top_point.Y = row[2] - 128 * y_mean
right_top_point.ID = 1
right_top_point.X = row[1] + 128 * x_mean
right_top_point.Y = row[2] - 128 * y_mean
right_bottom_point.ID = 2
right_bottom_point.X = row[1] + 128 * x_mean
right_bottom_point.Y = row[2] + 128 * y_mean
left_bottom_point.ID = 3
left_bottom_point.X = row[1] - 128 * x_mean
left_bottom_point.Y = row[2] + 128 * y_mean
#与第一个坐标相同,形成封闭的矩形
end_point.ID = 4
end_point.X = row[1] - 128 * x_mean
end_point.Y = row[2] - 128 * y_mean
point_array.add(left_top_point)
point_array.add(right_top_point)
point_array.add(right_bottom_point)
point_array.add(left_bottom_point)
#利用角点生成矩形(polygon)
feature_polygon = arcpy.Polygon(point_array)
feature_cursor = arcpy.da.InsertCursor(feature, ["SHAPE@"])
feature_cursor.insertRow([feature_polygon])
del feature_cursor
clip_feat = "D:/pycharm_workspace/tif5.tif"
outworkspace = "D:/pycharm_workspace/ccc"
outname = "D:/pycharm_workspace/ccc" + "//"+ "cilpimg" + "_" + str(row[0]) + ".bmp"
#arcpy.Clip_analysis(clip_feat, feature, outname)
#arcpy.Clip_analysis(clip_feat, feature, outname)
#arcpy.Clip_management(clip_feat, "left_top_point.X left_top_point.Y right_bottom_point.X right_bottom_point.Y",outname)
#利用生成的矩形进行裁剪
arcpy.Clip_management(clip_feat, "#", outname, feature, "#", "ClippingGeometry")
转载自:https://blog.csdn.net/weixin_40994552/article/details/85320331