I want to extract a variable called NVEL from the netCDF file which has five dimensions (i, j, tile, k, time) where i is longitude, j is latitude, k is the level of depths I want to extract NVEL(i, j, tile=3, k=1st level, time) the input file can be downloaded from here https://drive.google.com/file/d/12NQp_uLr_IZLLU6Fzr555gKGGJlrRE4H/view?usp=sharing
NVEL<- brick("NVEL_1992_01.nc", varname= "NVEL", lvar=1, nl=1)
NVEL <- NVEL[[which(getZ(NVEL) == 3)]]
This does not work. How to deal with a variable of 5 dimensions?
I see that this returns 50 (k) * 13 (tiles) * 1 (time) = 650 layers
library(terra)
f <- "NVEL_1992_01.nc"
x <- rast(f)
x
#class : SpatRaster
#dimensions : 90, 90, 650 (nrow, ncol, nlyr)
#resolution : 1, 1 (x, y)
#extent : -0.5, 89.5, -0.5, 89.5 (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=longlat +datum=WGS84 +no_defs
#data source : NVEL_1992_01.nc
#names : NVE_1, NVE_2, NVE_3, NVE_4, NVE_5, NVE_6, ...
The order is k-wise (and tile-wise within tiles). See (the rather lengthy) output from
terra::describe(f)
You can extract that information like this:
d <- describe(f, print=FALSE)
d <- unlist(strsplit(d, "\n"))
i <- grep("NETCDF_DIM_k=", d)
j <- grep("NETCDF_DIM_tile=", d)
k <- sapply(strsplit(d[i], "="), function(x) x[2])
tile <- sapply(strsplit(d[j], "="), function(x) x[2])
kt <- paste0("k", k, "_tile", tile)
names(x) <- kt
x
#class : SpatRaster
#dimensions : 90, 90, 650 (nrow, ncol, nlyr)
#resolution : 1, 1 (x, y)
#extent : -0.5, 89.5, -0.5, 89.5 (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=longlat +datum=WGS84 +no_defs
#data source : NVEL_1992_01.nc
#names : k0_tile0, k0_tile1, k0_tile2, k0_tile3, k0_tile4, k0_tile5, ...
This should happen automatigically in a future version. You can continue with terra (very similar to raster) or take the data back to a RasterBrick by doing
b <- brick(x*1)
(multiplying to get the values out of the file)