python实现ENVI功能:对遥感图像进行计算机分类

Author:v0_persistent
地理信息科学专业
把成果分享化为持续学习的动力!

前言

计算机遥感图像分类是通过模式识别理论,利用计算机将遥感图像分为若干地物类别的方法。提取信息类型包括分类、变化监测、物理量的提取、指标提取、特定地物和状态提取等,在多领域有广泛的应用。
本文将使用python进行遥感图像分类。

图像分类原理:依据和方法

统计模式识别的关键是提取待识别模式的一组统计特征值,然后按一定准则作出决策,对图像予以识别。
遥感图像分类的主要依据是地物的光谱特征,即地物电磁波辐射的多波段测量值,测量值可用作遥感图像分类的特征变量。python可以作为一种决策工具,将分类的准则和步骤写成代码令其执行,实现遥感图像分类。

本文将采用监督分类的最小距离分类法实现图像分类。

python实现图像分类

数据准备

已经过大气校正的遥感数据(本文以landset8数据为例)、ROI文件

获取训练区样本数据

(若已有ROI文件跳过此步)
在envi中,对已经进行过大气校正的遥感图像进行感兴趣区的选择,将ROI导出为dat文件
感兴趣区提取
将ROI保存为特征图层
不同类别地物对应不同的value值,此处用于后面的python分类。
在这里插入图片描述

进行最小距离分类

工具包的准备和数据引入

用到rasterio、matplotlib、numpy工具包,引入原遥感图像的7个波段数据和ROI图像7波段数据。

import rasterio as rio
import matplotlib.pyplot as plt
import numpy as np
import matplotlib
matplotlib.use('TkAgg')

# 原始图像7波段数据
data = rio.open('D:\\envi_practice\\遥感制图\\flaash_result\\image_rc_ac.dat')

band1 = data.read(1)
band2 = data.read(2)
band3 = data.read(3)
band4 = data.read(4)
band5 = data.read(5)
band6 = data.read(6)
band7 = data.read(7)

# 感兴趣区7波段数据
train_data = rio.open('D:\\envi_practice\\遥感制图\\tx_fenlei\\train_roi\\train_roi.dat')

# train_roi = train_data.read(1)
# plt.figure()
# plt.imshow(train_roi, cmap='summer', vmin=0, vmax=20)
# plt.show()

提取相应地物的栅格

ROI文件中不同的地物已经对应了不同的value值,此处借助ROI来提取原遥感图像的栅格。

# 水体
b1_water_mask = train_data == 3
# plt.imshow(b1_water_mask, cmap='grey')
b1_water = band1[train_roi == 3]
b2_water = band2[train_roi == 3]
b3_water = band3[train_roi == 3]
b4_water = band4[train_roi == 3]
b5_water = band5[train_roi == 3]
b6_water = band6[train_roi == 3]
b7_water = band7[train_roi == 3]

# 垂直堆叠 vertical stack
# vstack接受一个数组序列
train_water = np.vstack((b1_water, b2_water, b3_water, b4_water, b5_water, b6_water, b7_water))
train_water_center = train_water.mean(axis=1)

# 以此类推

# 建筑
b1_building_mask = train_data == 1
# plt.imshow(b1_building_mask, cmap='grey')
b1_building = band1[train_roi == 1]
b2_building = band2[train_roi == 1]
b3_building = band3[train_roi == 1]
b4_building = band4[train_roi == 1]
b5_building = band5[train_roi == 1]
b6_building = band6[train_roi == 1]
b7_building = band7[train_roi == 1]

train_building = np.vstack((b1_building, b2_building, b3_building, b4_building, b5_building, b6_building, b7_building))
train_building_center = train_building.mean(axis=1)

# 森林
b1_forest_mask = train_data == 2
# plt.imshow(b1_forest_mask, cmap='grey')
b1_forest = band2[train_roi == 2]
b2_forest = band2[train_roi == 2]
b3_forest = band3[train_roi == 2]
b4_forest = band4[train_roi == 2]
b5_forest = band5[train_roi == 2]
b6_forest = band6[train_roi == 2]
b7_forest = band7[train_roi == 2]

train_forest = np.vstack((b1_forest, b2_forest, b3_forest, b4_forest, b5_forest, b6_forest, b7_forest))
train_forest_center = train_forest.mean(axis=1)

# 渔场
b1_fishing_mask = train_data == 4
# plt.imshow(b1_fishing_mask, cmap='grey')
b1_fishing = band1[train_roi == 4]
b2_fishing = band2[train_roi == 4]
b3_fishing = band3[train_roi == 4]
b4_fishing = band4[train_roi == 4]
b5_fishing = band5[train_roi == 4]
b6_fishing = band6[train_roi == 4]
b7_fishing = band7[train_roi == 4]

train_fishing = np.vstack((b1_fishing, b2_fishing, b3_fishing, b4_fishing, b5_fishing, b6_fishing, b7_fishing))
train_fishing_center = train_fishing.mean(axis=1)

# 农地
b1_farmland_mask = train_data == 5
# plt.imshow(b1_farmland_mask, cmap='grey')
b1_farmland = band1[train_roi == 5]
b2_farmland = band2[train_roi == 5]
b3_farmland = band3[train_roi == 5]
b4_farmland = band4[train_roi == 5]
b5_farmland = band5[train_roi == 5]
b6_farmland = band6[train_roi == 5]
b7_farmland = band7[train_roi == 5]

train_farmland = np.vstack((b1_farmland, b2_farmland, b3_farmland, b4_farmland, b5_farmland, b6_farmland, b7_farmland))
train_farmland_center = train_farmland.mean(axis=1)

# 展示各类地物的光谱曲线
plt.plot(train_water_center, c='blue', label='water')
plt.plot(train_building_center, c='gray', label='building')
plt.plot(train_forest_center, c='green', label='forest')
plt.plot(train_fishing_center, label='fishing')
plt.plot(train_farmland_center, c='coral', label='farmland')

plt.title('the center of all class')
plt.legend()
plt.show()

各类地物光谱曲线如图所示
在这里插入图片描述

逐个像元分类生成分类图

row_sum = band1.shape[0]
col_sum = band1.shape[1]
# 创建一个与原遥感图像大小相同的空白数组,用于生成分类图像
image_class = np.zeros((row_sum, col_sum))
# 对逐个像元进行分类
for row in range(band1.shape[0]):
    for col in range(band1.shape[1]):
        # 取出每个像元的光谱
        r = np.array([band1[row, col], band2[row, col], band3[row, col], band4[row, col],
                      band5[row, col], band6[row, col], band7[row, col]])
        # 欧式距离求和
        dist_water = np.sum((r - train_water_center) ** 2)
        dist_building = np.sum((r - train_building_center) ** 2)
        dist_forest = np.sum((r - train_forest_center) ** 2)
        dist_fishing = np.sum((r - train_fishing_center) ** 2)
        dist_farmland = np.sum((r - train_farmland_center) ** 2)

        # 距离数组
        distance = np.array([dist_water, dist_building, dist_forest, dist_fishing, dist_farmland])
        # 找出数组中最小值出现的位置
        min_distance = distance.argmin()
        # 为像元分类
        image_class[row, col] = min_distance + 1

plt.figure()
plt.imshow(image_class, cmap='jet')
plt.show()

最终成果如图
在这里插入图片描述

点赞关注支持 -v- ❤️

若有错误,虚心请教

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值