Calculating the Distance Between Polygon and Point in R

Calculating the distance between polygon and point in R

You could use the rgeos package and the gDistance method. This will require you to prepare your geometries, creating spgeom objects from the data you have (I assume it is a data.frame or something similar). The rgeos documentation is very detailed (see the PDF manual of the package from the CRAN page), this is one relevant example from the gDistance documentation:

pt1 = readWKT("POINT(0.5 0.5)")
pt2 = readWKT("POINT(2 2)")
p1 = readWKT("POLYGON((0 0,1 0,1 1,0 1,0 0))")
p2 = readWKT("POLYGON((2 0,3 1,4 0,2 0))")
gDistance(pt1,pt2)
gDistance(p1,pt1)
gDistance(p1,pt2)
gDistance(p1,p2)

readWKT is included in rgeos as well.

Rgeos is based on the GEOS library, one of the de facto standards in geometric computing. If you don't feel like reinventing the wheel, this is a good way to go.

Distance of point feature to nearest polygon in R

Here I am using the gDistance function in the rgeos topology library. I am using a brute force double loop but it is surprisingly fast. It takes less than 2 seconds for 142 points and 10 polygons. I am sure that there is a more elgant way to perform the looping.

   require(rgeos)

# CREATE SOME DATA USING meuse DATASET
data(meuse)
coordinates(meuse) <- ~x+y
pts <- meuse[sample(1:dim(meuse)[1],142),]
data(meuse.grid)
coordinates(meuse.grid) = c("x", "y")
gridded(meuse.grid) <- TRUE
meuse.grid[["idist"]] = 1 - meuse.grid[["dist"]]
polys <- as(meuse.grid, "SpatialPolygonsDataFrame")
polys <- polys[sample(1:dim(polys)[1],10),]
plot(polys)
plot(pts,pch=19,cex=1.25,add=TRUE)

# LOOP USING gDistance, DISTANCES STORED IN LIST OBJECT
Fdist <- list()
for(i in 1:dim(pts)[1]) {
pDist <- vector()
for(j in 1:dim(polys)[1]) {
pDist <- append(pDist, gDistance(pts[i,],polys[j,]))
}
Fdist[[i]] <- pDist
}

# RETURN POLYGON (NUMBER) WITH THE SMALLEST DISTANCE FOR EACH POINT
( min.dist <- unlist(lapply(Fdist, FUN=function(x) which(x == min(x))[1])) )

# RETURN DISTANCE TO NEAREST POLYGON
( PolyDist <- unlist(lapply(Fdist, FUN=function(x) min(x)[1])) )

# CREATE POLYGON-ID AND MINIMUM DISTANCE COLUMNS IN POINT FEATURE CLASS
pts@data <- data.frame(pts@data, PolyID=min.dist, PDist=PolyDist)

# PLOT RESULTS
require(classInt)
( cuts <- classIntervals(pts@data$PDist, 10, style="quantile") )
plotclr <- colorRampPalette(c("cyan", "yellow", "red"))( 20 )
colcode <- findColours(cuts, plotclr)
plot(polys,col="black")
plot(pts, col=colcode, pch=19, add=TRUE)

The min.dist vector represents the row number of the polygon. For instance you could subset the nearest polygons by using this vector as such.

near.polys <- polys[unique(min.dist),]

The PolyDist vector contain the actual Cartesian minimum distances in the projection units of the features.

Calculating the maximum/(or minimum) distance between each polygons within a shapefile

You are likely looking for the sf::st_distance() function.

It supports polygon to polygon distances, so it will give you the shortest distance between two county polygons in the NC shapefile as a matrix. In case of neighboring polygons this will be zero.

In order to make the matrix easier to work with you may want to set the rownames and colnames; these can be then used for subsetting.

Consider this code:

library(sf) 
counties <- st_read(system.file("shape/nc.shp", package = "sf"), quiet = T)

# matrix of county to county distances
distances <- st_distance(counties, counties)

# create colnames and rownames to make work easier
colnames(distances) <- counties$NAME
rownames(distances) <- counties$NAME

# distance from county Mecklenburg / 100 values, inc. zeroes
distances["Mecklenburg",]

# counties bordering Mecklenburg cnty / 6 counties - just the zeroes from example above
distances[distances["Mecklenburg", ] == units::as_units(0, "meter"), "Mecklenburg"]

Finding the nearest distances between a vector of points and a polygon in R

I think you get only one distance value because you overwrite g in every iteration of your for-loop. I do however not not know if this is the only problem because I cannot reproduce your issue without suitable data.
Try changing the last loop to this:

g = rep(NA, dim(spdf)[1])
for(i in 1:dim(spdf)[1]){
g[i] = gDistance(spdf[i,],land_poly)
}

How to get the max distance from a point within a polygon in R


Solution using spatstat and the built-in dataset chorley:

library(spatstat)
W <- Window(chorley) # Polygonal window of the choley dataset
p <- list(x = c(350, 355), y = c(415, 425)) # Two points in polygon
plot(W, main = "")
points(p, col = c("red", "blue"), cex = 1.5)

Sample Image

v <- vertices(W) # Polygon vertices
d <- crossdist(v$x, v$y, p$x, p$y) # 2-column matrix of cross distances
i1 <- which.max(d[,1]) # Index of max dist for first (red) point
i2 <- which.max(d[,2]) # Index of max dist for second (blue) point

plot(W, main = "")
points(p, col = c("red", "blue"), cex = 1.5)
points(v$x[c(i1,i2)], v$y[c(i1,i2)], col = c("red", "blue"), cex = 1.5)

Sample Image

d[i1,1] # Max dist for first (red) point
#> [1] 21.35535
d[i2,2] # Max dist for second (blue) point
#> [1] 15.88226

Calculate minimum distance from polygon to spatial point

There are two issues with your current approach.
1.) The coordinates assignment isn't quite correct. It should be long/lat, but you're assigning it as lat/long.
2.) Directly setting the CRS in the current way will not actually alter the points in the necessary way. You need to assign an appropriate long/lat CRS first, and then perform a spTransform action.

#Open data on tin mines and define as spatial points
tin.01 <- read.csv("random/Tin Mine Location_01.csv")
coordinates(tin.01) <- c("Longitude","Latitude") ## Should be `9:8` if you wanted to stick with indexes, but using the names here is generally lower risk.
proj4string(tin.01) <- CRS(proj4string(Kinta)) ## Setting the initial projection to the same one the polygons are using. You should change this if your original data source uses some other known long/lat projection.

tin.km <- spTransform(tin.01,CRS(EPSG.3375)) ## Creating a transformed set of points for the distance calculation.

#Find distance between district and mines
gDistance(Kinta.km,tin.km,byid=TRUE) ## '0' distance means the mine is inside the district.

59
1 194384.372
2 223773.999
3 0.000
4 36649.914
5 102944.361
6 0.000
7 0.000
8 6246.066
9 0.000


Related Topics



Leave a reply



Submit