使用Python计算Sentinel-2影像NDVI结果全为0的问题排查

核心错误：近红外波段变量赋值错误 #

在读取近红外波段（nir）的代码块里，你误将nir = nir.astype(np.float64)写成了nir = red.astype(np.float64)，这导致nir数组完全等于red数组，计算nir - red时结果全为0，最终NDVI自然全部是0。

次要错误：文件路径未加引号 #

red_f和nir_f的路径是字符串，需要用引号包裹，且要带上文件后缀（比如.tif），否则Python无法识别为合法的路径字符串。

修正后的完整代码 #

import os
import numpy as np
import rasterio as rio

# 抑制numpy的除零、无效值警告
np.seterr(divide='ignore', invalid='ignore')

# 修正：路径添加引号，补充文件后缀
red_f = "/absolute/path/to/band/4.tif"
nir_f = "/absolute/path/to/band/8.tif"

def calc_ndvi():
    with rio.open(red_f) as src:
        red = src.read()
        red = red.astype(np.float64)
    
    with rio.open(nir_f) as src:
        nir = src.read()
        # 修正：将red改为nir，正确转换近红外波段的数据类型
        nir = nir.astype(np.float64)

    ndvi = np.divide((nir - red), (nir + red))
    
    # 可选：将计算结果保存为TIFF文件，方便后续查看
    with rio.open(red_f) as src:
        # 复制原影像的元数据，修改数据类型为float64
        profile = src.profile
        profile.update(dtype=np.float64, count=1)
        with rio.open("ndvi_result.tif", 'w', **profile) as dst:
            dst.write(ndvi)
    
    return ndvi

ndvi = calc_ndvi()

额外注意事项 #

• 确保路径中的文件真实存在，且对应Sentinel-2产品的Band4（红波段）和Band8（近红外波段）。
• 如果使用的是Sentinel-2 L1C级产品，建议先做大气校正再计算NDVI，结果会更准确；L2A级产品已做过大气校正，可直接计算。
• 计算完成后，建议用QGIS等GIS软件打开保存的NDVI结果文件，验证数值是否合理（正常NDVI范围为-1到1）。

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