3

I have a netcdf file that contains 79 variables, and for each variable there are 365 bands (one for each day of the year). I want to read all bands of one variable (i.e., a raster with 365 ayers). With the 'raster' package this would work as follows:

dailyvalues <- brick('GLOBAL_2010_day.nc', varname ='WDEP_PREC')

Result is a RasterBrick with 365 layers:

> dailyvalues 
class      : RasterBrick 
dimensions : 180, 360, 64800, 365  (nrow, ncol, ncell, nlayers)
resolution : 1, 1  (x, y)
extent     : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs 
source     : GLOBAL_2010_day.nc 
names      : X2010.01.01, X2010.01.02, X2010.01.03, X2010.01.04, X2010.01.05, X2010.01.06, X2010.01.07, X2010.01.08, X2010.01.09, X2010.01.10, X2010.01.11, X2010.01.12, X2010.01.13, X2010.01.14, X2010.01.15, ... 
Date       : 2010-01-01, 2010-12-31 (min, max)
varname    : WDEP_PREC

But I haven't suceeded to do this with the 'terra' package. I tried

dailyvalues <- rast('GLOBAL_2010_day.nc')
> dailyvalues
class       : SpatRaster 
dimensions  : 180, 360, 28835  (nrow, ncol, nlyr)
resolution  : 1, 1  (x, y)
extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +no_defs 
sources     : GLOBAL_2010_day.nc:WDEP_PREC  (365 layers) 
              GLOBAL_2010_day.nc:WDEP_SOX  (365 layers) 
              GLOBAL_2010_day.nc:WDEP_OXN  (365 layers) 
              ... and 76 more source(s)
varnames    : WDEP_PREC (WDEP_PREC) 
              WDEP_SOX (WDEP_SOX) 
              WDEP_OXN (WDEP_OXN) 
              ...
names       : WDEP_PREC_1, WDEP_PREC_2, WDEP_PREC_3, WDEP_PREC_4, WDEP_PREC_5, WDEP_PREC_6, ... 
unit        :          mm,          mm,          mm,          mm,          mm,          mm, ... 
time        : 2010-01-01 18:00:00 to 2010-12-31 12:00:00

The resulting SpatRaster has 79 'sources', but what is the syntax to use if I want to extract one 'source'? Adding varname = 'WDEP_PREC' in the rast function doesn't work. I tried dailyvalues$... but that calls single layers only (as listed under 'names').

Lena
  • 311
  • 2
  • 10

2 Answers2

3

The argument to use is called subds in terra.

dailyvalues = rast("GLOBAL_2010_day.nc", subds="WDEP_PREC")
dww
  • 30,425
  • 5
  • 68
  • 111
2

You should be able to do (as @dww says)

dailyvalues = rast("GLOBAL_2010_day.nc", subds="WDEP_PREC")

Or you can make a SpatRasterDataset

x = sds("GLOBAL_2010_day.nc")

And then extract the sub-dataset you want with

r <- x["WDEP_PREC"]

or

r <- x[1]

Furthermore, you can do

r <- rast('NETCDF:"GLOBAL_2010_day.nc:WDEP_PREC":WDEP_PREC')

(for the above to work you may need to add the path to the filename, or set your working directory to where it is)

Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63
  • When I tried the function x["WDEP_PREC"] and I have multiple variables such as "WDEP_PREC" and "WDEP_PREC_2", then x["WDEP_PREC"] selects all those variables. How do I avoid selecting multiple variables? – Herman Toothrot Feb 08 '22 at 14:54
  • My guess is that these are not variables, but rather the time dimension. or in terra-speak "layers". You can subset these too, with `[[` – Robert Hijmans Feb 08 '22 at 17:12