前言

为了创建深度学习的样本，需要做一个512*512标准的框，框内具有要素的tif文件

1.以特定点为中心，做一个1024*1024的框，截取tif图像

import geopandas as gpd
import rasterio
from rasterio.plot import show


def crop_tif_around_centroid(tif_path, shp_path, output_tif_path, crop_size=1024):
    # 读取shapefile文件以获取中心点坐标
    gdf = gpd.read_file(shp_path)
    centroid = gdf.iloc[0]['geometry'].centroid

    # 读取tif文件
    with rasterio.open(tif_path) as src:
        # 获取tif文件的地理变换参数
        transform = src.transform

        # 将中心点的地理坐标转换为像素坐标
        x_center, y_center = rasterio.transform.xy(transform, centroid.x, centroid.y)

        # 计算裁剪窗口的左上角和右下角像素坐标
        x_min = x_center - (crop_size / 2)
        y_min = y_center - (crop_size / 2)
        x_max = x_center + (crop_size / 2)
        y_max = y_center + (crop_size / 2)

        # 确保裁剪窗口不超出tif图像的边界
        x_min = max(0, min(x_min, src.width - 1))
        y_min = max(0, min(y_min, src.height - 1))
        x_max = min(src.width, max(x_max, 0))
        y_max = min(src.height, max(y_max, 0))

        # 创建裁剪窗口
        window = ((int(y_min), int(y_max)), (int(x_min), int(x_max)))

        # 读取裁剪窗口内的图像数据
        out_img = src.read(window=window)

        # 更新元数据
        out_meta = src.meta.copy()
        out_meta.update({
            "height": y_max - y_min,
            "width": x_max - x_min,
            "transform": rasterio.transform.from_origin(x_min, y_max, transform[0], transform[4])
        })

        # 保存裁剪的tif文件
        with rasterio.open(output_tif_path, 'w', **out_meta) as dst:
            dst.write(out_img)

        show(out_img)

# 使用函数
tif_path = r'E:\jwztestnew\杭州滨江.tif'
shp_path = r'E:\jwztestnew\newshp_5584_centroid.shp'
output_tif_path = r'E:\jwztestnew\杭州滨江_cropped1.tif'
crop_tif_around_centroid(tif_path, shp_path, output_tif_path)

加强版：提取shp文件中各个多边形的中心点输出为shp，然后提取各个点位以该点位的“TBBM”值为索引，查找文件夹下与之对应的tif，对相应的tif进行裁剪，获得512*512的矩形框

第一段代码：实现中心点的提取

import geopandas as gpd

def get_all_centroids_with_attributes(input_shp_path, output_shp_path):
    """
    提取输入shapefile文件中所有有效多边形的中心点，并将“TBBM”列数据作为属性保存到新的shapefile中。

    参数:
    input_shp_path (str): 输入shapefile文件的路径。
    output_shp_path (str): 输出shapefile文件的路径。
    """
    # 读取shapefile文件
    gdf = gpd.read_file(input_shp_path)

    # 过滤掉几何对象为None的记录
    gdf_with_geometry = gdf.dropna(subset=['geometry'])

    # 计算所有有效多边形的中心点
    centroids = gdf_with_geometry['geometry'].head(10).apply(lambda x: x.centroid)

    # 为了将“TBBM”列数据与中心点关联，我们需要创建一个包含中心点和对应“TBBM”值的新字典
    centroids_with_attributes = [{"geometry": centroid, "TBBM": gdf_with_geometry.loc[gdf_with_geometry.index, 'TBBM'].iloc[i]} for i, centroid in enumerate(centroids)]

    # 创建一个新的GeoDataFrame来存储所有中心点和属性
    centroids_gdf = gpd.GeoDataFrame(centroids_with_attributes, crs=gdf.crs)

    # 输出所有中心点和属性到新的shapefile
    centroids_gdf.to_file(output_shp_path)

