使用Rasterio导出的GeoTIFF在QGIS中显示全黑的问题求助

我使用Python的rasterio（1.3.6）和earthpy开发卫星影像处理工作流，通过栅格代数计算生成Sabins Ratio合成栅格，将每个比值作为单独波段导出为GeoTIFF。结果在Python环境中可正常可视化，但导出后在QGIS中仅显示黑色方块（仅左上角有代表兴趣区域的红点）。相关代码如下：

#file path
multi_bands = glob.glob('.././GIS/landsat_8/LC08_L2SP_090079_20230225_20230301_02_T1/LC08_L2SP_090079_20230225_20230301_02_T1_SR_*B[2:3:4:5:6:7].tif')

multi_bands.sort()

img_list = []

#reading bands
for img in multi_bands:
    with rio.open(img, 'r') as img_file:
        img_list.append(img_file.read(1))

#multiband array    
arr_st = np.stack(img_list)

# visualize all bands
ep.plot_bands(arr_st,
              cbar=False,
              cmap = 'gist_earth')
plt.show()

#avoid numpy error of zero division
np.seterr(divide='ignore', invalid='ignore')

#perform band ratios
ironOxideR = arr_st[0]/arr_st[2]
ClayHyR = arr_st[4]/arr_st[5]
ferrousR = arr_st[4]/arr_st[3]

# Create the RGB composition
sabinsRatio = np.stack((ironOxideR, ClayHyR, ferrousR))

#sabinsRatio plot
ep.plot_rgb(sabinsRatio,
            rgb=(0,1,2),
            figsize=(10, 10),
            stretch = True,
            title = 'RGB composition of Iron Oxide, Clay Hydration, and Ferrous')
plt.show()

#trying to export

#exporting sabinsRatio

# Update the metadata
out_meta = {"driver": "GTiff",
            "height": sabinsRatio.shape[1],
            "width": sabinsRatio.shape[2],
            "crs": "EPSG:32756",
            "count":3,
            "dtype":rio.float32
            }
           


# Write the mosaic raster to disk
with rio.open('./output/sabinRation2.tif', "w", **out_meta) as dest:
    dest.write(sabinsRatio, [1,2,3])

解决方案 #

1. 补充地理变换（GeoTransform）参数 #

导出的TIFF缺失地理变换信息，导致QGIS无法正确定位和渲染栅格。需要从原始影像中读取transform参数并写入输出元数据：

# 从第一个原始影像读取transform和基础元数据
with rio.open(multi_bands[0], 'r') as src:
    out_meta = src.meta.copy()
    # 更新元数据为输出需求
    out_meta.update({
        "count": 3,
        "dtype": rio.float32
    })

或者手动在原有out_meta中添加：

with rio.open(multi_bands[0], 'r') as src:
    transform = src.transform

out_meta = {
    "driver": "GTiff",
    "height": sabinsRatio.shape[1],
    "width": sabinsRatio.shape[2],
    "crs": "EPSG:32756",
    "count":3,
    "dtype": rio.float32,
    "transform": transform  # 新增地理变换参数
}

2. 处理浮点栅格的显示范围 #

比值计算后的浮点数据动态范围与QGIS默认显示设置不匹配，导致渲染异常。可以：

• 在QGIS中右键栅格图层 → 属性 → 样式，选择“最小/最大拉伸”，点击“加载”获取实际像素值范围后应用；
• 或在导出前将数据归一化到0-1范围：

def normalize(arr):
    min_val = np.nanmin(arr)
    max_val = np.nanmax(arr)
    return (arr - min_val) / (max_val - min_val)

# 对每个波段进行归一化
sabinsRatio_normalized = np.stack([normalize(band) for band in sabinsRatio])

3. 替换NaN值 #

比值计算中产生的NaN（除零或无效值）会被QGIS默认渲染为黑色，可将其替换为特定值：

sabinsRatio = np.nan_to_num(sabinsRatio, nan=0.0)

完整修改后的导出代码示例 #

# 读取原始影像元数据
with rio.open(multi_bands[0], 'r') as src:
    out_meta = src.meta.copy()
    out_meta.update({
        "count": 3,
        "dtype": rio.float32
    })

# 处理NaN值
sabinsRatio = np.nan_to_num(sabinsRatio, nan=0.0)

# 写入栅格文件
with rio.open('./output/sabinRation2.tif', "w", **out_meta) as dest:
    dest.write(sabinsRatio, [1,2,3])

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