1

Here is a list of files I would like to create a loop with for a raster calculation

"J:\NAIP\temp5\m_4208902_nw_16_1m3m_sd.tif"
"J:\NAIP\temp5\m_4208902_nw_16_3m_sd.tif"
"J:\NAIP\temp5\m_4208902_se_16_1m3m_sd.tif"
"J:\NAIP\temp5\m_4208902_se_16_3m_sd.tif"
"J:\NAIP\temp5\m_4208902_sw_16_1m3m_sd.tif"
"J:\NAIP\temp5\m_4208902_sw_16_3m_sd.tif"

Here is my R code that is unable to pair the correct rasters together for the calculation. An error message of mismatched origins is returned. The calculation runs successfully on one image pair when the patterns = the full intended file name pair.

library(raster)
library(rgdal)

#path to directory with folders
path <- "J:/NAIP/temp5/"

#get file names using list.files() function
stack_list <- list.files(path, pattern="16_3m_sd.tif$", full.names=T)
stack_list3 <- list.files(path, pattern="16_1m3m_sd.tif$", full.names=T)

#Loop over every file in stack_list
for(i in 1:length(stack_list)) {  
#stack image i
img  <- stack(stack_list[i])

#Loop over every file in stack_list3
for(j in 1:length(stack_list3)) {  

#stack image j
img3  <- stack(stack_list3[j])

#calc function for image i and j
calc <- (((img-img3) / (img+img3))*100 + 100) 
       
outname <- sub(pattern     = ".tif",
        replacement = "_ndti_1m3m.tif", 
        x           = stack_list[i])

#export raster
writeRaster(calc, 
            filename  = outname,
            datatype='FLT8S', # FLT8S =float
            overwrite = T)
}
}      

Output:

Error in compareRaster(e1, e2, extent = FALSE, rowcol = FALSE, crs = TRUE,  : 
  different origin
> stack_list
[1] "J:/NAIP/temp5/m_4208902_nw_16_3m_sd.tif"
[2] "J:/NAIP/temp5/m_4208902_se_16_3m_sd.tif"
[3] "J:/NAIP/temp5/m_4208902_sw_16_3m_sd.tif"
> stack_list3
[1] "J:/NAIP/temp5/m_4208902_nw_16_1m3m_sd.tif"
[2] "J:/NAIP/temp5/m_4208902_se_16_1m3m_sd.tif"
[3] "J:/NAIP/temp5/m_4208902_sw_16_1m3m_sd.tif"
>

Sample rasters: https://drive.google.com/file/d/1qMiRluXZzZ-EDx1yj9_EkS9eWke1xNWW/view?usp=sharing

2
  • 1
    So your rasters have a different origin? How can we help with that? Is it comparing the right pair of rasters? Can we get your data? Have you looked at the pair that fail to see what's going on?
    – Spacedman
    Commented Sep 30, 2021 at 17:50
  • Sample raster download has been linked above. There are 6 rasters. 3 pairs with the same origin/extent. The code isn't making the calculation in the correct order as the are listed alphabetically.
    – mdn
    Commented Sep 30, 2021 at 18:03

1 Answer 1

2

You are looping over all values of j for all values of i:

for(i in 1:length(stack_list)) {  
#stack image i
img  <- stack(stack_list[i])

#Loop over every file in stack_list3
for(j in 1:length(stack_list3)) {  

#stack image j
img3  <- stack(stack_list3[j])

which means you'll go on to compute calc for all 9 possible pairs of stacks.

What you want (I think) is to not have the j loop at all, and get the i element of stack_list3 to go with the i element of stack_list:

img3  <- stack(stack_list3[i])

To elaborate, lets see what your loops do when looped over the names. This is your basic structure:

for(i in 1:length(stack_list)){
 iname = stack_list[i]
 for(j in 1:length(stack_list3)){
   jname = stack_list3[j]

   message("iname = ",iname, " jname = ",jname)
   # now you do something with the corresponding data...
} 
} 

which produces this list of pairs:

iname = ./m_4208902_nw_16_3m_sd.tif jname = ./m_4208902_nw_16_1m3m_sd.tif
iname = ./m_4208902_nw_16_3m_sd.tif jname = ./m_4208902_se_16_1m3m_sd.tif
iname = ./m_4208902_nw_16_3m_sd.tif jname = ./m_4208902_sw_16_1m3m_sd.tif
iname = ./m_4208902_se_16_3m_sd.tif jname = ./m_4208902_nw_16_1m3m_sd.tif
iname = ./m_4208902_se_16_3m_sd.tif jname = ./m_4208902_se_16_1m3m_sd.tif
iname = ./m_4208902_se_16_3m_sd.tif jname = ./m_4208902_sw_16_1m3m_sd.tif
iname = ./m_4208902_sw_16_3m_sd.tif jname = ./m_4208902_nw_16_1m3m_sd.tif
iname = ./m_4208902_sw_16_3m_sd.tif jname = ./m_4208902_se_16_1m3m_sd.tif
iname = ./m_4208902_sw_16_3m_sd.tif jname = ./m_4208902_sw_16_1m3m_sd.tif

ie you are going to process each member of the first list with each member of the second list.

If you want to process element 1 of the first list with element 1 of the second list, and element 2 of the first list with element 2 of the second list and so on, then you only need a single loop and to get the i-th element of each list within the single loop:

for(i in 1:length(stack_list)){
   iname = stack_list[i]
   jname = stack_list3[i] # note `i` here
   message("iname = ",iname, " jname = ",jname)
   # do something with those files...
 }

which loops over three corresponding pairs to work on:

iname = ./m_4208902_nw_16_3m_sd.tif jname = ./m_4208902_nw_16_1m3m_sd.tif
iname = ./m_4208902_se_16_3m_sd.tif jname = ./m_4208902_se_16_1m3m_sd.tif
iname = ./m_4208902_sw_16_3m_sd.tif jname = ./m_4208902_sw_16_1m3m_sd.tif

note this depends on the folder containing only matching pairs of files. If the file naming scheme is as you give it then they should sort into the same order since the only difference is the extra 1m in the filename. You could use sort on both file lists to ensure they are in the same order, and check the lengths of each to make sure they are the same depending on how confident you are on your folder having the correct data in.

1
  • 1
    Spacedman - incredible. Thank you so much for your detailed walk through!
    – mdn
    Commented Oct 1, 2021 at 15:49

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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