calculate the area of an alpha shape
Dear Andr?, You're right! I assumed something that wasn't true - that $x in the ashape.obj were the points on the hull. That's what you get for free help, I guess.... Fortunately this seems to work. Notice how I index $x by the alpha.extremes field to extract just the points I want for bds: library(sp) bds <- achull.obj$ashape.obj$x[achull.obj$ashape.obj$alpha.extremes,] bds <- rbind(bds, bds[1,]) # close the ring ashape <- Polygon(bds) # convert to sp Polygon print(ashape at area) # The following plots the points in ashape superimposed on the hulls plot(achull.obj, col = "blue") plot(achull.obj$ashape.obj, add = TRUE, col = "red") points(coordinates(ashape), pch=16, cex=2, col='gold')
On Monday, November 04, 2013 04:32:42 PM Proosdij, Andre van wrote:
Many thanks Ashton! It does help a lot, but I'm afraid, it's not completely correct yet. One would expect that the area of the alpha shape is larger than the area of the alpha-convex hull (straight line segments instead of arcs). However, applying your suggestion, the area is smaller! If I plot the points and the alpha shape, the plot is correct. The area of the alpha shape seems to be calculated a polygon that includes the point (1.5, 2.0). See the updated example below. I don't know what happens, but apparently, the list of points that defines the polygon is incorrect. Any thoughts of how to adjust this? Best, Andr? x <- c(1,1,1.5,2.8,3,3) # original set of points y <- c(1,3,2,3,4,1) # original set of points #x <- c(1,1,2.8,3,3) # points without (1.5, 2.0) #y <- c(1,3,3,4,1) # points without (1.5, 2.0) z <- cbind(x,y) plot(z) ### Calculate the alpha convex hull library(alphahull) alpha <- 3 # Define alpha > 0 achull.obj <- ahull(z, alpha = alpha) ### Plot the sampled points and the alpha-convex hull object plot(achull.obj, add = TRUE, col = "blue") ### Get the length of the alpha-convex hull achulllength <- achull.obj$length achulllength ### Get the area of the alpha-convex hull achullarea <- areaahull(achull.obj) achullarea ### Plot the alpha shape plot(achull.obj$ashape.obj, add = TRUE, col = "red") ### Get the length of the alpha shape ashapelength <- achull.obj$ashape.obj$length ashapelength ### Get the area of the alpha shape library(sp) bds <- achull.obj$ashape.obj$x # matrix of coordinates in ashape bds <- rbind(bds, bds[1,]) # close the ring ashape <- Polygon(bds) # convert to sp Polygon ashapearea <- ashape at area ashapearea Ir. A.S.J. van Proosdij PhD Student t +31 317 48 1198 e andre.vanproosdij at wur.nl Naturalis Biodiversity Center, section NHN Biosystematics Group, Wageningen University Gen. Foulkesweg 37, 6703 BL Wageningen, the Netherlands www.bis.wur.nl -----Original Message----- From: Ashton Shortridge [mailto:ashton at msu.edu] Sent: maandag 4 november 2013 14:35 To: r-sig-geo at r-project.org Cc: Proosdij, Andre van Subject: Re: [R-sig-Geo] calculate the area of an alpha shape Dear Andr?, There may be a more elegant solution, but this works fine: library(sp) bds <- achull.obj$ashape.obj$x # matrix of coordinates in ashape bds <- rbind(bds, bds[1,]) # close the ring ashape <- Polygon(bds) # convert to sp Polygon print(ashape at area) On Monday, November 04, 2013 11:58:26 AM Proosdij, Andre van wrote:
Dear colleagues, I wish to calculate the area covered by the alpha shape based on a sample of points in a 2D plane. For the alpha-convex hull this works fine using the function ahull in the alphahull package. However, I want to get the area of the slightly different alpha shape. How can I do this? Below a short example. Many thanks in advance for your answers and ideas! Andr?
----- Ashton Shortridge Associate Professor ashton at msu.edu Dept of Geography http://www.msu.edu/~ashton 235 Geography Building ph (517) 432-3561 Michigan State University fx (517) 432-1671