5

I would like to upscale a 1/360 degree resolution raster (64800 x 129600) to a 1/12 resolution raster (2160 x 4320) by assigning to each cell in the low-res raster the number of cells in the high-res raster that are contained in that cell AND that have a certain value. In order to use the 'zonal' function, I first need to assign zones to the grid points of the high-res raster. Every 30x30 block of grid points should be one zone. (I.e. the zone raster should equal '1' in the top left 30x30 grid points, and 2160*4320=9331200 in the bottom right 30x30 grid points). Here is the code I am using to generate this zone raster ('raster_high' is the original high-res raster.)

zones <- raster_high;
for(i in 1:64800)
{
  j <- floor(i/30);
  zones[i,] <- rep(4320*j+1:4320*(j+1), each=30);
}

However, even only the first step i=1 takes ages. Is there a quicker way to do this?

2

What about

library(raster)

# MAke tests raster with integer values
r0 <- raster()
r0[] <- floor(rnorm(ncell(r), 30, 5))

# Define custom aggregation function
fun <- function(x, val = 28, ...) {
  length(x[x == val])
}

# Aggregate raster using custom function for spatial aggregation
r1 <- aggregate(r0, 30, fun = fun)

If you're concerned with NAs, you can use the following aggregating function instead.

fun2 <- function(x, val = 28, na.rm) {
  sum(x[x == val], na.rm = na.rm)
}

fun(c(1,1,2,3,NA), 1)
# [1] 3

fun2(c(1,1,2,3,NA), 1, na.rm = TRUE)
# [1] 2
  • (+1) aggregate is the right tool. Consider using sum(x==val) for the custom function--it's clear and more efficient. In fact, applying aggregate (with the sum function) to the indicator grid r0==val can be quite a bit faster. I'm getting very bad scaling (with both methods) for huge grids, though, even when though they fit in RAM. – whuber Jan 16 '17 at 23:05

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.