rast function error with MOD09GA hdf images in R

Question

I have some MODIS images downloaded. When I try to create a SpatRaster from the hdf files, using terra rast function, it works perfectly for "MOD09A1" but it doesn´t work for "MOD09GA".

terra::rast("C:/Users/User/AppData/Local/Temp/_modis/MOD09GA.A2011025.h08v06.006.2015216102410.hdf")

Error: [rast] number of rows and/or columns do not match

What is the problem? Is there any other function I could use? Thanks!

Answer 1

The problem is that the file has subdatasets with different resolutions.

To get the file your are using

# remotes::install_github("rspatial/luna")
aoi <- c(-106, -105, 26, 27)
f <- luna::getModis("MOD09GA", "2011-01-25", "2011-01-26", aoi, download=TRUE, path=".", user="*", password="*")
f 
[1] "./MOD09GA.A2011025.h08v06.006.2015216102410.hdf"

To see the subdatasets:

library(terra)
describe_sds(f)
id   name                                                          desc                                                               nrow ncol nlyr
 1   HDF4_EOS:EOS_GRID:...MODIS_Grid_1km_2D:num_observations_1km   [1200x1200] num_observations_1km MODIS_Grid_1km_2D (8-bit integer) 1200 1200    1
(...)
12        HDF4_EOS:EOS_GRID:...MODIS_Grid_500m_2D:sur_refl_b01_1   [2400x2400] sur_refl_b01_1 MODIS_Grid_500m_2D (16-bit integer) 2400 2400    1
(...)

So you need to access the different subdatasets seperately, like this

b1 <- rast(f, 12)
b2 <- rast(f, 13)
b1
# class       : SpatRaster 
# dimensions  : 2400, 2400, 1  (nrow, ncol, nlyr)
# resolution  : 463.3127, 463.3127  (x, y)
# extent      : -11119505, -10007555, 2223901, 3335852  (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=sinu +lon_0=0 +x_0=0 +y_0=0 +R=6371007.181 +units=m +no_defs 
# data source : MOD09GA.A2011025.h08v06.006.2015216102410.hdf:MODIS_Grid_500m_2D:sur_refl_b01_1 
# names       : sur_refl_b01_1 

Perhaps followed by

bb <-c(b1, b2)

Or create a SpatDataSet like this

x <- sds(f, 12:22)

In the development version of terra you can then proceed and do

y <- collapse(x)

or something like the below to get a particular set of subdatasets (that have the same spatial resolution)

r  <- rast(f, 12:18)
r
#class       : SpatRaster 
#dimensions  : 2400, 2400, 7  (nrow, ncol, nlyr)
#resolution  : 463.3127, 463.3127  (x, y)
#extent      : -11119505, -10007555, 2223901, 3335852  (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=sinu +lon_0=0 +x_0=0 +y_0=0 +R=6371007.181 +units=m +no_defs 
#source(s)   : MOD09GA.A2011025.h08v06.006.2015216102410.hdf:MODIS_Grid_500m_2D:sur_refl_b01_1  
          MOD09GA.A2011025.h08v06.006.2015216102410.hdf:MODIS_Grid_500m_2D:sur_refl_b02_1  
          MOD09GA.A2011025.h08v06.006.2015216102410.hdf:MODIS_Grid_500m_2D:sur_refl_b03_1  
          ... and 4 more source(s)
#names       : sur_refl_b01_1, sur_refl_b02_1, sur_refl_b03_1, sur_refl_b04_1, sur_refl_b05_1, sur_refl_b06_1, ...