从Google搜索中并不能立即看出来,因此认为在此处发布答案可能会很有用。还有一个关于非连续多边形的问题。
即时简单的答案:使用命令:
centroids <- getSpPPolygonsLabptSlots(polys)
R,sp中用于总体空间包的SpatialPolygonsDataFrame R数据类的类说明)
br />在以下代码中,该代码应在任何R安装中都可以复制(尝试!)。重心(应该在运行代码后显示该图):
到目前为止一切顺利。但是,当您在QGIS中计算多边形质心(菜单:“矢量” |“几何” |“多边形质心”)时,非连续多边形的结果会略有不同:
是3的东西:
快速简便的答案
警告使用R来计算非连续多边形的质心的人
关于应如何解决的问题在R中完成以正确说明多部分(不连续)的多边形
#1 楼
首先,我找不到任何文档说coordinates
或getSpPPolygonsLabptSlots
返回质心。实际上,后者的功能现在显示为“已弃用”,并应发出警告。包。进行gCentroid
会发现这一点。评论
根据许多R的空间包的开发者Roger Bivand的说法,它的作用是:“是。?“ Polygons-class”的类文档没有说明是这种情况,因为其他点可能被有效地插入了标签点。默认构造函数使用Polygons对象中最大的非孔环的质心。” -说明不连续。 stat.ethz.ch/pipermail/r-help/2009-February/187436.html。确认:gCentroid(sids,byid = TRUE)确实可以解决问题。
–RobinLovelace
2012年12月9日在11:16
不适用于我...即使应用gCentroid(polygon,byid = TRUE),我的质心也位于两个多边形之间。.因此,我假设那些被视为多部分多边形?如何将它们分开? points(coordinates(SC.tracks),pch = 16,col =“ blue”,cex = 0.4),但是,Produce不会从多边形产生质心...谢谢!
– Maycca
16年5月12日在18:26
stat.ethz.ch的链接不再起作用。仅出于完整性考虑,我几乎可以肯定现在可以在这里找到答案:r.789695.n4.nabble.com/…
– Exocom
16年6月2日在6:51
#2 楼
这是一种使用科幻小说的方法。如我所展示的,来自sf :: st_centroid和rgeos :: gCentroid的结果相同。library(sf)
library(ggplot2)
# I transform to utm because st_centroid is not recommended for use on long/lat
nc <- st_read(system.file('shape/nc.shp', package = "sf")) %>%
st_transform(32617)
# using rgeos
sp_cent <- gCentroid(as(nc, "Spatial"), byid = TRUE)
# using sf
sf_cent <- st_centroid(nc)
# plot both together to confirm that they are equivalent
ggplot() +
geom_sf(data = nc, fill = 'white') +
geom_sf(data = sp_cent %>% st_as_sf, color = 'blue') +
geom_sf(data = sf_cent, color = 'red')
评论
谢谢。为了完整起见,stackoverflow.com/questions/49343958/…显示它们非常接近,但并不完全相同。也就是说,我不知道哪个可以被认为是“更好”。
– hmhensen
19/12/18在23:02
#3 楼
我为克服此问题所做的工作是生成一个对多边形进行负缓冲的函数,直到该多边形足够小以至于需要凸多边形。使用的功能是centroid(polygon)
#' find the center of mass / furthest away from any boundary
#'
#' Takes as input a spatial polygon
#' @param pol One or more polygons as input
#' @param ultimate optional Boolean, TRUE = find polygon furthest away from centroid. False = ordinary centroid
require(rgeos)
require(sp)
centroid <- function(pol,ultimate=TRUE,iterations=5,initial_width_step=10){
if (ultimate){
new_pol <- pol
# For every polygon do this:
for (i in 1:length(pol)){
width <- -initial_width_step
area <- gArea(pol[i,])
centr <- pol[i,]
wasNull <- FALSE
for (j in 1:iterations){
if (!wasNull){ # stop when buffer polygon was alread too small
centr_new <- gBuffer(centr,width=width)
# if the buffer has a negative size:
substract_width <- width/20
while (is.null(centr_new)){ #gradually decrease the buffer size until it has positive area
width <- width-substract_width
centr_new <- gBuffer(centr,width=width)
wasNull <- TRUE
}
# if (!(is.null(centr_new))){
# plot(centr_new,add=T)
# }
new_area <- gArea(centr_new)
#linear regression:
slope <- (new_area-area)/width
#aiming at quarter of the area for the new polygon
width <- (area/4-area)/slope
#preparing for next step:
area <- new_area
centr<- centr_new
}
}
#take the biggest polygon in case of multiple polygons:
d <- disaggregate(centr)
if (length(d)>1){
biggest_area <- gArea(d[1,])
which_pol <- 1
for (k in 2:length(d)){
if (gArea(d[k,]) > biggest_area){
biggest_area <- gArea(d[k,])
which_pol <- k
}
}
centr <- d[which_pol,]
}
#add to class polygons:
new_pol@polygons[[i]] <- remove.holes(new_pol@polygons[[i]])
new_pol@polygons[[i]]@Polygons[[1]]@coords <- centr@polygons[[1]]@Polygons[[1]]@coords
}
centroids <- gCentroid(new_pol,byid=TRUE)
}else{
centroids <- gCentroid(pol,byid=TRUE)
}
return(centroids)
}
#Given an object of class Polygons, returns
#a similar object with no holes
remove.holes <- function(Poly){
# remove holes
is.hole <- lapply(Poly@Polygons,function(P)P@hole)
is.hole <- unlist(is.hole)
polys <- Poly@Polygons[!is.hole]
Poly <- Polygons(polys,ID=Poly@ID)
# remove 'islands'
max_area <- largest_area(Poly)
is.sub <- lapply(Poly@Polygons,function(P)P@area<max_area)
is.sub <- unlist(is.sub)
polys <- Poly@Polygons[!is.sub]
Poly <- Polygons(polys,ID=Poly@ID)
Poly
}
largest_area <- function(Poly){
total_polygons <- length(Poly@Polygons)
max_area <- 0
for (i in 1:total_polygons){
max_area <- max(max_area,Poly@Polygons[[i]]@area)
}
max_area
}
评论
速度慢,但结果却很好。它居中,可以很好地放置标签
–巴斯蒂安
18年1月11日在13:35
这正是芝加哥复杂病区边界所需要的。包含在包装中,这将是惊人的! gist.github.com/geneorama/40a5fd67fed2b4a5db469ce998c693ed
–geneorama
19/12/31在22:10
评论
我需要知道如何引用上面解释的函数质心。谢谢欢迎使用GIS StackExchange!作为新用户,请参观。这似乎是一个新问题,而不是对该问题的答案。请发布新问题。