Python读取遥感影像?Rasterio怎么用?

作者： GIS研习社更新时间：2025-12-18 20:00:56 分类：GIS基础理论

为什么你用Rasterio读影像总报错？可能是没搞懂这三层结构

“明明路径没错，代码也没拼写错误，怎么一运行就报‘IndexError: band index out of range’？”——这是我后台收到最多的提问之一。别慌，这不是你的错，而是你还没理解遥感影像在Python眼中的“三维本质”。

我在参与某省耕地遥感监测项目时，曾因没注意波段索引从1开始而非0，导致整批NDVI计算结果异常——整整三天白干。后来才明白，Rasterio的“坑”，其实都是设计哲学。

遥感影像不是照片，是带坐标的“数据立方体”

很多人把GeoTIFF当成普通图片，直接用PIL或cv2去读——这就像拿菜刀解剖精密仪器。遥感影像是一个“数据立方体”：长 × 宽 × 波段数，外加地理坐标、投影信息、元数据等“说明书”。Rasterio就是专门为解析这个复杂结构而生的瑞士军刀。

举个生活化类比：想象你买了一盒多层蛋糕（比如提拉米苏），每一层口味不同（对应不同波段），盒子上还贴着生产日期、产地、食用说明（元数据）。Rasterio不仅能帮你一层层切开蛋糕，还能告诉你每层糖分含量和最佳食用温度。

三步上手Rasterio：打开→读取→关闭，一个都不能少

最基础也最容易出错的操作，就是忘记“上下文管理”。很多教程只教.open()，却不强调必须用with语句——这会导致文件句柄泄露，尤其在批量处理时直接卡死系统。

import rasterio

# ✅ 正确姿势：用 with 自动管理资源
with rasterio.open('sentinel2.tif') as src:
    print(f"影像尺寸: {src.width} x {src.height}")
    print(f"波段数量: {src.count}")
    print(f"坐标系: {src.crs}")
    
    # 读取第1波段（注意：Rasterio波段索引从1开始！）
    band1 = src.read(1)
    
# 文件在这里自动关闭，安全又优雅

这里的关键点：

src.count：返回波段总数，避免你硬编码读第4波段结果影像只有3个波段。
src.read(1)：读第一波段。如果你写read(0)，恭喜你喜提IndexError。
src.crs：获取坐标系，后续做空间分析必备。

实战：批量读取并计算NDVI——避开内存爆炸的陷阱

研究生最爱干的事：一口气读进10GB影像到内存，然后电脑蓝屏。Rasterio的.read()默认读全图，对大影像极其不友好。正确做法是“窗口读取”或“按块处理”。

import rasterio
import numpy as np

with rasterio.open('landsat8.tif') as src:
    # 只读取我们需要的红光和近红外波段（假设是B4和B5）
    red = src.read(4)  # 红光
    nir = src.read(5)  # 近红外
    
    # 计算NDVI，记得处理除零和无效值
    ndvi = (nir.astype(float) - red.astype(float)) / (nir + red + 1e-8)
    
    # 创建输出文件（保持原投影和分辨率）
    profile = src.profile
    profile.update(dtype=rasterio.float32, count=1)
    
    with rasterio.open('ndvi_output.tif', 'w', **profile) as dst:
        dst.write(ndvi.astype(rasterio.float32), 1)

这段代码我优化过三次：

第一次：没加+1e-8，遇到全黑区域直接NaN满天飞。
第二次：忘了更新输出文件的dtype，结果保存成uint8，NDVI值全被截断。
第三次：加上count=1，否则会默认创建和输入一样多的波段，浪费空间。

进阶技巧：用Window切片读取超大影像

当单景影像超过8GB（如国产高分系列），直接.read()会撑爆内存。这时候要用rasterio.windows.Window分块读取。你可以把它想象成“用显微镜一格格观察地图”。

from rasterio.windows import Window

with rasterio.open('huge_image.tif') as src:
    # 定义一个512x512的小窗口，从(1000, 2000)像素位置开始
    window = Window(col_off=1000, row_off=2000, width=512, height=512)
    subset = src.read(1, window=window)  # 只读第一波段的这个小区域
    
    # 后续可对subset做分类、滤波等操作，再写回大图对应位置

这个技巧在做深度学习样本裁剪时特别有用——我带的学生靠它处理了200+景哨兵2号数据，没崩过一次机。