0

I have maximum daily temperature data (raster format) from E-OBS for the period 1980-2019. I have reduced this to 1985-2014 and only kept the months April-October (30 years, 214 days). Now I want to calculate a 90th percentile rasterlayer for each day from 1st May until 30th September (153 days). So for each day calculate the 90th percentile across the 30 years and with a 15-day centered window.

I have called the rasterstack RP:

RP <- setZ(RP, as.Date(names(RP), format = "X%Y.%m.%d")) # Set the time dimension using the rasterlayer names
RP <- subset(RP, which(getZ(RP) >= '1985-01-01' & (getZ(RP) <= '2014-12-31'))) # Keep only 1985-2014
RP <- subset(RP, which(substr(getZ(RP),6,7) %in% c('04','05','06','07','08','09','10'))) # Keep only April-October

I have used a foreach loop, because the for loop took too much time.

P90 <- stack() # Create empty rasterstack for the 90th percentile layers

registerDoParallel(makeCluster(detectCores()-1))

P90 <- foreach (d = 1:153, .packages = "raster", .combine = addLayer, .init = P90, .inorder = TRUE) %dopar% {
  Raster_1 <- stack() # Create empty rasterstack to store all the layers needed for the calculation (30 years x 15 days = 450 layers)
  
  for (i in 0:29) { # Calculation is done for all 30 years
    
    for (j in -7:7) { # Moving window of 15 days is used, 7 days before and 7 days after the chosen day
      Raster_1 <- addLayer(Raster_1, RP[[30 + d + i*214 + j]])
    }
  }
  
  P90 <<- addLayer(P90, calc(x = Raster_1, fun = function(x) quantile(x, probs = c(0.90), na.rm = TRUE)), forceapply = TRUE) # 90th percentile calculation
}

stopCluster()

I am getting the following error:

Error in { : task 1 failed - "cannot use this function"

Could anyone please tell me what I am doing wrong?

The strange thing is also, using the exact same code on a different dataset (not from E-OBS, but from the ESGF data portal) the code works perfectly fine and I receive a P90 rasterstack with 153 layers.

0

I have found a solution that works. The problem was not in the foreach loop but in the quantile calculation itself.

The following code gives me exactly what I needed:

P90 <- stack() # Create empty rasterstack for the 90th percentile layers

registerDoParallel(makeCluster(detectCores()-1))

P90 <- foreach (d = 1:153, .packages = "raster", .combine = addLayer, .init = P90, .inorder = TRUE) %dopar% {
  Raster_1 <- stack() # Create empty rasterstack to store all the layers needed for the calculation (30 years x 15 days = 450 layers)
  
  for (i in 0:29) { # Calculation is done for the 30 years of the reference period
    
    for (j in -7:7) { # Moving window of 15 days is used, 7 days before and 7 days after the chosen day
      Raster_1 <- addLayer(Raster_1, RP[[30 + d + i*214 + j]])
    }
  }
  
  Q90 <- function(x) {quantile(x, probs = 0.90, na.rm = TRUE)}
  P90 <<- addLayer(P90, calc(Raster_1, Q90)) 

}
stopCluster()

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.