Snappy读取TIFF影像波段替代方案及rasterio报错解决

• 触发IndexError的直接原因是src.read()参数传值错误：rasterio的read()方法接收从1开始计数的整数波段索引作为参数，无参数传入时默认读取全部波段。你传入了文件路径字符串S1_S2_stack，方法会逐字符解析字符串作为索引，第一个字符S不在合法的1~31波段范围内，因此抛出索引越界错误。
• 后续src1.tags()的调用逻辑也存在错误：read()返回的是numpy数组类型，不存在tags()方法，波段元数据、标签存储在rasterio打开的数据集对象src中，而非读取后的数值数组里。

方案1：rasterio实现（推荐） #

用上下文管理器自动管理文件开闭，避免资源泄漏，直接读取全部31个波段并提取波段信息：

import rasterio
import numpy as np

S1_S2_stack = 'S1_S2_stack.tif'

with rasterio.open(S1_S2_stack) as src:
    # 读取所有波段，返回数组维度为(波段数, 影像高度, 影像宽度)
    img_stack = src.read()
    print(f"影像维度校验：{img_stack.shape}，第一个值应为31")

    # 提取所有波段的元信息
    bands_info = []
    for band_id in range(1, src.count + 1):
        single_band_info = {
            "band_id": band_id,
            "band_desc": src.descriptions[band_id-1],
            "data_type": src.dtypes[band_id-1],
            "no_data_value": src.nodata,
            "band_tags": src.tags(band_id)
        }
        bands_info.append(single_band_info)

# 按需打印波段信息
for info in bands_info:
    print(info)

如果需要读取单个指定波段，给read()传入对应整数索引即可，比如读取第10个波段就写src.read(10)。

方案2：GDAL实现（备选） #

如果习惯GDAL接口，可使用以下代码实现相同功能：

from osgeo import gdal
import numpy as np

S1_S2_stack = 'S1_S2_stack.tif'
# 只读模式打开影像
ds = gdal.Open(S1_S2_stack, gdal.GA_ReadOnly)
total_bands = ds.RasterCount
print(f"总波段数：{total_bands}")

img_stack = []
bands_info = []
for band_id in range(1, total_bands + 1):
    band = ds.GetRasterBand(band_id)
    # 读取单波段数组
    band_arr = band.ReadAsArray()
    img_stack.append(band_arr)
    # 收集波段信息
    single_band_info = {
        "band_id": band_id,
        "band_desc": band.GetDescription(),
        "data_type": gdal.GetDataTypeName(band.DataType),
        "no_data_value": band.GetNoDataValue(),
        "band_tags": band.GetMetadata()
    }
    bands_info.append(single_band_info)

# 转换为(波段数, 高度, 宽度)维度的标准数组
img_stack = np.array(img_stack)
# 释放数据集资源
ds = None

• 原代码中训练、验证矢量路径的注释写反了，training_points对应testing.shp、validation_points对应training.shp，后续提取样本值时注意修正，避免数据集混淆。
• 后续如果需要提取矢量点对应位置的像元值做地物分类，直接用rasterio的sample方法即可完成，无需依赖ESA Snappy模块。

内容的提问来源于stack exchange，提问作者Gulnihal