Landsat8单波段影像经纬度转像素行列号越界如何解决

输出行列号超出范围的核心原因是输入坐标和影像的坐标系不匹配：

• 你传入的(92.8, 26.6)是WGS84地理坐标系（EPSG:4326）下的经纬度，单位是度
• 原生Landsat-8影像采用分带UTM投影坐标系，单位是米，src.transform是基于投影坐标生成的仿射变换参数，直接传入经纬度计算，相当于把度当成米做偏移计算，必然得到完全错误的行列值，你得到的负列号、超9万的行号就是这个原因导致的。

坐标转换后再计算行列号，步骤如下：

1. 读取影像时先确认自身坐标系，直接打印src.crs即可查看，你这个经纬度（92.8°E，26.6°N）对应的Landsat8影像一般为UTM 46N带，对应EPSG:32646
2. 将输入的WGS84经纬度转换到影像自身的投影坐标系下
3. 用投影后的坐标代入rowcol计算行列号，额外校验行列号是否在影像尺寸范围内，再读取对应像素值

可直接运行的参考代码：

import numpy as np
import rasterio
from rasterio.warp import transform

# 替换成你的影像路径
with rasterio.open("LC08_L1TP_xxxx_B4.tif") as src:
    # 输入WGS84经纬度
    lons = np.array([92.8])
    lats = np.array([26.6])
    
    # 坐标转换：从WGS84地理坐标系转到影像自身投影坐标系
    xs_proj, ys_proj = transform("EPSG:4326", src.crs, lons, lats)
    
    # 计算行列号
    rows, cols = rasterio.transform.rowcol(src.transform, xs_proj, ys_proj)
    row, col = rows[0], cols[0]
    
    # 校验坐标是否落在影像范围内
    if 0 <= row < src.height and 0 <= col < src.width:
        # 单波段影像读取对应像素值，numpy数组索引直接用[row, col]
        pixel_val = src.read(1)[row, col]
        print(f"numpy数组索引：行{row}，列{col}")
        print(f"对应像素值：{pixel_val}")
    else:
        print("输入坐标不在当前影像覆盖范围内，请检查坐标或影像文件是否正确")

如果按上述方法修改后仍有问题，逐一检查以下项：

• 确认影像没有经过错误的裁剪/重投影操作：如果做过预处理，要保证处理后的影像同步更新了transform和crs属性，不要沿用原始影像的仿射变换参数
• 确认经纬度顺序没有搞反：rowcol第一个参数传入东向坐标（经度/投影东坐标），第二个传入北向坐标（纬度/投影北坐标），不要调换顺序
• 确认使用的影像确实覆盖目标坐标：可先查询对应景Landsat8的幅面范围，排除点在影像外的情况

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