Home > Software engineering >  8 Day MODIS Raster to Monthly sum in R
8 Day MODIS Raster to Monthly sum in R

Time:05-17

I have MODIS data which contains 8 days product. Some of month has 3 images and some of month has four imgaes. So how I can build a R scrip to calculate monthly sum from this data

I am using the MODIS MOD17A2HGF GPP data set which comes as a 500x500 meter grid of the 8 day sum of GPP. The first file ends in _001 which means that it starts on January 1st and goes until the 8th. The next file, _009, then covers the 8th to the 16th of January and on and on. enter image description here

Now I want to sum GPP for months(January to December).

In this folder, it have 920 files (20 years) which every year has 46 files.

In ArcGIS, I know how to sum with raster calculator.But it will do so much repetitive works if I sum every month.

enter image description here

In R, I tried to do it just like this

library(terra)
a <- rast('MOD17A2HGF_GPP_2001_001.tif')
b <- rast('MOD17A2HGF_GPP_2001_009.tif')
c <- rast('MOD17A2HGF_GPP_2001_017.tif')
d <- rast('MOD17A2HGF_GPP_2001_025.tif')
GPP_January <- a b c d
writeRaster(GPP_January,'F:/test/GPP_January_2001.tif')

It will also do so much repetitive works if I sum every month.I don't know how to select the specified files of month in R loop.

I know how to sum annual GPP

library(tidyverse)
library(terra)

GPP2001 <- list.files(pattern = '2001',full.names = T) %>% 
  rast() %>% sum()

But how can I sum monthly GPP per year in R or ArcGIS?

CodePudding user response:

Get all files/layers

library(terra)
ff <- list.files(pattern="\\.tif$", full.names = T) 
r <- rast(ff) 

For each file, extract the year and the day number

# ff <- c('MOD17A2HGF_GPP_2001_001.tif', 'MOD17A2HGF_GPP_2001_009.tif')
s <- gsub('MOD17A2HGF_GPP_', "", basename(ff))
s <- gsub('.tif', "", s)
s <- strsplit(s, "_")
s <- do.call(rbind, s)
s
#     [,1]   [,2] 
#[1,] "2001" "001"
#[2,] "2001" "009"

Get the month and combine with year

d <- as.Date(as.numeric(s[,2]), origin = paste(s[,1], "-01-01", sep = "")) - 1
m <- format(d, "%m")
ym <- paste0(s[,1], m)
ym
#[1] "200101" "200101"

Now you can calculate the monthly average like this

x <- tapp(x, ym, mean)

(to be more precise, the average for the dates in each month; that is not exactly the monthly average, but probably good enough).

  • Related