R help - Jun 2018 - rasterize SpatialPolygon object using a RasterBrick object

I am trying to rasterize a SpatialPolygon object by a RasterBrick object.
The documentation of the raster::rasterize function explicitly says this is
allowed. Here's what I am doing

# load the raster package
library("raster")
# create a raster brick object using the example from the brick
function documentation
b <- brick(system.file("external/rlogo.grd",
package="raster"))
b
nlayers(b)
names(b)
extract(b, 870)
# create a SpatialPolygon object using the example from the function
documentation
Sr1 = Polygon(10*cbind(c(2,4,4,1,2),10*c(2,3,5,4,2)))
Sr2 = Polygon(10*cbind(c(5,4,2,5),10*c(2,3,2,2)))
Sr3 = Polygon(10*cbind(c(4,4,5,10,4),c(5,3,2,5,5)))
Sr4 = Polygon(10*cbind(c(5,6,6,5,5),c(4,4,3,3,4)), hole = TRUE)

Srs1 = Polygons(list(Sr1), "s1")
Srs2 = Polygons(list(Sr2), "s2")
Srs3 = Polygons(list(Sr3, Sr4), "s3/4")
SpP = SpatialPolygons(list(Srs1,Srs2,Srs3), 1:3)

plot(b[[1]])
plot(SpP, add = T)
# crop
clip1 = crop(b, extent(SpP))
# rasterize returns an error, but documentation says it should return
a RasterBrick object
clip2 = rasterize(SpP, b, mask = T)
Error in v[, r] <- rrv :
     number of items to replace is not a multiple of replacement length
# however, if I used only one layer, all would be fine
clip2 = rasterize(SpP, b[[1]], mask = T)

Of course, I could loop over the brick's layers, but as I understand it,
that would defeat the purpose of a brick object.

I want to use clip2 to then get the histogram of pixel values in the
layers, like this:

vals = getValues(clip2)

Can anyone tell me why I am getting this error, and how to go around it
efficiently?

	[[alternative HTML version deleted]]

I see this problem in 2.6-7 (version on CRAN) but it's now fixed in dev on
RForge, you can try it out by installing from there, or from the Github
mirror with devtools::install_github("rforge/raster/pkg/raster").

There's an imminent release to CRAN some time soon.

Cheers, Mike.


On Sat, 2 Jun 2018 at 04:48 Dimitri Szerman <dimitrijoe at gmail.com>
wrote:

I am trying to rasterize a SpatialPolygon object by a RasterBrick object.
The documentation of the raster::rasterize function explicitly says this is
allowed. Here's what I am doing

# load the raster package
library("raster")
# create a raster brick object using the example from the brick
function documentation
b <- brick(system.file("external/rlogo.grd",
package="raster"))
b
nlayers(b)
names(b)
extract(b, 870)
# create a SpatialPolygon object using the example from the function
documentation
Sr1 = Polygon(10*cbind(c(2,4,4,1,2),10*c(2,3,5,4,2)))
Sr2 = Polygon(10*cbind(c(5,4,2,5),10*c(2,3,2,2)))
Sr3 = Polygon(10*cbind(c(4,4,5,10,4),c(5,3,2,5,5)))
Sr4 = Polygon(10*cbind(c(5,6,6,5,5),c(4,4,3,3,4)), hole = TRUE)

Srs1 = Polygons(list(Sr1), "s1")
Srs2 = Polygons(list(Sr2), "s2")
Srs3 = Polygons(list(Sr3, Sr4), "s3/4")
SpP = SpatialPolygons(list(Srs1,Srs2,Srs3), 1:3)

plot(b[[1]])
plot(SpP, add = T)
# crop
clip1 = crop(b, extent(SpP))
# rasterize returns an error, but documentation says it should return
a RasterBrick object
clip2 = rasterize(SpP, b, mask = T)
Error in v[, r] <- rrv :
number of items to replace is not a multiple of replacement length
# however, if I used only one layer, all would be fine
clip2 = rasterize(SpP, b[[1]], mask = T)

Of course, I could loop over the brick's layers, but as I understand it,
that would defeat the purpose of a brick object.

I want to use clip2 to then get the histogram of pixel values in the
layers, like this:

vals = getValues(clip2)

Can anyone tell me why I am getting this error, and how to go around it
efficiently?

[[alternative HTML version deleted]]

______________________________________________
R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
<http://www.r-project.org/posting-guide.html>
and provide commented, minimal, self-contained, reproducible code.
-- 
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

	[[alternative HTML version deleted]]