计算同一区域粗分辨率与对应细分辨率栅格像元的欧氏距离
计算粗/细分辨率栅格像元质心的欧氏距离(Terra包实现)
需求:现有同一区域的两个栅格图层(粗分辨率460m,细分辨率100m),需计算每个粗分辨率栅格像元质心与其内部所有细分辨率栅格像元质心的欧氏距离(非最近距离),已完成栅格转点要素的操作,以下是后续步骤:
步骤1:建立细点与粗栅格的归属关联
首先明确每个细分辨率点属于哪个粗分辨率栅格像元,通过将粗栅格转为面要素后进行空间连接实现:
# 加载Terra包(若未加载) library(terra) # 将粗分辨率栅格转为单个像元对应的面要素(不合并) cr_poly = as.polygons(cr, dissolve = FALSE, na.rm = TRUE) # 为每个粗栅格面添加唯一ID,用于后续关联 cr_poly$coarse_id = seq_len(nrow(cr_poly)) # 空间连接:给每个细点匹配所属的粗栅格ID fr_p = terra::extract(cr_poly, fr_p, bind = TRUE) # 过滤掉不在粗栅格范围内的细点(若存在) fr_p = fr_p[!is.na(fr_p$coarse_id), ]
步骤2:提取坐标并计算欧氏距离
提取粗、细点的坐标信息,按归属ID关联后计算每对质心的欧氏距离:
# 提取粗点的坐标与对应ID cr_coords = data.frame( coarse_id = cr_poly$coarse_id, cr_x = xFromPoints(cr_p), cr_y = yFromPoints(cr_p) ) # 提取细点的坐标、所属粗栅格ID及原始值 fr_coords = data.frame( coarse_id = fr_p$coarse_id, fr_x = xFromPoints(fr_p), fr_y = yFromPoints(fr_p), fr_value = fr_p$B10_median ) # 关联粗、细点的坐标数据 merged_coords = merge(fr_coords, cr_coords, by = "coarse_id") # 计算平面欧氏距离(当前坐标系EPSG:7767为平面坐标系,无需投影转换) merged_coords$euclidean_distance = sqrt( (merged_coords$fr_x - merged_coords$cr_x)^2 + (merged_coords$fr_y - merged_coords$cr_y)^2 )
步骤3:生成结果输出(匹配示例图需求)
若需要输出细分辨率栅格格式的结果(每个细像元值为到所属粗像元质心的距离),执行以下代码:
# 创建与细栅格结构一致的结果栅格 result_rast = fr # 将计算得到的距离值赋值给结果栅格(按细点顺序匹配) values(result_rast) = merged_coords$euclidean_distance # 保存结果栅格到本地 writeRaster(result_rast, "path/fine_distance_to_coarse_centroid.tif", overwrite = TRUE)
可选:按粗栅格汇总距离统计值
若需要对每个粗栅格内的细像元距离进行统计(如均值、最大值),可使用分组统计:
# 加载dplyr包用于分组统计 library(dplyr) distance_stats = merged_coords %>% group_by(coarse_id) %>% summarise( avg_distance = mean(euclidean_distance), max_distance = max(euclidean_distance), min_distance = min(euclidean_distance) )
内容的提问来源于stack exchange,提问作者Nikos
相关产品推荐
相关产品推荐

