# [R-sig-Geo] extract mean raster value from polygon

Alexander Shenkin ashenkin at ufl.edu
Fri Mar 11 11:07:54 CET 2016

```Hello All,

I've been working to be able to make elevation profiles from a DEM along
a swath, rather than just a line (thanks to Forrest Stevens for the help
so far).  To that end, I've made a function, create_perp_buffers, that
creates polygons perpendicular to, and along, a transect (see graphic
below).

Now, I would like to be able to intersect each polygon with the
SpatialGridDataFrame DEM, and get an average elevation for each
polygon.  I've tried using over(), but it just hangs forever (my actual
DEM is quite large, a mosaic from SRTM).  Can anyone suggest a good way
to do this?

Thanks,
Allie

create_perp_buffers <- function(x1, y1, x2, y2, grid_dist, slice_width,
proj4string = "+init=epsg:28992
+towgs84=565.237,50.0087,465.658,-0.406857,0.350733,-1.87035,4.0812") {
xdiff <- x2 - x1
ydiff <- y2 - y1
lineangle <- atan(ydiff/xdiff)
dx = cos(lineangle)
dy = sin(lineangle)
left_angle = lineangle - 90
right_angle = lineangle + 90
midpoints = data.frame(x = seq(from = x1, to = x2, by =
dx*grid_dist), y = seq(from = y1, to = y2, by = dy*grid_dist))
begin.coord = data.frame("x" = cos(left_angle)*slice_width +
midpoints\$x, "y" = sin(left_angle)*slice_width + midpoints\$y)
end.coord = data.frame("x" = cos(right_angle)*slice_width +
midpoints\$x, "y" = sin(right_angle)*slice_width + midpoints\$y)

l <- vector("list", nrow(begin.coord))

for (i in seq_along(l)) {
l[[i]] <- Lines(list(Line(rbind(begin.coord[i,],
end.coord[i,]))), as.character(i))
}

sl <- SpatialLines(l)

names(begin.coord) <- c("begin_x", "begin_y")
names(end.coord) <- c("end_x", "end_y")
sldf <- SpatialLinesDataFrame(sl,
data.frame("lineID"=1:nrow(begin.coord), begin.coord, end.coord))
proj4string(sldf) = proj4string

blpi <- gBuffer(sldf, byid=TRUE, id=sldf\$lineID, width = grid_dist/2)

return(blpi)
}

my_blpi <- create_perp_buffers(180000, 331500, 181000, 332500, 100, 100)

data(meuse.grid)
coordinates(meuse.grid) <- ~x+y
gridded(meuse.grid) <- TRUE
proj4string(meuse.grid) <-
CRS(paste("+init=epsg:28992","+towgs84=565.237,50.0087,465.658,-0.406857,0.350733,-1.87035,4.0812"))
plot(meuse.grid, axes = T)

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20160311/228a7824/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: igeciice.png
Type: image/png
Size: 8315 bytes
Desc: not available
URL: <https://stat.ethz.ch/pipermail/r-sig-geo/attachments/20160311/228a7824/attachment.png>
```