如何在R中重新编码netcdf文件的非标准日期维度？

Question

我正在寻找关于如何重新编码netcdf文件/栅格砖的时间/日期维度的提示。我提取了月度温度数据，日期维度不是标准日期格式。例如，1850.04166666667,1850.125，分别表示1850年1月和2月的数据，R无法理解它们。这是数据结构

> eco.sst
class       : RasterBrick 
dimensions  : 180, 360, 64800, 2030  (nrow, ncol, ncell, nlayers)
resolution  : 1, 1  (x, y)
extent      : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
data source : C:\Users\Admin\AppData\Local\Temp\Rtmp6vFWNG\raster\r_tmp_2019-04-18_183528_11516_79080.grd 
names       : X1850.04166666667,     X1850.125, X1850.20833333333, X1850.29166666667,     X1850.375, X1850.45833333333, X1850.54166666667,     X1850.625, X1850.70833333333, X1850.79166666667,     X1850.875, X1850.95833333333, X1851.04166666667,     X1851.125, X1851.20833333333, ... 
min values  :      -0.668686926,  -0.273689777,      -1.902773499,      -3.406341553,  -3.626811743,      -2.131400585,      -1.569969416,  -1.630665064,      -2.418994188,      -1.971702099,  -0.822018623,      -3.423746109,      -1.303600550,  -0.786648691,      -1.452626705, ... 
max values  :      3.926926e+00,  1.550823e+00,      1.151636e+00,      6.622851e-01,  9.673859e-01,      7.178870e-01,      9.010030e-01,  4.363060e-01,      5.231520e-01,      3.346115e-01,  2.156055e-01,      6.418970e-01,      2.259051e+00,  3.802529e+00,      2.077996e+00, ... 
time        : 1850.04166666667, 2019.125 (min, max)

我写了一个函数来替换小数点;

fun.repl = function (x, na.rm = TRUE) str_replace(x, ".04166666667", "-1") 

如此1850.04166666667，即1850年1月成为1850-1;

然后我将其应用于以下形式的栅格：

stackApply(eco.sst, "time", fun.repl)

eco.sst是栅格砖，而“时间”是参考时间维度。但功能不起作用。我怀疑我不能很好地理解数据结构以成功应用该功能。我是新手，也发现这是netcdf的一个罕见问题。

我期望1850.04166666667,1850.125,1850.20833333333,1850.29166666667 ...这意味着1850年1月到4月成为标准日期...... 1850-1,1850-2,1850-3,1850-4等

Answer 1

这是另一种获取日期的方法

todate <- function(x) {
    year <- trunc(x)
    month <- round(12 * ((x - year) + 1/24))
    as.Date(paste0(year, "-", month, "-15")) 
}

dts <- c(1850.04166666667, 1850.125)
d <- todate(dts)
d
#[1] "1850-01-15" "1850-02-15"

请注意，我每个月使用第15天才能获得有效日期 - 但这似乎是预期的。

现在您可能想要使用setZ，请参阅下面的示例。 （stackApply在这里没有意义，你使用的论点也无效）

library(raster)
r <- raster(ncol=10, nrow=10)
s <- stack(lapply(1:2, function(x) setValues(r, runif(ncell(r)))))
s <- setZ(s, d)
s
#class      : RasterStack 
#dimensions : 10, 10, 100, 2  (nrow, ncol, ncell, nlayers)
#resolution : 36, 18  (x, y)
#extent     : -180, 180, -90, 90  (xmin, xmax, ymin, ymax)
#crs        : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
#names      :    layer.1,    layer.2 
#min values : 0.01477963, 0.01178438 
#max values :  0.9980334,  0.9939610 
#time        : 1850-01-15, 1850-02-15