1

I am working with series of raster data.

I am checking temporal outlier of each pixel by using image before and after, then replace current pixel value (t) by average of image before (t-1) and after (t+1). I am able to run with test of 2 rows x 2 columns, but with the real data (NA value in background of image). It can be run but very slow. It cost more than one week to complete one image (979x526 pixels)

I am a beginner R user.

Following is the script

# rm(list=ls())
# setwd ("D:\\5Validation\\92NDVIValid\\GapTestNDVI")
# I've commented these out because these could catch people out

Files <- list.files(pattern="*.tif") 

library(raster)
library(stats)

#### Need to difine variables in loop
Ndvit       <- as.matrix(raster(Files[2]))                # The second image

Ndvit_1     <- as.matrix(raster(Files[1]))                 # The first image

Ndvit_3     <- as.matrix(raster(Files[3]))                 # The third image

Con1        <- (Ndvit-Ndvit_1) <(-0.01*Ndvit_1)           # Condition 1  
Con2        <- (Ndvit-Ndvit_3) <(-0.01*Ndvit_3)           # Condition 2

Con         <-  Con1 & Con2                               # Conbine both Condition 1&2 (T==1, F==0)
summary(Con1)

## Dimension of image in series
NRow       <- nrow(raster(Files[1]))

NCol       <- ncol(raster(Files[1]))

#Create blank Table to merge all output  
NdviReli  <- matrix(nrow = 979, ncol = 526)

Out       <- matrix(nrow = 979, ncol = 526)


## Need image before and after to fill, 
## So it will process from second image to lenghth of image - 1)

for (t in 2:(length(Files)-1)) {   
  for (i in 1:NRow) { 
    for (j in 1:NCol) { 
        x <- t-1         ## Set variabe t-1
        y <- t+1         ## ## Set variabe t+1

        ## Assigned input image as matrix
        Ndvit       <- as.matrix(raster(Files[t]))
        Ndvit_x     <- as.matrix(raster(Files[x]))
        Ndvit_y     <- as.matrix(raster(Files[y]))

        ## Condition 
        Con1     <- (Ndvit-Ndvit_x) <(-0.01*Ndvit_x) 
        Con2     <- (Ndvit-Ndvit_y) <(-0.01*Ndvit_y)  
        Con      <- Con1 & Con2
        if(is.na(Con[i,j])) {

           NdviReli[i,j] <- NA

          } else if  (Con[i,j] == 1) {
          ## If Condition satify
          NdviReli[i,j] <- (Ndvit_x[i,j] + Ndvit_y[i,j])/2

          } else  {
          ## If Condition satify
          NdviReli[i,j] <- Ndvit[i,j]
          }

          Reli <- print(paste(NdviReli[i,j])) ## Using this if want to check the value one by one
          Out <- raster(NdviReli)                                       ## Assign matrix as raster
          Ex       <- extent(c(487124.6, 618624.6, 1902496, 2147246))   ## Assign coordination, projection base on image
          extent(Out) <- Ex    
          projection(Out) <- CRS("+proj=utm +zone=47 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0")
          writeRaster(Out,filename= paste("Reli",Files[t]),format="GTiff", overwrite=T)
    } 
  }
}
  • It's the three nested for loops. Can you provide a minimal reproducible example (i.e. an example image) that we can test alternatives with? – Phil Sep 23 '16 at 9:37
  • @ Midavalo: Thank you very much. Please, help me on this issue, I am new R user. I am sharing link of sample data as below. drive.google.com/drive/folders/… – Bang Tam Sep 24 '16 at 14:55
  • The Google Drive link has restricted permissions; I can't download any files. – Phil Sep 26 '16 at 20:50
  • @ Midavalo: Could you please to check the the link as below 4shared.com/folder/O6XyaGoE/TestingData.html – Bang Tam Sep 28 '16 at 6:59
  • It says I need to sign up for an account to download the files, which I'm not going to do. You could just change the permissions on your Google Drive version to public. – Phil Sep 28 '16 at 8:59

Your Answer

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

Browse other questions tagged or ask your own question.