2

I have a Data.frame consisting of Longitude and Latitude called Points1, with which I use Voronoi polygon. I also have another Data.frame of Longitude and Latitude called Points2.

I then want to count the number of occurences of Points2 in each Voronoi Polygon. Regarding the plot, I want to count the occurences of triangles in each Polygon. Sample

Sample code:

library("deldir")
Points1 <- data.frame(x=runif(10), y=runif(10))
Points2 <- data.frame(x=runif(100), y=runif(100))
x=Points1$x
y=Points1$y

vtess <- deldir(x, y)

plot(x, y, type="n", asp=1)
points(x, y, pch=20, col="red", cex=0.5)
points(Points2$x,Points2$y,pch=2)
plot(vtess, wlines="tess", wpoints="none", number=FALSE, add=TRUE, lty=1)
Improve this question
1
  • 1
    Note that voronoi/delaunay/dirichlet operations from deldir (or spatstat) with lat-long coordinates assume the lat-longs to be a cartesian coordinate system, not spherical coordinates. In a small area, convert to a local metric coordinate system.
    – Spacedman
    Commented Aug 6, 2018 at 21:44

1 Answer 1

Reset to default
2

You can do this using dirichlet from spatstat package and then quadratcount:

> library(spatstat)
> voro = dirichlet(ppp(Points1$x, Points1$y))
> quadratcount(ppp(Points2$x, Points2$y), tess=voro)
tile
 1  2  3  4  5  6  7  8  9 10 
16  7  7 16  6 16 12  4 13  3
Improve this answer
2
  • @Squeezie Beware though! you should create proper ppp obejcts with the correct "window" or points that fall outside the bounding box of your tesselation won't get counted. Although in theory voronoi tiles go off to infinity, spatstat has to crop them finitely.
    – Spacedman
    Commented Aug 7, 2018 at 7:39
  • Thanks for the reminder. I also had a problem with ppp() ignoring a duplicate points which made the indices of the points "corrupt". Solved by using duplicated()!
    – Squeezie
    Commented Aug 7, 2018 at 11:41

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.