python遥感图像裁剪成深度学习样本

您所在的位置:网站首页 python裁剪 python遥感图像裁剪成深度学习样本

python遥感图像裁剪成深度学习样本

2023-04-01 07:35| 来源: 网络整理| 查看: 265

前言

如果将图像直接输入到深度学习网络中,会导致内存溢出,因此需要将图像裁剪成图像块输入到网络中。裁剪方法包括规则格网裁剪滑动窗口裁剪以及随机裁剪

规则格网裁剪滑动窗口裁剪随机裁剪正文

规则格网裁剪属于重复率为0的滑动窗口裁剪,滑动窗口裁剪代码为:

import os 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 ''' 滑动窗口裁剪函数 TifPath 影像路径 SavePath 裁剪后保存目录 CropSize 裁剪尺寸 RepetitionRate 重复率 ''' def TifCrop(TifPath, SavePath, CropSize, RepetitionRate): 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)) + 1 # 裁剪图片,重复率为RepetitionRate for i in range(int((height - CropSize * RepetitionRate) / (CropSize * (1 - RepetitionRate)))): for j in range(int((width - CropSize * RepetitionRate) / (CropSize * (1 - RepetitionRate)))): # 如果图像是单波段 if(len(img.shape) == 2): cropped = img[int(i * CropSize * (1 - RepetitionRate)) : int(i * CropSize * (1 - RepetitionRate)) + CropSize, int(j * CropSize * (1 - RepetitionRate)) : int(j * CropSize * (1 - RepetitionRate)) + CropSize] # 如果图像是多波段 else: cropped = img[:, int(i * CropSize * (1 - RepetitionRate)) : int(i * CropSize * (1 - RepetitionRate)) + CropSize, int(j * CropSize * (1 - RepetitionRate)) : int(j * CropSize * (1 - RepetitionRate)) + CropSize] # 写图像 writeTiff(cropped, geotrans, proj, SavePath + "/%d.tif"%new_name) # 文件名 + 1 new_name = new_name + 1 # 向前裁剪最后一列 for i in range(int((height-CropSize*RepetitionRate)/(CropSize*(1-RepetitionRate)))): if(len(img.shape) == 2): cropped = img[int(i * CropSize * (1 - RepetitionRate)) : int(i * CropSize * (1 - RepetitionRate)) + CropSize, (width - CropSize) : width] else: cropped = img[:, int(i * CropSize * (1 - RepetitionRate)) : int(i * CropSize * (1 - RepetitionRate)) + CropSize, (width - CropSize) : width] # 写图像 writeTiff(cropped, geotrans, proj, SavePath + "/%d.tif"%new_name) new_name = new_name + 1 # 向前裁剪最后一行 for j in range(int((width - CropSize * RepetitionRate) / (CropSize * (1 - RepetitionRate)))): if(len(img.shape) == 2): cropped = img[(height - CropSize) : height, int(j * CropSize * (1 - RepetitionRate)) : int(j * CropSize * (1 - RepetitionRate)) + CropSize] else: cropped = img[:, (height - CropSize) : height, int(j * CropSize * (1 - RepetitionRate)) : int(j * CropSize * (1 - RepetitionRate)) + CropSize] writeTiff(cropped, geotrans, proj, SavePath + "/%d.tif"%new_name) # 文件名 + 1 new_name = new_name + 1 # 裁剪右下角 if(len(img.shape) == 2): cropped = img[(height - CropSize) : height, (width - CropSize) : width] else: cropped = img[:, (height - CropSize) : height, (width - CropSize) : width] writeTiff(cropped, geotrans, proj, SavePath + "/%d.tif"%new_name) new_name = new_name + 1 # 将影像1裁剪为重复率为0.1的256×256的数据集 TifCrop(r"Data\data2\tif\data2.tif", r"Data\train\image1", 256, 0.1) TifCrop(r"Data\data2\label\label.tif", r"data\train\label1", 256, 0.1)

随机裁剪的话,只需要随机生成裁剪图像的左上角坐标,然后以此为基准取特定大小的矩阵块就可以了。代码:

import random import gdal import numpy as np import os # 读取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 ''' 随机裁剪函数 ImagePath 原始影像路径 LabelPath 标签影像路径 IamgeSavePath 原始影像裁剪后保存目录 LabelSavePath 标签影像裁剪后保存目录 CropSize 裁剪尺寸 CutNum 裁剪数量 ''' def RandomCrop(ImagePath, LabelPath, IamgeSavePath, LabelSavePath, CropSize, CutNum): dataset_img = readTif(ImagePath) width = dataset_img.RasterXSize height = dataset_img.RasterYSize proj = dataset_img.GetProjection() geotrans = dataset_img.GetGeoTransform() img = dataset_img.ReadAsArray(0,0,width,height)#获取哟昂数据 dataset_label = readTif(LabelPath) label = dataset_label.ReadAsArray(0,0,width,height)#获取标签数据 # 获取当前文件夹的文件个数len,并以len+1命名即将裁剪得到的图像 fileNum = len(os.listdir(IamgeSavePath)) new_name = fileNum + 1 while(new_name < CutNum + fileNum + 1): # 生成剪切图像的左上角XY坐标 UpperLeftX = random.randint(0, height - CropSize) UpperLeftY = random.randint(0, width - CropSize) if(len(img.shape) == 2): imgCrop = img[UpperLeftX : UpperLeftX + CropSize, UpperLeftY : UpperLeftY + CropSize] else: imgCrop = img[:, UpperLeftX : UpperLeftX + CropSize, UpperLeftY : UpperLeftY + CropSize] if(len(label.shape) == 2): labelCrop = label[UpperLeftX : UpperLeftX + CropSize, UpperLeftY : UpperLeftY + CropSize] else: labelCrop = label[:, UpperLeftX : UpperLeftX + CropSize, UpperLeftY : UpperLeftY + CropSize] writeTiff(imgCrop, geotrans, proj, IamgeSavePath + "/%d.tif"%new_name) writeTiff(labelCrop, geotrans, proj, LabelSavePath + "/%d.tif"%new_name) new_name = new_name + 1 # 裁剪得到300张256*256大小的训练集 RandomCrop(r"Data\data2\tif\data2.tif", r"Data\data2\label\label.tif", r"Data\train\image1", r"Data\train\label1", 256,300)后记

