rgeos – Update

[This article was first published on GSoC 2010 R, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

So I have been remiss in posting updates of progress on rgeos as my summer has gotten busier but the project is progressing smoothly and everything is on schedule. We have just recently passed an important personal milestone, parsing and running of all the JTS unit tests are now working and the rgeos code passes all but a handful of tests. (Failing tests appear mostly to do with issues surrounding dimension collapse which, I do not believe, is currently handled by geos). We are nearly feature complete and will now be focusing on testing and writing up documentation / examples.

Below are a couple of examples that came up as answers to questions posted to the R-sig-geo mailing list. Note that the naming convention of the geos functions has changed to use a single g (for geos) as a prefix instead of using the awkward RGEOS per Hadley’s earlier suggestion. Old function names are retained for the time being but will throw a depricated warning.

Pruning a Delaunay tessalation using a polygon:

library(rgeos)
library(spdep)

#bounding polygon
bound <- SpatialPolygons(list(Polygons(list(Polygon(cbind( 
x=c(0.34063939, 0.56754974, 0.95361248, 0.96464284, 0.60694389, 
0.58173163, 0.58330740, 0.91421832, 1.00403700, 0.96464284, 
0.50294332, 0.39263968, 0.22560845, 0.02548613, 0.25397225, 
-0.01233226, 0.34063939), 
y=c(1.02171515, 0.70486571, 0.92207697, 0.84834471, 0.62116963, 
0.53747356, 0.47569788, 0.35812482, 0.02134774, -0.07231216, 
0.15287015, 0.14489909, -0.03444965, 0.09707276, 0.56138672, 
0.58729265, 1.02171515) ) )),ID="bound"))) 

# points that form the basis of the Delaunay tessalation
dat <- data.frame( 
x=c(0.34433527, 0.08592805, 0.55564179, 0.03938242, 0.98044051, 
0.19835405, 0.94186612, 0.56208017, 0.31093811, 0.54341230, 
0.93508703, 0.38160473, 0.89435383, 0.55457428, 0.22406338), 
y=c(0.997803213, 0.094993232, 0.509774611, 0.615526709, 0.004314438, 
0.676662461, 0.026060349, 0.165807011, 0.596449068, 0.469553704, 
0.888788079, 0.163129754, 0.340335312, 0.621845245, 0.019412254) 
) 


dat.nb = tri2nb(dat) 


filternb = function(nb, coords, bound) { 
        sl = nb2lines(nb,as.list(rep(NA,length(nb))), coords) 
        
        # any line that crosses the polygon should be removed
        out=!gCrosses(sl,bound,byid=TRUE) 
        
        k=1 
        for(i in 1:length(nb)) { 
                l = length(nb[[i]]) 
                nb[[i]] = nb[[i]][ out[k:(k+l-1)] ] 
                k=k+l 
        } 
        
        return(nb) 
} 

dat.nb.new = filternb(dat.nb, dat, bound) 

par(mfrow=c(1,2))

plot(bound,col='red') 
plot(dat.nb,dat,add=T) 

plot(bound,col='red') 
plot(dat.nb.new,dat,add=T) 

IDing polygons containing points:

This can be done with existing functions in sp (see the overlay function) but this is how you could do it with rgeos. (I just realized that this is pretty much exactly the same as the cities example from my last post)

library(rgeos) 
gt <- GridTopology(c(0.05,0.05), c(.1,.1), c(10,10)) 
grd <- SpatialGrid(gt) 
spi <- as(grd, "SpatialPixels") 
spol <- as(spi, "SpatialPolygons") 

set.seed(1) 
x=runif(25) 
y=runif(25) 

pts = SpatialPoints(cbind(x,y)) 

ct = gContains( spol,pts,byid=c(TRUE,TRUE) ) 
colsub = apply(ct, 2,any) 

plot(spol) 
plot(spol[colsub,],add=T,col='green') 
plot(pts,add=T,pch=16) 


To leave a comment for the author, please follow the link and comment on their blog: GSoC 2010 R.

R-bloggers.com offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

Never miss an update!
Subscribe to R-bloggers to receive
e-mails with the latest R posts.
(You will not see this message again.)

Click here to close (This popup will not appear again)