Terra中负缓冲区处理复杂多面状要素时出现异常结果

环境信息 #

• terra版本：1.7-55（近期更新）
• 操作系统：Windows 11
• 硬件：13代Intel(R) Core(TM) i7-1355U 1.70 GHz、16GB内存

在使用terra包的buffer函数生成负缓冲区时，遇到三类异常问题，以下是问题描述、临时解决办法及复现代码：

问题1：聚合面生成超出原要素范围的异常负缓冲区 #

现象：对聚合后的面状要素执行负缓冲区操作，结果出现位于原要素外的异常面。
临时修复：

1. 对原始要素执行disagg拆分后再计算缓冲区；
2. 调整聚合前要素的几何范围，避免包含跨度极大的面。

复现代码：

library(terra)

# 问题1：聚合面生成异常负缓冲区
# 生成三个面状要素的SpatVector
p1<-vect("POLYGON ((-66.940374 45.008097, -66.971762 45, -66.986546 45.008447, -66.940267 45.026992, -66.940374 45.008097))")
p2<-vect("POLYGON ((-66.828 45.032592, -66.829 45.035, -66.824 45.031, -66.8301 45.0275, -66.834 45.027, -66.828 45.032592))")
p3<-vect("POLYGON ((-67 46, -66.2 46, -66.2 45.65, -66.37 45.45, -66.37 45.105, -66.459 45.05, -66.46 45.17,
-67 45.15, -67 46))")

wkt.p<-rbind(p1, p2, p3)
crs(wkt.p)<-"+proj=longlat +datum=WGS84"

hol_buff<-buffer(aggregate(wkt.p), -2000)
plot(wkt.p, ylim=c(45, 46))
plot(hol_buff, add=T, col="green")

# 拆分后的面无异常负缓冲区
hol_buff<-buffer(wkt.p, -2000)
plot(hol_buff, add=T, col="red")

# 修改p3上限后聚合面无异常
p1<-vect("POLYGON ((-66.940374 45.008097, -66.971762 45, -66.986546 45.008447, -66.940267 45.026992, -66.940374 45.008097))")
p2<-vect("POLYGON ((-66.828 45.032592, -66.829 45.035, -66.824 45.031, -66.8301 45.0275, -66.834 45.027, -66.828 45.032592))")
p3<-vect("POLYGON ((-67 45.8, -66.2 45.8, -66.2 45.65, -66.37 45.45, -66.37 45.105, -66.459 45.05, -66.46 45.17,
-67 45.15, -67 45.8))")

wkt.p<-rbind(p1, p2, p3)
crs(wkt.p)<-"+proj=longlat +datum=WGS84"

hol_buff<-buffer(aggregate(wkt.p), -2000)
plot(wkt.p, ylim=c(45, 46))
plot(hol_buff, add=T, col="blue")

问题2：带内环的面状要素负缓冲区异常 #

现象：带孔洞（内环）的面执行负缓冲区操作时，缓冲区会在孔洞外生成，甚至超出外围面边界，与预期效果不符。
临时修复：使用fillHoles移除面要素的内环后，再进行缓冲区计算。

复现代码：

# 问题2：带孔面的负缓冲区超出边界
par(mfrow=c(1,2))
# 带内环的面
wkt.y<-"POLYGON ((-66.902 45.013, -66.901 45.01, -66.906 45.0114,  -66.9045 45.0155, -66.902 45.0155, -66.902 45.013),
(-66.903371 45.012734, -66.903372 45.012914, -66.903626 45.012913, -66.903624 45.012733, -66.903371 45.012734))"
y<-vect(wkt.y,crs="+proj=longlat +datum=WGS84")
plot(y, main="terra output")
plot(buffer(y, -40), add=T, col="red")

# 无内环的面
wkt.w<-"POLYGON ((-66.902 45.013, -66.901 45.01, -66.906 45.0114,  -66.9045 45.0155, -66.902 45.0155, -66.902 45.013))"
w<-vect(wkt.w,crs="+proj=longlat +datum=WGS84")
plot(y, main="desired output")
plot(buffer(w, -40), add=T, col="red")# 理想的负缓冲区效果
plot(y, add=T)# 叠加原始面
# 可通过fillHoles移除内环
w2<-fillHoles(y)

# 内环位于边缘的面
par(mfrow=c(1,2))
wkt.y<-"POLYGON ((-66.902 45.012,  -66.906 45.0114,  -66.906 45.0155, -66.902 45.0155, -66.902 45.012),
(-66.903371 45.012734, -66.903372 45.012914, -66.903626 45.012913, -66.903624 45.012733, -66.903371 45.012734))"
y<-vect(wkt.y,crs="+proj=longlat +datum=WGS84")
plot(y, xlim=c(-66.907, -66.901), ylim=c(45.01, 45.017), main="terra output")
plot(buffer(y, -200), add=T, col="red", border="forestgreen", lwd=2)
plot(y, add=T)# 叠加原始面

# 无内环的面
wkt.w<-"POLYGON ((-66.902 45.012,  -66.906 45.0114,  -66.906 45.0155, -66.902 45.0155, -66.902 45.012))"
w<-vect(wkt.w,crs="+proj=longlat +datum=WGS84")
plot(y, main="desired output")
plot(buffer(w, -200), add=T, col="red")# 理想的负缓冲区效果
plot(y, add=T)# 叠加原始面

问题3：R无规律频繁崩溃 #

现象：执行缓冲区操作过程中，R进程无规律频繁崩溃。
临时修复：替换输入文件的原始投影CRS为非投影坐标系。

内容的提问来源于stack exchange，提问作者kl-higgins