arcgis—python网格化批量裁剪遥感影像

一、生成网格shp

在arcgis中，找到创建渔网，具体位置在arctool box->数据管理工具->采样->创建渔网

依次输入

二、python批量裁剪代码

这是之前我看到一个博主的代码，具体是谁，我也不清楚了，只需要修改shp_path，img_path，save_path三个地方 ，改成你对应的文件夹路径，就能使用了。

# -*- coding: utf-8 -*-
# @Time : 2020/1/8 15:04
# @Author : Zhao HL
# @File : extract by mask.py
import arcpy
import os, time

shp_path = r'D:\postgraduate\test.shp'
img_path = r'D:\postgraduate\zg.tif'
save_path = r'E:\postgraduate\layer'

arcpy.CheckOutExtension("Spatial")


def clear_folder(root_path):
    '''
    clear all files in root_path(include files in subfolders, but not folders)
    :param root_path:
    :return:
    '''
    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)

    print('clear ' + root_path)


def clip_img(in_raster, clip_shp, workspace):
    '''
    according to features in the Shp, extract the raters one by one;
    all result rasters are saved at wordspace according to the FID field in Shp
    :param in_raster:
    :param clip_shp:
    :param workspace:
    :return:
    '''
    if arcpy.GetParameterAsText(0) != '':
        in_raster = arcpy.GetParameterAsText(0)
    if arcpy.GetParameterAsText(1) != '':
        clip_shp = arcpy.GetParameterAsText(1)
    if arcpy.GetParameterAsText(2) != '':
        workspace = arcpy.GetParameterAsText(2)

    clear_folder(workspace)
    arcpy.env.workspace = workspace

    t1 = time.time()
    for row in arcpy.SearchCursor(clip_shp):
        mask = row.getValue("Shape")
        FID = int(row.getValue("FID"))
        FID_name = str(FID).zfill(5) + '.tif'
        img_8_path = os.path.join(workspace, FID_name)

        # region method 1: slow but steady
        mask_raster = arcpy.sa.ExtractByMask(in_raster, mask)
        arcpy.CopyRaster_management(in_raster=mask_raster,
                                    out_rasterdataset=img_8_path,
                                    config_keyword="DEFAULTS",
                                    background_value=255,
                                    nodata_value=255,
                                    onebit_to_eightbit="",
                                    colormap_to_RGB="",
                                    pixel_type="8_BIT_UNSIGNED",
                                    scale_pixel_value=None,
                                    RGB_to_Colormap=None)
        # endregion

        # region method 2: fast but probably get null result
        # arcpy.Clip_management(in_raster, '#', img_8_path, mask, 0, "ClippingGeometry")
        # endregion

        t2 = time.time() - t1
        arcpy.AddMessage(FID_name + ' is generated successfully,  total time:' + str(round(t2, 2)))


if __name__ == "__main__":
    pass
    clip_img(img_path, shp_path, save_path)

三、运行代码

在arcgis中工具箱->我的工具箱->工具箱.tbx，右键工具箱，创建脚本，选择前面准备好的代码，其它设置均可保持默认。

四、其它方法

这个代码，没有使用arcpy，而是采用numpy来进行遥感的读取，所以裁剪出来的图像与上述的方法得到的结果不同，cut_photo用于裁剪图像，paste_photo用于拼接裁剪后的图像。

'''
coding: utf-8 -*-
@Time : 2023/9/13
@Author : Dengzhiyong
@File : extracy remote sesening
'''
import os
from PIL import Image
import cv2
import numpy as np
Image.MAX_IMAGE_PIXELS = None


def cut_photo(img_path):
    # 创建临时文件夹
    tmp = './VOC2012/DEM_zc/'
    if not os.path.exists(tmp): os.mkdir(tmp)

    # 输入图片
    # predict_img = Image.open(img_path)
    predict_img =cv2.imread(img_path, cv2.IMREAD_GRAYSCALE)
    p_size = predict_img.shape[0:2]  # 获取图像大小
    print(p_size[0], p_size[1])
    p_width, p_height = p_size[1], p_size[0]  # 获取图像宽和高
    print("宽:", p_width, "高:", p_height)
    # 需要把图片切割成若干个[256,256]的部分，以减少预测过程中的显存或内存占用
    # 计算X和Y对应切割后的数量
    w_remain, h_remain = p_width % 82, p_height % 82
    w_num = p_width // 82 if w_remain == 0 else p_width // 82 + 1
    h_num = p_height // 82 if h_remain == 0 else p_height // 82 + 1
    # 切割图片1
    right_max, bottom_max = p_width, p_height
    for y in range(0, h_num):
        y_top, y_bottom = y * 82, (y + 1) * 82
        y_bottom = y_bottom if y_bottom <= bottom_max else bottom_max
        for x in range(0, w_num):
    #         # 裁剪坐标为[y0:y1, x0:x1]dd
            x_left, x_right = x * 82, (x + 1) * 82  # 每一张图像左端边界至右端边界
            x_right = x_right if x_right <= right_max else right_max
            print(x_left,x_right,y_top,y_bottom)
    #         # (left, upper, right, lower)
    #         cut_small_img = predict_img.crop((x_left, y_top, x_right, y_bottom))
            print(y_top, y_bottom, x_left, x_right)
            cut_small_img = predict_img[y_top:y_bottom, x_left:x_right]
            # cut_small_img.imwrite(os.path.join(tmp, '{}_{}.png'.format(x, y)),)
            cv2.imwrite(os.path.join(tmp, '{}_{}.tif'.format(y, x)),cut_small_img)
    return p_width, p_height


seg_photo_size = cut_photo('zg_PODU_123.tif')
#
#
# def paste_photo(p_width, p_height):
#     # 新建一张图像
#     New_image = Image.new('RGB', (p_width, p_height))
#     # 计算循环
#     w_remain, h_remain = p_width % 512, p_height % 512
#     w_num = p_width // 512 if w_remain == 0 else p_width // 512 + 1
#     h_num = p_height // 512 if h_remain == 0 else p_height // 512 + 1
#
#     for x in range(0, h_num):
#         x_left = x * 512
#         for y in range(0, w_num):
#             y_top = y * 512
#             small_img = Image.open(
#                 os.path.join(r'E:/project/Swin-Transformer-Semantic-Segmentation-main/tools/rs128/Predict_mask/', '{}_{}.png'.format(x, y)))
#             # im.paste(ims,(0,0))
#             New_image.paste(small_img, (y_top, x_left))
#
#     return New_image
# #
# #
# img = paste_photo(seg_photo_size[0], seg_photo_size[1])
# img.save('deeplabv3+_zg_9_10.tif')