提取坐标至不同土地利用类型距离时遇terra::distance范围未知问题

问题 #

需求：提取家庭坐标到各类土地利用类型（水域、耕地、建筑等）的距离，最终为每个坐标生成对应距离列。
数据来源：SENTINEL-2 10M土地利用/覆盖数据（需合并两个栅格得到完整研究区）、带10km缓冲区的研究区栅格、家庭坐标数据。
操作流程：读取并投影数据→合并土地利用栅格→裁剪至研究区→提取水域类栅格后调用terra::distance计算距离，此时出现警告：[distance] unknown extent，寻求解决方法。

尝试代码：

# 读取家庭坐标数据框
df <- import('clean_R/clean_gps.csv') 

# 读取研究区范围栅格并投影
extent <- rast("data/areas/extent.tif")
extent <- project(extent, "EPSG:32628")

# 读取家庭点数据
points <- st_read('data/points/points.shp') 

## 预处理2020年土地利用数据
lu_20a <- rast("2020/28P_20200101-20210101.tif") # 导入土地利用栅格
lu_20b <- rast("2020/28Q_20200101-20210101.tif") # 导入土地利用栅格
lu_20a <- project(lu_20a, "EPSG:32628") # 重投影栅格
lu_20b <- project(lu_20b, "EPSG:32628") # 重投影栅格
lu_20 <- merge(lu_20a, lu_20b) # 合并两个栅格
plot(lu_20) # 检查栅格范围
plot(st_geometry(study_area), col = "dark sea green", add =TRUE) # 叠加索马里边界
crs(lu_20)  # 检查坐标系

remask = resample(lu_20, extent) # 确保栅格坐标一致
lu_20_2 <- crop(lu_20, remask)  # 裁剪至研究区范围
plot(lu_20_2)

# 提取2020年水域数据
water20 <- lu_20  # 创建新栅格用于提取水域（水域分类值为1）
coltab(water20)
water20[water20>1] <- NA # 保留仅水域值
water20[water20<1] <- NA # 保留仅水域值
water20[water20 == 0] <- NA # 保留仅水域值

plot(water20) # 检查是否仅保留水域
water20_dis <- distance(water20) # 计算每个像素到最近水域的距离
df$lu20_water_points <- extract(water20_dis, points) # 提取每个坐标点的距离值

解决建议 #

1. 检查水域栅格的有效性 #

[distance] unknown extent警告通常指向两个核心问题：要么输入栅格丢失了有效范围信息，要么栅格中没有任何非NA的水域像素。

• 先运行以下代码验证：

  print(water20) # 查看栅格范围、分辨率等元信息
  sum(!is.na(values(water20))) # 统计非NA像素数量
  

  如果非NA像素数为0，说明水域提取逻辑错误：
  • 确认SENTINEL-2土地利用数据中水域的分类值是否真的是1，可通过unique(values(lu_20))查看所有分类值；
  • 简化提取逻辑，直接用water20 <- lu_20 == 1，避免多次赋值导致的逻辑冲突。

2. 修正栅格使用逻辑 #

代码中存在明显逻辑矛盾：裁剪得到lu_20_2后，提取水域却用了原始的lu_20，这会导致栅格范围与研究区不匹配，甚至丢失有效范围信息。修正如下：

# 改用裁剪后的栅格提取水域
water20 <- lu_20_2
water20[water20 != 1] <- NA # 直接保留分类值为1的像素

3. 确保栅格范围与投影完整 #

• 运行crs(water20)确认坐标系是否为EPSG:32628，投影过程出错会导致范围信息丢失；
• 如果extent.tif是研究区掩码，建议用mask()替代单纯的crop()，确保栅格仅保留研究区内的有效像素：

  lu_20_2 <- mask(crop(lu_20, extent), extent)
  

4. 匹配点数据与栅格的坐标系 #

确保家庭点数据的坐标系与水域栅格完全一致，否则提取距离时会出错：

points <- st_transform(points, crs(water20))

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