QGIS 栅格 Rescale 结果异常解决方案

在做多因子叠加分析、多指标综合评价的时候，小编最常干的一件事，就是把一堆来自不同来源的栅格指标统一到 0 到 1 的标准化范围，方便后续加权叠加。理论上来说，只要按照最小值和最大值做一个线性缩放，就应该能拿到 0 到 1 的结果。

但在实际使用 QGIS 的栅格 Rescale 工具、栅格计算器，或者 GDAL 的 gdal_translate 命令时，不少同学会发现：结果栅格的最小值和最大值总是“差一点点”，比如期望是 0 到 1，实际却是 0.006 多到 0.88 多，或者 0 到 0.99999999999999，看着就让人不放心。今天小编就结合一条来自 GIS Stack Exchange 的问答，聊聊问题背后的原因和一个更靠谱的解决思路。

原问答链接：https://gis.stackexchange.com/questions/500062/qgis-raster-rescale-gives-bad-results，GDAL 参数说明可参考官方文档 https://gdal.org/en/stable/programs/gdal_translate.html#cmdoption-gdal_translate-scale。

问题原因

表面上是“重标定之后的栅格最小值和最大值对不上”，本质上主要有两点原因。

第一，线性缩放依赖于正确的源栅格最小值和最大值。如果 Rescale 或者 gdal_translate 使用的最小值、最大值不是栅格真实统计出来的值，而是手工输入或者来自图层属性中的近似统计值，那么换算之后的结果自然会偏离预期。这也是为什么回答者第一步就建议使用 gdalinfo 重新检查原始栅格的统计信息。

第二，栅格运算是在浮点数空间里完成的，浮点数天生就存在精度误差。也就是说，即便你从理论上给出了完全精确的最小值和最大值，运算后的结果在机器里也很难刚好等于零或者一，往往会得到一个极小的负数或者 0.99999999999999 这样的值。从数学角度看，这些值已经足够接近 0 和 1，但在 GIS 软件的属性窗口里却显得有些“别扭”，于是就容易让人误以为缩放出了问题。

解决方法

小编按照问答中的思路，总结了一套更稳妥的操作流程，既适用于命令行环境，也可以迁移到 QGIS 的图形界面中。

第一步，使用 gdalinfo 精确获取源栅格的真实最小值和最大值。

示例命令：

gdalinfo REM__IDW.tif -stats

在输出结果的元数据部分，可以看到类似这样的信息：

STATISTICS_MAXIMUM=53.368576049805
STATISTICS_MINIMUM=-1.4455108642578

这里的 STATISTICS_MINIMUM 和 STATISTICS_MAXIMUM 就是后续缩放要用的真实数值。小编的经验是，尽量直接复制这两个值，避免自己抄写时出现位数或者小数点错误。

第二步，使用 gdal_translate 做线性缩放，将栅格拉伸到 0 到 1 范围。

示例命令：

gdal_translate -scale -1.4455108642578 53.368576049805 0 1 \
  REM__IDW.tif scaled.tif

这个命令的含义，可以理解为：把原始栅格中最小值约为 -1.4455 的像元映射为 0，把最大值约为 53.3686 的像元映射为 1，其余像元按线性比例映射到 0 到 1 之间。注意 -scale 后面直接跟四个数字，分别是源最小值、源最大值、目标最小值和目标最大值。

第三步，检查结果栅格的统计信息，确认缩放效果。

示例命令：

gdalinfo scaled.tif -stats

测试结果得到的最小值约为 -2.2551405187698 乘以 10 的负 16 次方，最大值约为 0.99999999999999。这个最小值虽然不是严格意义的零，但从数量级来看已经非常接近零，可以认为是浮点数计算引入的极小误差。对于绝大多数 GIS 分析和制图场景，这样的误差完全可以忽略。

第四步，如果希望结果范围更加贴近 0 到 1，可以使用 -exponent 参数做微调。

在 gdal_translate 的文档中，-exponent 可以对缩放后的值再做一次指数变换，用得恰当时可以“挤压”一下数值范围，使最小值更接近零。示例命令：

gdal_translate -scale -1.4455108642578 53.368576049805 0 1 \
  REM__IDW.tif scaled2.tif -exponent 1

再用 gdalinfo 查看统计信息，最小值就是 0，最大值仍然是 0.99999999999999。这里的 0.99999999999999 对于计算机和实际应用来说已经可以视为 1。

如果你更习惯在 QGIS 图形界面中操作，也可以在处理工具箱中找到 GDAL 的 Translate 工具，把上面的参数对应填入输入文件、输出文件和 scale 相关参数中，效果是一致的。

第五步，如果只是在 QGIS 内部进行简单验证，也可以使用栅格计算器手动写出线性缩放表达式，例如：

("REM__IDW@1" - (-1.4455108642578)) / (53.368576049805 - (-1.4455108642578))

这里 REM__IDW@1 表示输入栅格的第一个波段。表达式的含义是：先把每个像元减去最小值，再除以极差，最后得到 0 到 1 范围内的结果。计算完成后，可以在图层属性里查看统计信息，同样会看到非常接近 0 和 1 的浮点数。

其他方法

除了直接使用 gdal_translate 的 -scale 选项之外，小编平时还会用到两种思路：

一种是完全在 QGIS 的栅格计算器里完成标准化，把上面的线性表达式封装成模型或脚本，方便批量复用。对于喜欢图形界面的同学，这种方式上手更快，也便于和后续的加权叠加、重分类流程串联在一起。

另一种是结合数据类型控制输出结果。比如在做完 0 到 1 的缩放后，可以选择使用 Float32 作为输出数据类型，既能保留足够的精度，又能避免 Float64 带来的文件体积膨胀问题。如果后续只需要分级制图，还可以继续对栅格做整数化处理。

不管使用哪种方法，有两点共通的经验值得记住：一是尽量基于真实统计值做缩放，而不是肉眼估计；二是正确理解浮点数精度误差，不要被 0.99999999999999 这样的数字吓到。

总结

QGIS 栅格 Rescale 结果“对不上”的问题，看似是工具不准，实则是源最小值和最大值的选择加上浮点数精度的双重影响。通过先用 gdalinfo 获取真实统计值，再用 gdal_translate 的 -scale 加上可选的 -exponent 参数进行缩放，就能得到一个数值范围非常稳定的 0 到 1 标准化栅格。

对于做多因子叠加、适宜性评价的 GISer 来说，理解这套缩放思路不仅能避免“结果看着不顺眼”的焦虑，也能在后续的栅格建模中做到心中有数。如果你在 QGIS 中有更好用的栅格标准化技巧，或者有其他处理浮点误差的经验，欢迎留言讨论。