[R-sig-Geo] Obtaining 3d surface area using polygon layer and DEM grid

Barry Rowlingson b.rowlingson at lancaster.ac.uk
Wed Apr 14 17:12:55 CEST 2010


On Tue, Apr 13, 2010 at 5:46 PM, Barry Rowlingson
<b.rowlingson at lancaster.ac.uk> wrote:

>  I don't know of any R implementation of a surface area algorithm for
> R at the moment.

 Well, I do now because I've just written it. It's an implementation
of the Jenness algorithm. He doesn't say how he copes with the edges,
so what I do is extend the grid by 1 pixel all round and fill the edge
with the nearest value from the grid.

 It's written in C and does a 2000x2000 grid in a couple of seconds.

 > b=matrix(runif(2000*2000),2000,2000)
 > system.time(print(sacalc(b)))
 [1] 4592546
    user  system elapsed
   1.868   0.076   1.941

 now currently it doesn't take into account the grid spacing - it
assumes the cells are 1 unit apart and that's the same scale as the
values. Simple enough to add scaling but fiddly to integrate it into
SpatialGrid objects. Matrices are so much simpler :)

 Also, it's one C file and one R file, no documentation, and I've not
tested it hardly, and I can't build it for Windows.

Barry

-- 
blog: http://geospaced.blogspot.com/
web: http://www.maths.lancs.ac.uk/~rowlings
web: http://www.rowlingson.com/
twitter: http://twitter.com/geospacedman
pics: http://www.flickr.com/photos/spacedman



More information about the R-sig-Geo mailing list