为何terra::mask与terra::crop在touches=TRUE时结果不同？哪个更可靠？

问题背景 #

在测试terra::mask的示例时发现，touches参数在terra::mask和terra::crop(..., mask=TRUE)中的作用存在细微但明确的差异：当touches=FALSE时，两者的栅格统计结果完全一致；但当touches=TRUE时，结果出现明显差异。可复现代码如下：

library(terra)
#> terra 1.8.50

f <- system.file("ex/elev.tif", package = "terra")
r <- terra::rast(f)

f <- system.file("ex/lux.shp", package = "terra")
v <- terra::vect(f)[1, ]

# With touches=FALSE
bench::mark(
  mask = terra::mask(r, v, touches = FALSE) |>
    terra::global("sum", na.rm = TRUE),
  crop = terra::crop(r, v, mask = TRUE, touches = FALSE) |>
    terra::global("sum", na.rm = TRUE),
  check = TRUE
)
#> # A tibble: 2 × 6
#>   expression      min   median `itr/sec` mem_alloc `gc/sec`
#>   <bch:expr> <bch:tm> <bch:tm>     <dbl> <bch:byt>    <dbl>
#> 1 mask         2.48ms   2.53ms      393.     143KB     18.4
#> 2 crop         2.43ms   2.49ms      399.     773KB     22.6

# With touches=TRUE
bench::mark(
  mask = terra::mask(r, v, touches = TRUE) |>
    terra::global("sum", na.rm = TRUE),
  crop = terra::crop(r, v, mask = TRUE, touches = TRUE) |>
    terra::global("sum", na.rm = TRUE),
  check = TRUE
)
#> Error: Each result must equal the first result:
#> `mask` does not equal `crop`

1. 差异产生的原因 #

两者的差异源于处理栅格与矢量边界相交时的逻辑顺序和底层实现不同：

• terra::mask保留原始栅格的完整范围，直接对每个原始栅格单元判断是否与矢量要素相交（当touches=TRUE时，只要单元边界与矢量边界接触就会被保留）。
• terra::crop(..., mask=TRUE)先将栅格裁剪到矢量要素的外接矩形范围，再在裁剪后的栅格上执行mask操作。裁剪步骤会改变边缘栅格单元的空间范围，导致这些单元与矢量边界的相交判断结果和原始栅格上的判断出现偏差。

简言之，mask是全局范围内的原始单元判断，而crop+mask是先缩小范围再做判断，边缘单元的判定逻辑在两种场景下不一致。

2. 算法可靠性分析 #

两者的可靠性取决于你的业务需求：

• 若需要严格基于原始栅格的单元边界进行相交判断，terra::mask更可靠——它直接在原始栅格的所有单元上执行判断，没有经过裁剪带来的范围变更，结果完全贴合原始栅格的空间属性。
• 若你的目标是先裁剪到矢量的外接矩形再做mask，terra::crop(..., mask=TRUE)的结果符合预期，但需注意其边缘单元的判断是基于裁剪后的栅格范围，而非原始栅格。

从空间逻辑的严谨性角度，terra::mask的判断逻辑更直接对应touches参数的定义，而crop+mask的结果受裁剪步骤影响，可能偏离原始栅格的单元属性。

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