python矩阵堆叠-实现遥感影像波段叠加

2021/5/22 14:55:15

本文主要是介绍python矩阵堆叠-实现遥感影像波段叠加,对大家解决编程问题具有一定的参考价值,需要的程序猿们随着小编来一起学习吧!

前文读取landsat8数据,计算NDVI并保存成tiff文件和python批量读取landsat8的波段分别介绍了对landsat8数据的不同读取方式,以及对利用波段数据进行简单的指数计算。本文在此基础上实现波段的叠加。

准备工作

这一部分参考python批量读取landsat8的波段,利用文中的函数,我们对landsat8的波段路径进行了一个整理,存放到了bands列表中并返回。
看一下格式:第一张影像的前5个波段

from read_landsat8 import read_landsat8_bands
import numpy as np
from osgeo import gdal_array

base_path = 'data'
bands = read_landsat8_bands(base_path)

bands[0][:5]

[‘data\LC08_L1TP_121034_20130402_20200912_02_T1\LC08_L1TP_121034_20130402_20200912_02_T1_B1.TIF’, ‘data\LC08_L1TP_121034_20130402_20200912_02_T1\LC08_L1TP_121034_20130402_20200912_02_T1_B2.TIF’, ‘data\LC08_L1TP_121034_20130402_20200912_02_T1\LC08_L1TP_121034_20130402_20200912_02_T1_B3.TIF’, ‘data\LC08_L1TP_121034_20130402_20200912_02_T1\LC08_L1TP_121034_20130402_20200912_02_T1_B4.TIF’, ‘data\LC08_L1TP_121034_20130402_20200912_02_T1\LC08_L1TP_121034_20130402_20200912_02_T1_B5.TIF’]

读取波段1、2、3:

# 读取波段
B1 = bands[0][0]
B2 = bands[0][1]
B3 = bands[0][2]
B1_gdal = gdal_array.LoadFile(B1)
B2_gdal = gdal_array.LoadFile(B2)
B3_gdal = gdal_array.LoadFile(B3)

print(B1_gdal.shape)
print(B2_gdal.shape)
print(B3_gdal.shape)  # 维度相等

(7301, 7341)
(7301, 7341)
(7301, 7341)

波段叠加:

本文只是简单演示,这里我们只对前面读取的123波段进行叠加
先用numpy转化成ndarray形式

# 转化成ndarray形式
B1_np = np.array(B1_gdal)
B2_np = np.array(B2_gdal)
B3_np = np.array(B3_gdal)
print(B1_np.shape)
print(B2_np.shape)
print(B3_np.shape)

维度仍保持

(7301, 7341)
(7301, 7341)
(7301, 7341)

用构建数组的方式叠加

B123 = np.array([B1_np, B2_np, B3_np])
print(B123.shape)

输出可以看成是channel_first形式

(3, 7301, 7341)

用stack方式叠加

B123 = np.stack([B1_np, B2_np, B3_np])
print(B123.shape)

输出维度与上面相同

(3, 7301, 7341)

stack叠加可以按照不同的轴实现不同的维度叠加,默认axis=0

B123 = np.stack([B1_np, B2_np, B3_np], axis = 0)
print(B123.shape)

(3, 7301, 7341)

B123 = np.stack([B1_np, B2_np, B3_np], axis = 1)
print(B123.shape)

(7301, 3, 7341)

B123 = np.stack([B1_np, B2_np, B3_np], axis = 2)
print(B123.shape)

(7301, 7341,3)

用concatenate方式叠加

B123 = np.concatenate([B1_np, B2_np, B3_np])
print(B123.shape)

可以看到,它是按照行叠加了,这种叠加方式也有一些用处,不过在这里不满足我们的要求。

(21903, 7341)



这篇关于python矩阵堆叠-实现遥感影像波段叠加的文章就介绍到这儿,希望我们推荐的文章对大家有所帮助,也希望大家多多支持为之网!


扫一扫关注最新编程教程