使用Rasterio执行栅格乘法后无法写入输出栅格的问题排查

错误原因分析 #

你遇到的NotGeoreferencedWarning由以下问题导致：

• 写入输出栅格时设置transform=None，导致栅格缺失地理变换信息，无法关联空间坐标系
• CRS参数错误：输入栅格使用EPSG:3035，但代码中误写为epsg:3050
• Nodata设置为-0无实际意义，应使用合理空值（如输入栅格自带的nodata或np.nan）
• 手动指定width/height易出错，复用输入栅格的元数据更可靠

修正后的代码 #

import rasterio
import numpy as np

# 输入栅格文件路径
raster_woody_path = 'C:/Users/Spatial datasets/European countries/Input_datasets/Forest_Biomass_Map/JRC_forest_yield_upsampled_to_1km.tif'
raster2_path = 'C:/Users/Throughput_datasets/eu_nuts0_for_rents_prices_rast.tif'

# 打开输入栅格并保留元数据
with rasterio.open(raster_woody_path) as src1, rasterio.open(raster2_path) as src2:
    # 读取栅格数据
    rast1 = src1.read()
    rast2 = src2.read()
    
    # 执行乘法计算，同时处理nodata
    result = rast1 * rast2
    
    # 提取并统一nodata值
    nodata_src1 = src1.nodata if src1.nodata is not None else np.nan
    nodata_src2 = src2.nodata if src2.nodata is not None else np.nan
    result[(rast1 == nodata_src1) | (rast2 == nodata_src2)] = nodata_src1

    # 复制源栅格元数据作为输出基础
    out_meta = src1.meta.copy()
    # 更新输出元数据配置
    out_meta.update({
        "dtype": "float32",
        "nodata": nodata_src1,
        "compress": "LZW"  # 可选，添加压缩减少文件体积
    })

    # 输出文件路径
    output_path = "C:/Users/Throughput_datasets/eu_forestry_monetary_yields_unmasked.tif"

    # 写入结果栅格
    with rasterio.open(output_path, 'w', **out_meta) as dst:
        dst.write(result.astype(out_meta["dtype"]))

关键修改说明 #

1. 复用源栅格元数据：直接从输入栅格复制meta信息，确保width、height、transform、CRS等空间参数完全匹配输入栅格
2. 修正地理变换参数：元数据中自带正确的transform，无需手动设置为None
3. 合理处理Nodata：继承输入栅格的nodata值，同时将输入中nodata的位置在结果中标记为nodata
4. 数据类型对齐：将计算结果转换为输出指定的float32类型，避免类型不兼容问题
5. 可选压缩配置：添加compress: LZW参数，减少输出文件的存储空间

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