time units in years miscoded when saving with writeCDF

[dramanica]

When time is codes as (years), writing with writeCDF leads to a nc without units (and it looks like the time is coded as seconds, which might lead to issues if we have dates in the deep past). If we try to read the file we just generated back with rast, time is then lost. Here is a a simple reprex:

a <- rast(ncols=5, nrows=5, nl=50)
values(a) <- 1:prod(dim(a))
time(a, tstep="years") <- (1:50) * 20
writeCDF(a,"b.nc",varname = "bio01",overwrite=TRUE)
# read back the file
b <- rast("b.nc")
# but b has no time, as we just get NA with
time(b)

Here is an illustration of a simple fix (that could be implemented in writeCDF):

# to fix it, we remove the units before saving (so that we save the raw years, without conversion into other units)
time(a, tstep="") <- time(a)
writeCDF(a,"c.nc",varname = "bio01",overwrite=TRUE)
# and then fix the units in the netcdf file
nc_in <- ncdf4::nc_open("c.nc", write=TRUE)
ncdf4::ncatt_put(nc_in,varid="time", 
                 attname = "units",
                 attval = "years since 0-01-01 00:00:00.0")
#not strictly necessary, but it helps with some other libraries to recognise the time axis, so it would be a nice addition:
ncdf4::ncatt_put(nc_in, varid="time", attname="axis", attval = "T")
ncdf4::nc_close(nc_in)
# now when we read the file, we have time coded as years and with the same values as in the original raster
c <- rast("c.nc")
time(a)==time(c)

So, I would suggest writing the time in years, and using "years since 0-01-01 00:00:00.0" as the units (thus closely mirroring how time is represented in the original SpatRaster. Thanks!

[rhijmans]

Thanks for reporting this. The netcdf writing had not kept up with the expansion of the time steps supported by SpatRaster. I now get:

library(terra)
a <- rast(ncols=5, nrows=5, nl=50)
values(a) <- 1:prod(dim(a))
time(a, tstep="years") <- (1:50) * 20
writeCDF(a,"b.nc",varname = "bio01",overwrite=TRUE)
b <- rast("b.nc")

b
#class       : SpatRaster 
#dimensions  : 5, 5, 50  (nrow, ncol, nlyr)
#resolution  : 72, 36  (x, y)
#extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
#coord. ref. : lon/lat WGS 84 
#source      : b.nc 
#varname     : bio01 
#names       : bio01_1, bio01_2, bio01_3, bio01_4, bio01_5, bio01_6, ... 
#time (years): 20 to 1000 

time(b)
# [1]   20   40   60   80  100  120  140  160  180  200  220  240  260  280  300
#[16]  320  340  360  380  400  420  440  460  480  500  520  540  560  580  600
#[31]  620  640  660  680  700  720  740  760  780  800  820  840  860  880  900
#[46]  920  940  960  980 1000

[dramanica]

Thank you @rhijmans for the super fast fixes!