使用Shapefile提取多光谱栅格像素值存为DataFrame时遇索引错误

错误原因 #

• rasterio.mask.mask()返回的是长度为2的元组：(masked_array, out_transform)，仅第一个元素是掩码后的栅格数组（形状为(波段数, 行数, 列数)），你错误地尝试取元组的第2、3个位置（[1]、[2]等），这些位置不存在，直接触发索引越界错误。
• 重复调用rasterio.mask.mask()完全冗余，会重复执行掩码计算，大幅降低效率。

修正后的代码 #

import rasterio
import geopandas as gpd
import pandas as pd
import numpy as np

# 加载掩码矢量文件
mask = gpd.read_file('/disco02/emurillo/Documentos/metodologia_usuev/calval_v2g.shp')

# 加载栅格并一次性完成掩码提取
with rasterio.open('rep_pncm_2018_lluvia_C2_b2b7_int16.tif') as src:
    # 执行掩码操作，获取所有波段的掩码后数组
    masked_data, _ = rasterio.mask.mask(src, mask.geometry, crop=True, nodata=src.nodata)
    # 获取NoData值用于过滤无效像素
    nodata_val = src.nodata

# 将每个波段的二维数组展平为一维，同时过滤NoData值
band_data = {}
for i in range(masked_data.shape[0]):
    band_name = f'band{i+1}'
    # 展平并过滤无效值
    flattened = masked_data[i].flatten()
    valid_vals = flattened[flattened != nodata_val]
    band_data[band_name] = valid_vals

# 构建DataFrame
df = pd.DataFrame(band_data)

# 合并Shapefile属性：单个多边形直接赋值；多多边形需额外处理像素归属
if len(mask) == 1:
    for col in mask.columns:
        df[col] = mask[col].values[0]
else:
    print("注意：多多边形场景需匹配每个像素对应的多边形属性，当前代码未处理此逻辑")

# 查看结果示例
print(df.head())

关键说明 #

1. 一次性执行掩码操作，避免重复计算，同时正确获取掩码后的数组（masked_data）。
2. 展平数组时过滤NoData值，避免无效像素进入DataFrame。
3. 多多边形场景下，需要额外处理每个像素对应的多边形属性（比如用rasterio.features.rasterize标记像素归属），当前代码默认处理单个多边形的情况。

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