裁剪好了就可以语义分割了

欢迎留言交流~

添加地理坐标

本来如果裁剪只是做样本的话是用不到坐标信息的,因为深度学习只是学习的像素,学习不了坐标。但是好多同学留言加上坐标信息,其实加上坐标信息只需要利用仿射变换矩阵计算裁剪的样本的左上角的坐标,然后更新新的仿射变换矩阵就可以了。代码:

import os try: import gdal except : 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 CoordTransf(Xpixel, Ypixel, GeoTransform): XGeo = GeoTransform[0]+GeoTransform[1]*Xpixel+Ypixel*GeoTransform[2]; YGeo = GeoTransform[3]+GeoTransform[4]*Xpixel+Ypixel*GeoTransform[5]; return XGeo, YGeo ''' 滑动窗口裁剪函数 TifPath 影像路径 SavePath 裁剪后保存目录 CropSize 裁剪尺寸 RepetitionRate 重复率 ''' def TifCrop(TifPath, SavePath, CropSize, RepetitionRate): if not os.path.exists(SavePath): os.makedirs(SavePath) 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)) + 1 # 裁剪图片,重复率为RepetitionRate for i in range(int((height - CropSize * RepetitionRate) / (CropSize * (1 - RepetitionRate)))): for j in range(int((width - CropSize * RepetitionRate) / (CropSize * (1 - RepetitionRate)))): # 如果图像是单波段 if(len(img.shape) == 2): cropped = img[int(i * CropSize * (1 - RepetitionRate)) : int(i * CropSize * (1 - RepetitionRate)) + CropSize, int(j * CropSize * (1 - RepetitionRate)) : int(j * CropSize * (1 - RepetitionRate)) + CropSize] # 如果图像是多波段 else: cropped = img[:, int(i * CropSize * (1 - RepetitionRate)) : int(i * CropSize * (1 - RepetitionRate)) + CropSize, int(j * CropSize * (1 - RepetitionRate)) : int(j * CropSize * (1 - RepetitionRate)) + CropSize] XGeo, YGeo = CoordTransf(int(j * CropSize * (1 - RepetitionRate)), int(i * CropSize * (1 - RepetitionRate)), geotrans) crop_geotrans = (XGeo, geotrans[1], geotrans[2], YGeo, geotrans[4], geotrans[5]) # 写图像 writeTiff(cropped, crop_geotrans, proj, SavePath + "/%d.tif"%new_name) # 文件名 + 1 new_name = new_name + 1 # 向前裁剪最后一列 for i in range(int((height-CropSize*RepetitionRate)/(CropSize*(1-RepetitionRate)))): if(len(img.shape) == 2): cropped = img[int(i * CropSize * (1 - RepetitionRate)) : int(i * CropSize * (1 - RepetitionRate)) + CropSize, (width - CropSize) : width] else: cropped = img[:, int(i * CropSize * (1 - RepetitionRate)) : int(i * CropSize * (1 - RepetitionRate)) + CropSize, (width - CropSize) : width] XGeo, YGeo = CoordTransf(width - CropSize, int(i * CropSize * (1 - RepetitionRate)), geotrans) crop_geotrans = (XGeo, geotrans[1], geotrans[2], YGeo, geotrans[4], geotrans[5]) # 写图像 writeTiff(cropped, crop_geotrans, proj, SavePath + "/%d.tif"%new_name) new_name = new_name + 1 # 向前裁剪最后一行 for j in range(int((width - CropSize * RepetitionRate) / (CropSize * (1 - RepetitionRate)))): if(len(img.shape) == 2): cropped = img[(height - CropSize) : height, int(j * CropSize * (1 - RepetitionRate)) : int(j * CropSize * (1 - RepetitionRate)) + CropSize] else: cropped = img[:, (height - CropSize) : height, int(j * CropSize * (1 - RepetitionRate)) : int(j * CropSize * (1 - RepetitionRate)) + CropSize] XGeo, YGeo = CoordTransf(int(j * CropSize * (1 - RepetitionRate)), height - CropSize, geotrans) crop_geotrans = (XGeo, geotrans[1], geotrans[2], YGeo, geotrans[4], geotrans[5]) writeTiff(cropped, crop_geotrans, proj, SavePath + "/%d.tif"%new_name) # 文件名 + 1 new_name = new_name + 1 # 裁剪右下角 if(len(img.shape) == 2): cropped = img[(height - CropSize) : height, (width - CropSize) : width] else: cropped = img[:, (height - CropSize) : height, (width - CropSize) : width] XGeo, YGeo = CoordTransf(width - CropSize, height - CropSize, geotrans) crop_geotrans = (XGeo, geotrans[1], geotrans[2], YGeo, geotrans[4], geotrans[5]) writeTiff(cropped, crop_geotrans, proj, SavePath + "/%d.tif"%new_name) new_name = new_name + 1引用:

Wang, Z.; Zhou, Y.; Wang, S.; Wang, F.; Xu, Z. House building extraction from high resolution remote sensing image based on IEU-Net. J. Remote Sens. 2021, in press. [CrossRef]



【本文地址】


今日新闻


推荐新闻


CopyRight 2018-2019 办公设备维修网 版权所有 豫ICP备15022753号-3