多波段遥感影像和单波段影像裁剪代码,保留坐标信息

import os
from osgeo import gdal
import numpy as np

# 读取tif数据集
def readTif(fileName):
    dataset = gdal.Open(fileName)
    if dataset == None:
        print(fileName + "文件无法打开")
    return dataset

# 保存tif文件函数
def writeTiff(im_data, im_geotrans, im_proj, path):
    if 'int8' in im_data.dtype.name:
        datatype = gdal.GDT_Byte
    elif 'int16' in im_data.dtype.name:
        datatype = gdal.GDT_UInt16
    else:
        datatype = gdal.GDT_Float32

    if len(im_data.shape) == 3:
        im_bands, im_height, im_width = im_data.shape
    elif len(im_data.shape) == 2:
        im_data = np.array([im_data])
        im_bands, im_height, im_width = im_data.shape

    # 创建文件
    driver = gdal.GetDriverByName("GTiff")
    dataset = driver.Create(path, int(im_width), int(im_height), int(im_bands), datatype)
    if (dataset != None):
        dataset.SetGeoTransform(im_geotrans)  # 写入仿射变换参数
        dataset.SetProjection(im_proj)  # 写入投影
    for i in range(im_bands):
        dataset.GetRasterBand(i + 1).WriteArray(im_data[i])
    del dataset

# 滑动窗口裁剪函数
def TifCrop(TifDir, SavePath, CropSize):
    for filename in os.listdir(TifDir):
        if filename.endswith(".tif"):
            TifPath = os.path.join(TifDir, filename)
            dataset_img = readTif(TifPath)
            width = dataset_img.RasterXSize
            height = dataset_img.RasterYSize
            proj = dataset_img.GetProjection()
            geotrans = dataset_img.GetGeoTransform()
            img = dataset_img.ReadAsArray(0, 0, width, height)  # 获取数据

            # 获取当前文件夹的文件个数len,并以len+1命名即将裁剪得到的图像
            new_name = len(os.listdir(SavePath))

            # 对图像进行填充
            pad_height = CropSize - (height % CropSize)
            pad_width = CropSize - (width % CropSize)
            img = np.pad(img, ((0, 0), (0, pad_height), (0, pad_width)), mode='constant')

            # 裁剪图片
            for i in range(int(height / CropSize)):
                for j in range(int(width / CropSize)):
                    # 如果图像是单波段
                    if (len(img.shape) == 2):
                        cropped = img[i * CropSize: (i + 1) * CropSize, j * CropSize: (j + 1) * CropSize]
                    # 如果图像是多波段
                    else:
                        cropped = img[:, i * CropSize: (i + 1) * CropSize, j * CropSize: (j + 1) * CropSize]
                    # 写图像
                    writeTiff(cropped, geotrans, proj, SavePath + "/%d.tif" % new_name)
                    # 文件名 + 1
                    new_name = new_name + 1

TifCrop(r"文件夹地址/",
        r"裁剪存放位置/", 256)

上面的裁剪代码是用于裁剪标签和多波段影像的代码,只需要改一下这行代码中img = np.pad(img, ((0, 0), (0, pad_height), (0, pad_width)), mode='constant')的(0,0),有(0,0)代表多波段,无则反之。

保留坐标信息的代码是这一行: writeTiff(cropped, geotrans, proj, SavePath + "/%d.tif" % new_name)。

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值