遥感影像语义分割过程记录
遥感影像语义分割
步骤包括前期原始影像，提取样本，制作标签，转换成tensor数据集，选择模型（deeplab），数据喂入训练，新的影像预测。
原始影像
样本获取
主要采用gdal库，在原始影像中随机生成坐标点采用固定尺寸进行裁剪。代码如下：
import gdal
import random
def clip_raster(in_put,out_put,xsize=300,ysize=300):
# 读取要切的原图
in_ds = gdal.Open(in_put)
# print("open tif file succeed")
# 读取原图中的每个波段
in_band1 = in_ds.GetRasterBand(1)
in_band2 = in_ds.GetRasterBand(2)
in_band3 = in_ds.GetRasterBand(3)
im_width = in_ds.RasterXSize # 栅格矩阵的列数
im_height = in_ds.RasterYSize # 栅格矩阵的行数
# 定义切图的起始点坐标(相比原点的横坐标和纵坐标偏移量)
offset_x = random.randint(0,im_width-xsize) # 这里是随机取点，可根据自己的实际需要设置
offset_y = random.randint(0,im_height-ysize)
# 定义切图的大小（矩形框）
block_xsize = xsize # 行
block_ysize = ysize # 列
# 从每个波段中切需要的矩形框内的数据(注意读取的矩形框不能超过原图大小)
out_band1 = in_band1.ReadAsArray(offset_x, offset_y, block_xsize, block_ysize)
out_band2 = in_band2.ReadAsArray(offset_x, offset_y, block_xsize, block_ysize)
out_band3 = in_band3.ReadAsArray(offset_x, offset_y, block_xsize, block_ysize)
# 获取Tif的驱动，为创建切出来的图文件做准备
gtif_driver = gdal.GetDriverByName("GTiff")
# 创建切出来的要存的文件（3代表3个不都按，最后一个参数为数据类型，跟原文件一致）
out_ds = gtif_driver.Create(out_put, block_xsize, block_ysize, 3, in_band1.DataType)
# print("create new tif file succeed")
# 获取原图的原点坐标信息
ori_transform = in_ds.GetGeoTransform()
if ori_transform:
print (ori_transform)
print("Origin = ({}, {})".format(ori_transform[0], ori_transform[3]))
print("Pixel Size = ({}, {})".format(ori_transform[1], ori_transform[5]))
# 读取原图仿射变换参数值
top_left_x = ori_transform[0] # 左上角x坐标
w_e_pixel_resolution = ori_transform[1] # 东西方向像素分辨率
top_left_y = ori_transform[3] # 左上角y坐标
n_s_pixel_resolution = ori_transform[5] # 南北方向像素分辨率
# 根据反射变换参数计算新图的原点坐标
top_left_x = top_left_x + offset_x * w_e_pixel_resolution
top_left_y = top_left_y + offset_y * n_s_pixel_resolution
# 将计算后的值组装为一个元组，以方便设置
dst_transform = (top_left_x, ori_transform[1], ori_transform[2], top_left_y, ori_transform[4], ori_transform[5])
# 设置裁剪出来图的原点坐标
out_ds.SetGeoTransform(dst_transform)
# 设置SRS属性（投影信息）
out_ds.SetProjection(in_ds.GetProjection())
# 写入目标文件
out_ds.GetRasterBand(1).WriteArray(out_band1)
out_ds.GetRasterBand(2).WriteArray(out_band2)
out_ds.GetRasterBand(3).WriteArray(out_band3)
# 将缓存写入磁盘
out_ds.FlushCache()
print("FlushCache succeed")
# 计算统计值
# for i in range(1, 3):
# out_ds.GetRasterBand(i).ComputeStatistics(False)
# print("ComputeStatistics succeed")
if __name__ == "__main__":
# 生成100张图
in_put = "。。。。。\\TL_2020_10M_02_Clip1.tif"
for i in range(100):
out_put = '。。。。。\\clip'+str(i)+".tif"
clip_raster(in_put,out_put)
print("End!")将所有随机裁剪的结果加载至arcgis中如下：
样本制作
1、labelme 制作样本
labelme是很常见的用来制作样本的工具，但是目前不能处理.tif格式的影像，如果使用这个工具需要将图像转jpg或者png，我们希望不通过转换格式来制作样本，因此采用下面方法。
2、arcgis 制作样本
简单来说就是创建一个shp文件，画出我们分类的实体，并在属性表中标记类别，然后在工具箱里的【Conversion Tools】->【To Raster】->【Feature to Raster】将shp转换为栅格文件，特别需要注意的是，需要在环境设置里的处理范围设置为与原影像一致。
主要参考： https://zhuanlan.zhihu.com/p/163353715?utm_source=qq
样本读入tensorflow
目前仍然存在问题，无法读取tif格式。
原始的方法是通过tf.gfile.FastGFile实现对图片的读取read成二进制格式。如下为其中一部分二进制。
修改后的方法是采用gdal读取数据，然后转换为numpy数组，再用np的方法直接转换为二进制。如下。
但问题是对应的二进制文件却与之前的不同。
因此，先使用tif转png\jpg的方式。
tif转png
1、直接将gdal读取的tif转换成np数组，再用PIL读取。
import numpy as np
from PIL import Image
import gdal
tif = gdal.Open("C:\\Users\\94320\\Desktop\\image\\clip.tif")
print("打开tif：",tif)
# print(label2)
cols = tif.RasterXSize # 栅格矩阵的列数
rows = tif.RasterYSize # 栅格矩阵的行数
im_bands = tif.RasterCount # 波段数
label2 = tif.ReadAsArray(0, 0, cols,rows)
data=np.empty([rows,cols,3],dtype = float)
###########################
# 直接将gdal读取的tif转换成np数组，再用PIL读取。
###########################
for i in range(im_bands):
band=tif.GetRasterBand(i+1)
data1 = band.ReadAsArray()
data1 = data1
print("data1", data1)
data[:, :, i]=data1
print(Image.fromarray(np.uint8(data))) # 将16位的数值直接用np.uint8转换成8位的
Image.fromarray(np.uint8(data)).save("C:\\Users\\94320\\Desktop\\image\\tiff.png")得到的结果如图所示：显示不符合实际情况。
2、以百分比截断进行拉伸
import numpy as np
from PIL import Image
import gdal
def tif2png(input_tif,output_dir):
cols = input_tif.RasterXSize # 栅格矩阵的列数
rows = input_tif.RasterYSize # 栅格矩阵的行数
im_bands = input_tif.RasterCount # 波段数
data=np.empty([rows,cols,3],dtype = float)
out = np.zeros_like(data, dtype=np.uint8) # 创建一个空的图
for i in range(im_bands):
band = input_tif.GetRasterBand(i + 1)
data1 = band.ReadAsArray()
data[:, :, i] = data1
# 这里以百分比截断进行拉伸
a = 0
b = 255
c = np.percentile(data[:, :,i], 0.7)
d = np.percentile(data[:, :,i], 99.9)
t = a + (data[:, :,i] - 1) * (b - a) / (d - c)
t[t < a] = a
t[t > b] = b
out[:, :,i] = t
outImg=Image.fromarray(np.uint8(out))
outImg.save(output_dir)
tif = gdal.Open("C:\\Users\\94320\\Desktop\\image\\clip.tif")
output_dir = "C:\\Users\\94320\\Desktop\\image\\tiff1.png"
print("打开tif：",tif)
tif2png(tif,output_dir)结果如下：与原tif直观差异不大。
稍作修整，批量处理
import numpy as np
from PIL import Image
import gdal
import os
def tif2png(input_dir,output_dir):
input_tif = gdal.Open(input_dir)
cols = input_tif.RasterXSize # 栅格矩阵的列数
rows = input_tif.RasterYSize # 栅格矩阵的行数
im_bands = input_tif.RasterCount # 波段数
data=np.empty([rows,cols,3],dtype = float)
out = np.zeros_like(data, dtype=np.uint8) # 创建一个空的图
for i in range(im_bands):
band = input_tif.GetRasterBand(i + 1)
data1 = band.ReadAsArray()
data[:, :, i] = data1
# 这里以百分比截断进行拉伸
a = 0
b = 255
c = np.percentile(data[:, :,i], 0.7)
d = np.percentile(data[:, :,i], 99.9)
t = a + (data[:, :,i] - 1) * (b - a) / (d - c)
t[t < a] = a
t[t > b] = b
out[:, :,i] = t
outImg=Image.fromarray(np.uint8(out))
outImg.save(output_dir)
if __name__ == "__main__":
tif_path = "C:\\Users\\94320\\Desktop\\image\\" # 输入tif文件夹路径
output_dir = "C:\\Users\\94320\\Desktop\\image_jpg\\" # 输出png文件夹路径
tif = os.listdir(tif_path)
for i, img in enumerate(tif):
img_name = os.path.join(tif_path, img)
save_file = os.path.join(output_dir, img.strip('.tif') + '.png')
tif2png(img_name,save_file)
本文章为转载内容，我们尊重原作者对文章享有的著作权。如有内容错误或侵权问题，欢迎原作者联系我们进行内容更正或删除文章。
