Skip to main content

Working with spatial objects

loading packages

library(rgdal)
library(sp)
library(raster)

create raster

df <- data.frame( x = rep( 0:4, each=5 ),y = rep( 0:4,  5),l = 1:25)

coordinates(df) <- ~ x + y # coerce to spatial points data frame by assigning x and y values as coordinates

gridded(df) <- TRUE # coerce to SpatialPixelsDataFrame
r <- raster(df) # coerce to raster

create spatial polygon

coords = matrix(c(2,1.5,  2,3.5,   3.5,3.5,   3.5,1.5,  2,1.5 ), ncol = 2, byrow = TRUE) # coordinates

p1 = Polygon(coords) 
pl = SpatialPolygons(list(Polygons(list(p1), ID = "a")))

plotting

plot(r)
text(r)
plot(pl, add = TRUE)

Simple operations:

taking spatial average of raster values inside polygon

r.l = extract(r,pl,weights = T) # taking area weights
r.l
## [[1]]
##      value    weight
## [1,]    14 0.1666667
## [2,]    19 0.3333333
## [3,]    13 0.1666667
## [4,]    18 0.3333333
v = sum(r.l[[1]][,1] * r.l[[1]][,2])
v
## [1] 16.83333
r.l2 = extract(r,pl)
r.l2
## [[1]]
## [1] 14 19 13 18
lapply(r.l2,FUN = mean) # without taking area weights
## [[1]]
## [1] 16

Comments