2

I have multiple variables in a timeseries which I want to save as NetCDF. But when I combine all variables with terra::sds() the time information disappears. I tried the following:

library("terra") ts <- c("20210101", "20210301", "20210501", "20210701", "20210901", "20211101") ts_dt <- as.Date(ts, "%Y%m%d") r1 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(r1) <- 1:ncell(r1) r2 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(r2) <- 1:ncell(r2) r3 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(r3) <- 1:ncell(r3) r4 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(r4) <- 1:ncell(r4) r5 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(r5) <- 1:ncell(r5) r6 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(r6) <- 1:ncell(r6) r_p <- rast(sds(r1, r2, r3, r4, r5, r6)) names(r_p) <- ts_dt terra::time(r_p) <- ts_dt 

Everything is fine when I check the SpatRaster information:

class : SpatRaster
dimensions : 18, 36, 6 (nrow, ncol, nlyr)
resolution : 10, 10 (x, y)
extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84
sources : memory
memory
memory
... and 3 more source(s)
names : 2021-01-01, 2021-03-01, 2021-05-01, 2021-07-01, 2021-09-01, 2021-11-01
min values : 1, 1, 1, 1, 1, 1
max values : 648, 648, 648, 648, 648, 648
time : 2021-01-01 to 2021-11-01

So I create another variable:

s1 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(s1) <- 1:ncell(s1) s2 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(s2) <- 1:ncell(s2) s3 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(s3) <- 1:ncell(s3) s4 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(s4) <- 1:ncell(s4) s5 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(s5) <- 1:ncell(s5) s6 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(s6) <- 1:ncell(s6) s_p <- rast(sds(s1, s2, s3, s4, s5, s6)) names(s_p) <- ts_dt terra::time(s_p) <- ts_dt 

In a last step I combine both variables with sds():

p <- sds(r_p, s_p) 

But p hasn't any time variable:

class : SpatRasterDataset
subdatasets : 2
dimensions : 18, 36 (nrow, ncol)
nlyr : 6, 6
resolution : 10, 10 (x, y)
extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
coord. ref. : lon/lat WGS 84
source(s) : memory
names : Mean Temp, Mean Prec

When I write p to NetCDF, the resulting .nc file, as expected, hasn't any time information. So my question is as the headline says: How can I write multiple timeseries variables to one NetCDF file while preserving information about time?

1 Answer 1

4

A SpatRasterDataset does not have a time dimension. The SpatRasters it contains can have that dimension, and it may be different for each SpatRaster.

Here is a simplified script that suggests that all is good (you did not show what you get, so I cannot compare)

library("terra") #terra 1.5.36 ts <- c("20210101", "20210301", "20210501", "20210701", "20210901", "20211101") ts_dt <- as.Date(ts, "%Y%m%d") r1 <- rast(xmin=-180, xmax=180,ymin=-90, ymax=90, nrows=18, ncols=36) values(r1) <- 1:ncell(r1) # use "c()" instead of "rast(sds())" r_p <- c(r1, r1, r1, r1, r1, r1) terra::time(r_p) <- ts_dt # not setting layer names as that concept does not exist in NetCDF s_p <- c(r1, r1, r1, r1, r1, r1) terra::time(s_p) <- ts_dt p <- sds(r_p, s_p) names(p) <- c("first", "second") p #class : SpatRasterDataset #subdatasets : 2 #dimensions : 18, 36 (nrow, ncol) #nlyr : 6, 6 #resolution : 10, 10 (x, y) #extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax) #coord. ref. : lon/lat WGS 84 #source(s) : memory #names : first, second 

You can see that time is preserved in the SpatRasterDataset

p[1] #or p["first"] #class : SpatRaster #dimensions : 18, 36, 6 (nrow, ncol, nlyr) #resolution : 10, 10 (x, y) #extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax) #coord. ref. : lon/lat WGS 84 #sources : memory # memory # memory # ... and 3 more source(s) #names : lyr.1, lyr.1, lyr.1, lyr.1, lyr.1, lyr.1 #min values : 1, 1, 1, 1, 1, 1 #max values : 648, 648, 648, 648, 648, 648 #time : 2021-01-01 to 2021-11-01 

The the NetCDF file also preserved the time stamps:

writeCDF(p, "test.nc", overwrite=TRUE) # Read back pp <- sds("test.nc") pp[1] #class : SpatRaster #dimensions : 18, 36, 6 (nrow, ncol, nlyr) #resolution : 10, 10 (x, y) #extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax) #coord. ref. : lon/lat WGS 84 #source : test.nc:first #varname : first #names : first_1, first_2, first_3, first_4, first_5, first_6 #time : 2021-01-01 to 2021-11-01 # Or read a single sub-dataset rast("test.nc", 1) #class : SpatRaster #dimensions : 18, 36, 6 (nrow, ncol, nlyr) #resolution : 10, 10 (x, y) #extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax) #coord. ref. : lon/lat WGS 84 #source : test.nc:first #varname : first #names : first_1, first_2, first_3, first_4, first_5, first_6 #time : 2021-01-01 to 2021-11-01 # Or read both sub-datasets as a single SpatRaster rast("test.nc") #class : SpatRaster #dimensions : 18, 36, 12 (nrow, ncol, nlyr) #resolution : 10, 10 (x, y) #extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax) #coord. ref. : lon/lat WGS 84 #sources : test.nc:first (6 layers) # test.nc:second (6 layers) #varnames : first # second #names : first_1, first_2, first_3, first_4, first_5, first_6, ... #time : 2021-01-01 to 2021-11-01 

If you get different results, please use the development version of terra instead: install.packages('terra', repos='https://rspatial.r-universe.dev')

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.