0

I want to get proper GPM IMERG data with variables including time in R and work with it. I used "https://disc.gsfc.nasa.gov/" and got the link list in txt file. I don't know what to do from here. Could someone give an idea?

2
  • First, download the data. IMERG is distributed with different time resolutions. Each file contains one timestep, so isn't a multidimensional time file. Use terra library to open the file in R
    – aldo_tapia
    Sep 16 at 13:19
  • Could you help me walk through this? I have downloaded the nc4 files using wget, 265 in total. I'm unable to use terra maybe because I'm still quite new and unfamiliar with using new pacakges. my idea is to combine these files into a single so I can compare with different sets of observations.
    – hat6ytrs
    Sep 17 at 16:25

1 Answer 1

0

Terra can open IMERG rasters, but these files don't have CRS and are flipped, so you must correct them before computing or implementing any analysis. In this case, I list the files inside a folder, then I iterate through files opening, correcting the X Y dimension (transposing the file), adding the CRS, and saving the result to a list):

library(terra)

fnames <- list.files(pattern = 'IMERG')

rlist <- list()

for (i in seq_along(fnames)) {
  r = rast(fnames[i])
  r = trans(r)
  ext(r) = c(-180,180,-90,90)
  crs(r) = "epsg:4326"
  rlist[[i]] = r
}

Then, you can perform any arithmetic operation or zonal stats or whatever you want to achieve.

If you need only one band, the process is faster:

rlist <- list()

for (i in seq_along(fnames)) {
  r = rast(fnames[i])['precipitationCal']
  r = trans(r)
  ext(r) = c(-180,180,-90,90)
  crs(r) = "epsg:4326"
  rlist[[i]] = r
}

plot(app(rast(rlist), sum), type = "interval",
     breaks = c(0.1,5,10,20,50,100))

enter image description here

I also post the same process with python using a custom function:

import rioxarray as rioxr

def load_imerg(file, var):
  da = rioxr.open_rasterio(file, engine='h5netcdf')
  da = da[0][var]
  return da.assign_coords({"x": (da.x / 10) - 90}) \
      .assign_coords({"y": (da.y / 10) - 180}) \
      .transpose('band', 'x', 'y') \
      .rio.write_crs(4326) \
      .rename({'x': 'y', 'y': 'x'})

You can loop using this function saving the results to a list as well.

I'm working with 30" product, so I found Python solution faster than R for doing the same task over thousands files.

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.