# 使用函数
input_shp_path = r'E:\jwztestnew\newshp_5584(1)\newshp_5584.shp'
output_shp_path = r'E:\jwztestnew\center\newshp_5584_10_centroids_with_bttm.shp'
get_all_centroids_with_attributes(input_shp_path, output_shp_path)

第二段代码：根据提取到的中心点做缓冲区矩形

import geopandas as gpd
import rasterio
import os
from rasterio.windows import Window
from rasterio.transform import Affine


def crop_tif_around_centroid_for_each_tbbm(tif_dir, shp_path, output_dir, crop_size=512):
    # 读取shapefile文件以获取中心点坐标和TBBM列的值
    gdf = gpd.read_file(shp_path)

    # 创建输出目录，如果不存在
    if not os.path.exists(output_dir):
        os.makedirs(output_dir)

    # 遍历shapefile中的每个点
    for index, row in gdf.iterrows():
        tbbm = row['TBBM']
        centroid = row['geometry'].centroid
        print(f'处理TBBM值为 {tbbm} 的点 {index + 1}。')

        # 在tif目录下查找名称与TBBM值一致的tif文件
        tif_files = [f for f in os.listdir(tif_dir) if f.startswith(tbbm) and f.endswith('.tif')]
        if not tif_files:
            print(f'未找到TBBM值为 {tbbm} 的tif文件。')
            continue

        # 处理每个找到的tif文件
        for tif_file in tif_files:
            print(f'开始处理文件：{tif_file}。')
            tif_path = os.path.join(tif_dir, tif_file)
            output_tif_path = os.path.join(output_dir, tif_file)

            with rasterio.open(tif_path) as src:
                # 获取tif文件的地理变换参数
                transform = src.transform

                # 将中心点的地理坐标转换为像素坐标
                y_center, x_center = rasterio.transform.rowcol(transform, centroid.x, centroid.y)

                # 计算裁剪窗口的左上角和右下角像素坐标
                x_min = x_center - (crop_size / 2)
                y_min = y_center - (crop_size / 2)
                x_max = x_center + (crop_size / 2)
                y_max = y_center + (crop_size / 2)

                # 确保裁剪窗口不超出tif图像的边界
                x_min = max(0, min(x_min, src.width - 1))
                y_min = max(0, min(y_min, src.height - 1))
                x_max = min(src.width, max(x_max, 0))
                y_max = min(src.height, max(y_max, 0))
                width = x_max - x_min
                height = y_max - y_min

                # 创建裁剪窗口
                window = Window(x_min, y_min, width, height)

                # 读取裁剪窗口内的图像数据
                out_img = src.read(window=window)

                # 新的仿射变换系数
                coord_ltx, coord_lty = rasterio.transform.xy(transform, y_min, x_min)
                window_affine = Affine(
                    src.transform[0],  # X 像素尺寸
                    src.transform[1],  # X 方向旋转
                    coord_ltx,  # X 偏移量
                    src.transform[3],  # Y 方向旋转
                    src.transform[4],  # Y 像素尺寸
                    coord_lty  # Y 偏移量
                )
                new_transform = window_affine

                # 更新元数据
                out_meta = src.meta.copy()
                out_meta.update({
                    "height": y_max - y_min,
                    "width": x_max - x_min,
                    "transform": new_transform
                })

            # 保存裁剪的tif文件
            with rasterio.open(output_tif_path, 'w', **out_meta) as dst:
                dst.write(out_img)

            print(f'文件 {tif_file} 处理完成，已保存到 {output_tif_path}。')


# 使用函数
tif_dir = 'E:\\jwztestnew\\20240329_Sample_8bit'
shp_path = 'E:\\jwztestnew\\center\\newshp_5584_all_centroids_with_bttm.shp'
output_dir = 'Z:\\YuXuening\\ProcessData\\20240410_Sample_8bit_512'
crop_tif_around_centroid_for_each_tbbm(tif_dir, shp_path, output_dir)