如何为240m分辨率栅格图像设置唯一ID并分解至30m分辨率？

原代码的核心问题 #

1. 自定义分解函数冗余且错误：R内置raster::disaggregate已能高效处理栅格分辨率提升，自定义函数里的y-1:rows是语法错误，且xyFromCell用法完全错误
2. ID生成存在精度风险：用paste(x,y)拼接浮点坐标，可能因精度丢失导致不同像素生成相同ID
3. 循环变量与逻辑错误：str_replace(file, ...)里的file未定义，disaggregfate_raster是拼写错误，ID关联的left_join键完全反向，逻辑不成立

优化后的完整代码（适配时间序列） #

先加载必要包：

library(raster)
library(dplyr)

1. 定义生成唯一ID栅格的函数

用单元格索引做ID，完全避免浮点精度问题：

create_id_raster <- function(r) {
  # 为每个像素生成唯一整数ID（基于栅格的单元格顺序）
  id_raster <- setValues(r, 1:ncell(r))
  names(id_raster) <- "pixel_id"
  return(id_raster)
}

2. 单栅格处理示例

# 读取240m分辨率的NDVI栅格
r <- raster("your_ndvi_file.tif")
target_res <- 30
# 计算分解倍数：240/30=8
fact <- res(r)[1] / target_res

# 生成ID栅格
id_raster <- create_id_raster(r)

# 分解NDVI和ID栅格：NDVI用双线性插值，ID用最近邻保证原始ID不变
disagg_ndvi <- disaggregate(r, fact = fact, method = "bilinear")
disagg_id <- disaggregate(id_raster, fact = fact, method = "ngb")

# 合并分解后的结果并保存
disagg_stack <- stack(disagg_ndvi, disagg_id)
writeRaster(disagg_stack, "ndvi_disagg_with_id.tif", overwrite = TRUE)

# 可选：转成CSV查看
disagg_df <- as.data.frame(disagg_stack, xy = TRUE)
write.csv(disagg_df, "ndvi_disagg_with_id.csv", row.names = FALSE)

3. 时间序列批量处理（多文件）

# 获取所有NDVI文件（带完整路径）
ndvi_files <- list.files(pattern = glob2rx("*ndvi*.tif$"), full.names = TRUE)
target_res <- 30

for (file in ndvi_files) {
  # 读取当前栅格
  r <- raster(file)
  fact <- res(r)[1] / target_res
  
  # 生成并分解ID栅格
  id_raster <- create_id_raster(r)
  disagg_ndvi <- disaggregate(r, fact = fact, method = "bilinear")
  disagg_id <- disaggregate(id_raster, fact = fact, method = "ngb")
  
  # 生成输出文件名
  base_name <- basename(file)
  out_raster <- gsub(".tif", "_disagg_with_id.tif", base_name)
  out_csv <- gsub(".tif", "_disagg_with_id.csv", base_name)
  
  # 保存结果
  disagg_stack <- stack(disagg_ndvi, disagg_id)
  writeRaster(disagg_stack, out_raster, overwrite = TRUE)
  
  # 可选保存CSV
  disagg_df <- as.data.frame(disagg_stack, xy = TRUE)
  write.csv(disagg_df, out_csv, row.names = FALSE)
  
  cat("处理完成：", base_name, "\n")
}

4. 多波段时间序列栅格处理（单文件多波段）

如果你的时间序列是一个多波段栅格文件，直接处理整个栈更高效：

# 读取多波段时间序列栅格
ndvi_stack <- stack("ndvi_time_series.tif")
target_res <- 30
fact <- res(ndvi_stack)[1] / target_res

# 基于第一个波段生成ID栅格（所有波段的单元格位置一致）
id_raster <- create_id_raster(ndvi_stack[[1]])
disagg_id <- disaggregate(id_raster, fact = fact, method = "ngb")

# 分解整个时间序列栈
disagg_ndvi_stack <- disaggregate(ndvi_stack, fact = fact, method = "bilinear")

# 添加ID波段并保存
disagg_full_stack <- addLayer(disagg_ndvi_stack, disagg_id)
writeRaster(disagg_full_stack, "ndvi_time_series_disagg_with_id.tif", overwrite = TRUE)

关键说明 #

• ID的可靠性：用1:ncell(r)生成的整数ID是绝对唯一的，不会出现浮点坐标拼接的精度问题
• 分解效率：内置disaggregate函数底层是C实现，比自定义R代码快很多，且自动处理投影和地理信息
• 时间序列兼容性：不管是多文件还是多波段格式，流程都能直接适配，保证每个时间点的分解像素都关联原始240m像素的ID

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