在SpatialPolygon中计算Centroid WITHIN/INSIDE

lok*_*oki 7 gis r polygon geospatial centroid

在像ArcMap这样的软件中,可以为多边形的多边形创建质心.在如下所示的情况下,这是必要的.

R其中可以计算空间多边形的质心rgeos::gCentroid().但是,无法强制计算多边形内的质心.

library(rgdal)
library(rgeos)

x <- readWKT("POLYGON ((1441727.5096940901130438 6550163.0046194596216083, 
             1150685.2609429201111197 6669225.7427449300885201, 
             975398.4520359700545669 6603079.7771196700632572, 
             866257.6087542800232768 6401334.5819626096636057, 
             836491.9242229099618271 6106985.0349301798269153, 
             972091.1537546999752522 5835786.5758665995672345, 
             1547561.0546945100650191 5782869.8033663900569081, 
             1408654.5268814601004124 5600968.3978968998417258, 
             720736.4843787000281736 5663807.0652409195899963, 
             598366.4479719599476084 6001151.4899297598749399, 
             654590.5187534400029108 6341803.2128998702391982, 
             869564.9070355399744585 6784981.1825891500338912, 
             1451649.4045378800947219 6788288.4808704098686576, 
             1441727.5096940901130438 6550163.0046194596216083))")
plot(x)
Run Code Online (Sandbox Code Playgroud)

这是多边形 x

gCentroid()创建一个质心,在这个特定情况下,它位于多边形之外.尽管几何上是正确的,但某些应用程序需要多边形内的质心,因为它们可以由ArcMap计算.

xCent <- gCentroid(x, byid = TRUE)
points(xCent, col = "red", pch = 16)
Run Code Online (Sandbox Code Playgroud)

在此输入图像描述

所需的输出(来自ArcMap)如下所示:

在此输入图像描述

是否有可能在R中生成像这样的质心?

编辑:

经过一番挖掘,结果发现ArcMap在Polygon中选择了一个随机点:

"对于输入多边形:输出点将位于多边形内."

因此问题必须是:是否有一个函数可以在多边形的任意位置创建一个点?

lok*_*oki 7

sf解决方案

随着软件包的出现sf,事情变得容易了一些。只需使用:

library(sf)
y <- st_as_sf(x) # only necessary when you don't already have an sf object
st_point_on_surface(y)
Run Code Online (Sandbox Code Playgroud)

“返回一个保证位于(多)曲面上的点”。

sp解决方案

正如问题更新中所指出的,ArcMap 似乎只是将一个点放置在多边形内的随机位置。gPointsOnSurface(..., n = 1, type = 'random')这也可以通过以下方式实现。

xCent2 <- gPointOnSurface(x, byid = T)
points(xCent2, col = "blue", pch = 16)
Run Code Online (Sandbox Code Playgroud)

我编写了这个函数,它首先找到质心,如果它不在内部(即它不与多边形重叠/相交),则将其替换为表面上的点。此外,它返回一个新列,指示一个点是否是真正的质心。

gCentroidWithin <- function(pol) {
  require(rgeos)

  pol$.tmpID <- 1:length(pol)
  # initially create centroid points with gCentroid
  initialCents <- gCentroid(pol, byid = T)

  # add data of the polygons to the centroids
  centsDF <- SpatialPointsDataFrame(initialCents, pol@data)
  centsDF$isCentroid <- TRUE

  # check whether the centroids are actually INSIDE their polygon
  centsInOwnPoly <- sapply(1:length(pol), function(x) {
    gIntersects(pol[x,], centsDF[x, ])
  })

  if(all(centsInOwnPoly) == TRUE){
        return(centsDF)
    }

  else {
    # substitue outside centroids with points INSIDE the polygon
    newPoints <- SpatialPointsDataFrame(gPointOnSurface(pol[!centsInOwnPoly, ], 
                                                        byid = T), 
                                        pol@data[!centsInOwnPoly,])
    newPoints$isCentroid <- FALSE
    centsDF <- rbind(centsDF[centsInOwnPoly,], newPoints)

    # order the points like their polygon counterpart based on `.tmpID`
    centsDF <- centsDF[order(centsDF$.tmpID),]

    # remove `.tmpID` column
    centsDF@data <- centsDF@data[, - which(names(centsDF@data) == ".tmpID")]

    cat(paste(length(pol), "polygons;", sum(centsInOwnPoly), "actual centroids;", 
              sum(!centsInOwnPoly), "Points corrected \n"))

    return(centsDF)
  }
Run Code Online (Sandbox Code Playgroud)