[R-sig-Geo] voronoi (thiessen) polygons inside irregular area

[Adrian Baddeley]

Josh Firth asks about Dirichlet-Voronoi-Thiessen polygons in an irregular area.

The Dirichlet tile of each data point is the set that is closer to this point than to 
any of the other data points. 'Closer' is measured by Euclidean distance. 
Basically you're saying that you want to replace Euclidean distance by 
another measure of distance - say, the shortest-path distance inside the non-convex window.

Tinkering with the results of the Dirichlet tessellation is not going to do this;
you'd need a completely different algorithm.

Rolf Turner writes:
> I think it suffices to eliminate from each tile any polygons that do not
> contain any points of the pattern.  

This is not a solution because it removes part of the window. 

The code in Rolf's solution is not going to work if the window has holes 
(since some of the $bdry[[i]] could be holes rather than positive boundaries).
  
Also, an idiom like
    ttt$tiles <- xxx
is dangerous, and in this case, will lead to problems, because it is no longer true that
the union of tt$tiles is equal to tt$window.  The correct way to create a changed tessellation
is something like
     tess(tlles=xxx, window=www).


Prof Adrian Baddeley FAA
University of Western Australia
________________________________________
From: Rolf Turner [r.turner at auckland.ac.nz]
Sent: Wednesday, 22 January 2014 6:53 AM
To: buzzfuzz006
Cc: r-sig-geo at r-project.org; Adrian Baddeley
Subject: Re: [R-sig-Geo] voronoi (thiessen) polygons inside irregular area

The problem is that the tile corresponding to point 8 consists of
two discontiguous polygons one of which is empty of points of the
pattern.  The empty polygon is juxtaposed with one of the polygons
comprising the tile for point 10, which makes points 10 and 8 neighbours.

I think it suffices to eliminate from each tile any polygons that do not
contain any points of the pattern.  The following function does
that:

bar <- function(X) {
require(spatstat)
ttt <- dirichlet(X)
sss <- ttt$tiles
f <- function(x,pat){
     www <- lapply(x$bdry,function(y){owin(poly=y)})
     gp  <- sapply(www,function(w){
                any(inside.owin(pat$x,pat$y,w))})
     www[[seq(along=gp)[gp]]]
}
rrr <- lapply(sss,f,pat=X)
ttt$tiles <- rrr
ttt
}

Try:

junk <- bar(poly.w.points)
junk.sp <- as(junk,"SpatialPolygons")
nays <- poly2nb(junk.sp)
nays[10] # Gives 3  7 11 14.  OMMMMMMM!!!!

Rolf Turner

On 21/01/14 20:10, buzzfuzz006 wrote:

> Hi,
> Is there a simple way to draw voronoi/thiessen polygons inside an irregular
> polygon similar to this, but where the polygons created do not extend over
> the boundaries of the overlaying polygon. E.g
> Using the packages
>
> library(spatstat);library(spdep)
>
> If we had set points:
>
> x<-c(0.9,1.7,2.4,2.9,4.83, 0.73, 2.31, 3.69, 4.23, 2.86, 1.91, 4.32, 4.60,
> 1.82)
> y<-c(1.9,0.9,2.8,1.9,1.81, 1.66, 4.54, 5.66, 1.99, 4.03, 4.32, 5.98, 5.56,
> 3.41)
>
> within the set irregular polygon:
>
> x.p<-c(0.1, 6.0, 6.0, 1.0, 5.0, 3.0, 3.5, 0.1)
> y.p<-c(0.1, 1.0, 6.5, 5.5, 5.0, 1.0, 4.8, 5.0)
>
> now I set the polygon to needed dirichlet format
>
> poly.l<-list(x.p,y.p)
> names(poly.l)<-c("x","y")
>
> and create the spatial point pattern
>
> poly.w.points<-ppp(x=x,y=y,poly=poly.l)
>
> then carry out the voronoi tessellation
>
> spatstat.options("gpclib"=T)
> tess<-dirichlet(poly.w.points)
>
> unfortunately if we are to plot this it looks like this
>
> plot(tess)
> text(x,y,1:length(x))
>
> and if we check the neighbours we see that the voronoi polygons extend
> outside the polygon area
> tess.sp<-as(tess,"SpatialPolygons")
> neighs<-poly2nb(tess.sp)
>
> e.g.
> neighs[10]
> shows that point 10 is neighbours (shares boundary) with 8,9 and 13 – whilst
> I want a result whereby 10 would only be neighbours with 7,11,14 and 3.