I am using MOD09Q1 now and I would like to calculating NDVI whole year data. One year, I have 46 band 1 images and 46 band 2 images.

  1. list.files

I creating list of band 1 and list of band 2.

Band1list  <- list.files(path = ".", pattern = "*.b01.tif", all.files = FALSE,
                         full.names = FALSE, recursive = FALSE,
                         ignore.case = FALSE)  

Band2list  <- list.files(path = ".", pattern = "*.b02.tif", all.files = FALSE,
                         full.names = FALSE, recursive = FALSE,
                         ignore.case = FALSE)
  1. raster::stack

Stacking all 46 band 1 (band 2) images to 'StackBand1' ('StackBand2')

StackBand1 <- stack(Band1list)

StackBand2 <- stack(Band2list)
  1. NDVI calculation via for loop

I try to stack images of same date (e.g. "2012.01.01.Smb01.tif" and "2012.01.01.Smb02.tif") and calculate NDVI

for (i in 1:length(Band1list))
  StackBand1[[i]] == i

  StackBand2[[i]] == i 

  Stack <- brick(StackBand1[[i]],StackBand2[[i]])

In here, I could not execute loop for all 46 images. It can be run all and result is showed for the last couple only.

  • 1
    Welcome to GIS SE. Please format your code for readability using edit > code sample. – Aaron Sep 19 '15 at 6:34
  • 3
    You are overwriting the variable Stack in every loop – Iris Sep 21 '15 at 9:24
  • 2
    In addition to the observation made by @Iris , you aren't really calculating anything in your loop. No NDVI calculations or anything. – Mikkel Lydholm Rasmussen Sep 21 '15 at 9:52
  • Why not use their ready NDVI or EVI products? modis.gsfc.nasa.gov/data/dataprod/… – reima Sep 21 '15 at 11:09
  • Thank you for your comment. I am try to use 250 met for spatial resolution and 8days for temporal resolution – Bang Tam Oct 26 '15 at 8:23

There's no need for a for (or any other kind of) loop since such stuff is entirely included in the raster package. If you want to calculate the NDVI, then basically all you need is

## required package

## calculate ndvi from red (band 1) and near-infrared (band 2) channel
ndvi <- overlay(StackBand1, StackBand2, fun = function(x, y) {
  (y-x) / (y+x)
  • Thank you. This question refers to a file with different bands. But what if I have different file directory/folder for NIR and Red? How do you do the NDVI? – user2543 Mar 20 '18 at 1:44

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.