Geo-sp
Geo-sp

Reputation: 1704

Fill the gaps in raster stack using subsequent layers

Suppose I have a raster stack that each layer has data gap in. I want to use the next two layers in the stack to fill the gaps of each layer:

library(raster)
r1 <- raster(ncol=20,nrow=20, xmn=0, xmx=20, ymn=0,ymx=20)
r1[] <- 1:20
r2 <- r3 <- r4 <- r5 <- r1
set.seed(0)
r1[sample(1:ncell(r1), size = 20)] <- NA
r2[sample(1:ncell(r2), size = 30)] <- NA
r3[sample(1:ncell(r3), size = 10)] <- NA
r4[sample(1:ncell(r4), size = 18)] <- NA
r5[sample(1:ncell(r5), size = 18)] <- NA

s <- stack(r1, r2, r3, r4, r5)

In this case r2 and r3 will be used to fill the gaps in r1 and so forth.

Upvotes: 1

Views: 280

Answers (3)

sermomon
sermomon

Reputation: 317

Here is an adaptation for those who want to fill gaps with a moving window. This is using the average of the pixel values before and after the gap. The code also makes sure to fill in the gaps in the beginning and at the end of the RasterStack:

fill_gaps <- function(IMG, WINDOW=2){
  
  # Arguments:
  # IMG: RasterStack or RasterBrick containing the datacube with gaps. The empty pixels must be NA value, not 0.
  # WINDOW: (integer) the number of images before and after the gap used to fill the gap. By default WINDOW=2.
  
  # Note that the NA values will be replaced by the mean value of the images in ''window''.
  # If your empty pixels are stored as 0 use: img <- img[img < 0] <- NA before applying the function.
  
  if (!require("raster")) install.packages("raster")
  library(raster)
  
  s <- raster::stack(IMG)
  
  for(i in 1:nlayers(s)){
    if(i <= WINDOW) {
      s[[i]] <- cover(s[[i]], mean(s[[WINDOW+1]]:s[[WINDOW*2]]), na.rm=TRUE)}
    if(i >= (nlayers(s)-WINDOW)) {
      s[[i]] <- cover(s[[i]], mean(s[[nlayers(s)-(WINDOW*2+1)]]:s[[nlayers(s)-WINDOW]]), na.rm=TRUE)}
    else {
      s[[i]] <- cover( s[[i]], mean(s[[(i-WINDOW):(i+WINDOW)]], na.rm=TRUE))}
  }

  return(s)

}

Upvotes: 1

MikeJewski
MikeJewski

Reputation: 357

Not sure if this is what you want, but it will give you a start. I'm still new to R, so there is probably another way to do it.

library(raster)
r1 <- raster(ncol=20,nrow=20, xmn=0, xmx=20, ymn=0,ymx=20)
r1[] <- 1:20
r2 <- r3 <- r4 <- r5 <- r1
set.seed(0)
r1[sample(1:ncell(r1), size = 20)] <- NA
r2[sample(1:ncell(r2), size = 30)] <- NA
r3[sample(1:ncell(r3), size = 10)] <- NA
r4[sample(1:ncell(r4), size = 18)] <- NA
r5[sample(1:ncell(r5), size = 18)] <- NA

s <- stack(r1, r2, r3, r4, r5)

for(i in 1:(nlayers(s) - 2) ){
    s[[i]] <- merge( s[[i]], mask( s[[(i+1)]], s[[i]], inverse = TRUE))
    s[[i]] <- merge( s[[i]], mask( s[[(i+2)]], s[[i]], inverse = TRUE))
}

Upvotes: 2

Robert Hijmans
Robert Hijmans

Reputation: 47536

MikeJewski's solution might work, but the cover function is designed for this and more direct. It is not clear how you want to use the next two layers. The mean:

for(i in 1:(nlayers(s) - 2) ){
    s[[i]] <- cover( s[[i]], mean( s[[(i+1):(i+2)]], na.rm=TRUE))
}

Or first the closest (as MikeJewski assumed):

for(i in 1:(nlayers(s) - 2) ){
    s[[i]] <- cover( s[[i]], s[[(i+1)]])
    s[[i]] <- cover( s[[i]], s[[(i+2)]])
}

This would be another, but probably inefficient, approach:

f <- function(x) {
    for(i in 1:((ncol(x)-2)) ){
        x[is.na(x[,i]),i]  <- x[is.na(x[,i]),i+1]  
        x[is.na(x[,i]),i]  <- x[is.na(x[,i]),i+2]  
    }
    x
}

ss <- calc(s, f)

Upvotes: 2

Related Questions