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

[Barry Rowlingson]

On Wed, Apr 14, 2010 at 5:00 PM, Sharon Baruch-Mordo
<sharonb_m at yahoo.com> wrote:
> Thanks you Barry and Edzer for your replies!
>
> I should have mentioned that I'm a relatively new to spatial analysis in
> R, and unfortunately I'm not familiar with programing in C. But I'm ready to
> learn! In trying to understand the code below, I see that you are generating
> a 2000x2000 matrix with a random uniform number then using the function
> sacalc, but how is it calculating the 3d surface (I tried the web for
> documentation on sacalc but didn't get too far... )? Would that be the
> equivalent of getting a DEM clipped by the polygon of interest and then
> running it through sacalc?
>

 Ah, sorry, sacalc() is the function I wrote, and it's not included. The line:

   b=matrix(runif(2000*2000),2000,2000)

 just creates a 2000x2000 matrix of random numbers (representing
heights) between 0 and 1. It won't be a very smooth surface! The
resulting output being about 2000x2000 shows I'm getting in the right
ballpark. In fact if I change the values in "b" to be smaller, meaning
a relatively flatter surface compared to the grid size, the result
tends towards 2000x2000, which is comforting.

 sacalc just works on a 2d matrix of values, and treats them as
heights in the DEM. It can't clip it to a polygon, but Edzer's
comments on NA (missing data) values has given me some ideas...

Barry