summing rasters with a condition given by other rasters
Hi Martin,
On May 10, 2015, at 4:42 PM, Martin Brandt <martin.brandt at mailbox.org> wrote:
Hi Ben, thanks again, I see, but i'm not sure if i understand it correctly. The thing is that each pixel in the rasters sos and eos has different values, i.e. the start of season and end of season is not the same for each pixel. Thus, I cannot really define eos and sos but need to use the given rasters containing the values..
Oooooo. I think I see what you are after. I'll bet there are better ways, but my brute force and ignorance approach would be to prepend your start-of-season raster (sos) and end-of-season raster (eos) to your raster data. I would still use calc but modify the sum_segment function. After the computation I would then drop the sos and eos layer form the data. Note that sos and eos still have to be numeric indices into the layers of the data, b.
library(raster)
# create a multiple layer brick
b <- brick(system.file("external/rlogo.grd", package="raster"))
# add layers
b <- addLayer(b,b,b)
b
# class : RasterStack
# dimensions : 77, 101, 7777, 9 (nrow, ncol, ncell, nlayers)
# resolution : 1, 1 (x, y)
# extent : 0, 101, 0, 77 (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=merc
# names : red.1, green.1, blue.1, red.2, green.2, blue.2, red.3, green.3, blue.3
# min values : 0, 0, 0, 0, 0, 0, 0, 0, 0
# max values : 255, 255, 255, 255, 255, 255, 255, 255, 255
# define the sos and eos by indexed position
# sos will range from 1 to 4
sos <- raster(matrix(sample(1:4, ncell(b), replace = TRUE), nrow = nrow(b), ncol = ncol(b)),
template = b)
#eos will range from 6 to 9
eos <- raster(matrix(sample(6:nlayers(b), ncell(b), replace = TRUE), nrow = nrow(b), ncol = ncol(b)),
template = b)
# prepend sos and eos to your data
b <- addLayer(sos, eos, b)
# define the summing function - add 2 to account for sos and eos at the start
sum_segment <- function(x, ...) {
sum(x[(x[1]+2):(x[2] + 2)],...)
}
# run the computation
s <- calc(b, sum_segment)
#restore b
b <- dropLayer(b, c(1,2))
s
#class : RasterLayer
#dimensions : 77, 101, 7777 (nrow, ncol, ncell)
#resolution : 1, 1 (x, y)
#extent : 0, 101, 0, 77 (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=merc
#data source : in memory
#names : layer
#values : 0, 2295 (min, max)
b
#class : RasterStack
#dimensions : 77, 101, 7777, 9 (nrow, ncol, ncell, nlayers)
#resolution : 1, 1 (x, y)
#extent : 0, 101, 0, 77 (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=merc
#names : red.1, green.1, blue.1, red.2, green.2, blue.2, red.3, green.3, blue.3
#min values : 0, 0, 0, 0, 0, 0, 0, 0, 0
#max values : 255, 255, 255, 255, 255, 255, 255, 255, 255
You might be able to do something similar with the overlay(). If that doesn't do it then we'll both have to wait for someone to come to our rescue!
Cheers,
Ben
cheers, Martin -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/summing-rasters-with-a-condition-given-by-other-rasters-tp7588222p7588233.html Sent from the R-sig-geo mailing list archive at Nabble.com.
_______________________________________________ R-sig-Geo mailing list R-sig-Geo at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Ben Tupper Bigelow Laboratory for Ocean Sciences 60 Bigelow Drive, P.O. Box 380 East Boothbay, Maine 04544 http://www.bigelow.org