与计算机视觉中的目标检测不同，遥感图像目标检测的检测范围更大，例如机场内的飞机检测、港口内的舰船检测，根据影像分辨率的不同，检测所需图像数据相比自然图像更多、范围更大，如下图1。以高分系列卫星数据为例，高分二号空间分辨率最高达到0.8m，成像幅宽45km；高分六号空间分辨率为2m，成像幅宽800km。
图1 (目标在城市中)
在给定大范围区域快速对感兴趣目标进行搜索任务时，一般最常见的思路就是“化整为块”，将大幅图像分块为不同的小块，块之间可以不重叠，也可以重叠，块的大小满足目标检测算法的输入要求，例如Faster RCNN的输入大小一般设置为1024*1024，对每个块进行检测，然后将结果进行拼接，输出检测信息及包含目标的图像块，经纬度坐标等等。
我自己采用Python库中有cv2或者PIL等读取图像，能够处理1万像素附近的图像，但是对于较大宽幅图像，尤其是超过几万个像素的图像读取，总是报错(这里还在查找原因，要么是本身cv2和PIL不适合读取大范围图像，要么是读取设置有问题)
下面是在Python下，基于GDAL进行大宽幅图像的分块代码，旨在实现大宽幅图像的快速分块。
一、python下的GDAL安装
第一步添加清华源：sudo add-apt-repository -y ppa:ubuntugis/ubuntugis-unstable
第二步：sudo apt update
第三步：sudo apt upgrade
第四步：sudo apt install gdal-bin python-gdal python3-gdal
如果之前安装过gdal，在第三步更新就可以；如果没有装过，则在第四步进行安装
二、图像分块
待处理的图像大小为43648*53759像素，空间分辨率为4.5m，覆盖面积约10559平方公里。分块的python代码如下(这里只考虑无重叠分块，有重叠后续再更新)：
# *_*coding: utf-8 *_*
# author --liming--
import os
import numpy as np
from osgeo import gdal
# 定义读取和保存图像的类
class GRID:
def load_image(self, filename):
image = gdal.Open(filename)
img_width = image.RasterXSize
img_height = image.RasterYSize
img_geotrans = image.GetGeoTransform()
img_proj = image.GetProjection()
img_data = image.ReadAsArray(0, 0, img_width, img_height)
del image
return img_proj, img_geotrans, img_data
def write_image(self, filename, img_proj, img_geotrans, img_data):
# 判断栅格数据类型
if 'int8' in img_data.dtype.name:
datatype = gdal.GDT_Byte
elif 'int16' in img_data.dtype.name:
datatype = gdal.GDT_UInt16
else:
datatype = gdal.GDT_Float32
# 判断数组维度
if len(img_data.shape) == 3:
img_bands, img_height, img_width = img_data.shape
else:
img_bands, (img_height, img_width) = 1, img_data.shape
# 创建文件
driver = gdal.GetDriverByName('GTiff')
image = driver.Create(filename, img_width, img_height, img_bands, datatype)
image.SetGeoTransform(img_geotrans)
image.SetProjection(img_proj)
if img_bands == 1:
image.GetRasterBand(1).WriteArray(img_data)
else:
for i in range(img_bands):
image.GetRasterBand(i+1).WriteArray(img_data[i])
del image # 删除变量,保留数据
if __name__ == '__main__':
import time
import argparse
parser = argparse.ArgumentParser(description='load remote sensing image and split to patch')
parser.add_argument('--image_path',
default='/media/lm/1E7FBDC6EEE168BC/RS_Dataset/test_image/',
help='remote sensing image path')
parser.add_argument('--patch_size',
default=1000,
help='patch size')
parser.add_argument('--patch_save',
default='/media/lm/1E7FBDC6EEE168BC/RS_Dataset/patch_save/',
help='save path of patch image')
args = parser.parse_args()
print('待处理图像路径为:{}'.format(args.image_path))
print('分块大小为:{}'.format(args.patch_size))
print('分块图像保存路径:{}'.format(args.patch_save))
image_path = args.image_path
image_list = os.listdir(image_path)
#image_list.sort(key=lambda x: int(x[:-4])) # 对文件夹中的图像进行排序
image_num = len(image_list)
t_start = time.time()
for k in range(image_num):
time_start = time.time()
img_name = image_path + image_list[k]
proj, geotrans, data = GRID().load_image(img_name)
# 图像分块
patch_size = args.patch_size
patch_save = args.patch_save
channel, width, height = data.shape
num = 0
for i in range(width//patch_size):
for j in range(height//patch_size):
num += 1
sub_image = data[:, i*patch_size:(i+1)*patch_size, j*patch_size:(j+1)*patch_size]
GRID().write_image(patch_save + '{}.tif'.format(num), proj, geotrans, sub_image)
time_end = time.time()
print('第{}张图像分块完毕, 耗时:{}秒'.format(k+1, round((time_end-time_start), 4)))
t_end = time.time()
print('所有图像处理完毕,耗时:{}秒'.format(round((t_end-t_start), 4)))分块的结果截图如下：
共分为2279个块，总耗时为43.2258秒。每个块的大小为1000*1000，即每个块的覆盖面积为4.5平方公里。
图2
图3
完毕！
