使用R进行转换为NetCDF时保留栅格变量名称

Question

使用R进行转换为NetCDF时保留栅格变量名称

5

使用以下代码将多年份的每月温度数据的栅格文件读入，文件名可通过以下格式访问names(object)：'Jan.1981'，'Feb.1981'等（适用于下面代码的两个年份示例文件 here - 添加所有文件会使其过大。

使用以下代码将其读入并写入到NetCDF中：

#Load Packages
library(raster)
library(ncdf4)

#Read in temperature files
r1 <- brick('TavgM_1981.grd')
r2 <- brick('TavgM_1982.grd')

#stack them together 
TempStack = stack(r1, r2)

#set the coordinate system (as it was missing)
crs(TempStack) <- ('+proj=lcc +lat_1=53.5 +lat_2=53.5 +lat_0=46.834 +lon_0=5 +x_0=1488375 +y_0=-203375 +datum=WGS84 +to_meter=2500 +no_defs +ellps=WGS84 +towgs84=0,0,0')

#reproject to get in lat/lon instead of meters
TempStack<-projectRaster(TempStack, crs=CRS("+init=epsg:4326"))

#Extract monthly data names to assign to netCDf later
names <- names(TempStack)

#write the raster file to NetCDF
writeRaster(TempStack, "Temp.nc", overwrite=TRUE, format="CDF",     varname="Temperature", varunit="degC", 
        longname="Temperature -- raster stack to netCDF, monthly average", xname="Longitude",   yname="Latitude", zname='Time', zunit=names)

当我将数据写入NetCDF并绘制月度数据时，它按照从第1个月到第24个月的顺序进行组织，但我希望它显示为“1981年1月”，“1981年2月”等。

我认为在writeRaster中添加zunit参数会起作用，但实际上并没有，数字仍然是1-24，而不是Jan，Feb等。

- Pad

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- dww · Accepted Answer

你的例子中存在一些误解。首先，你应该意识到netcdf维度中的值必须是数字。它们不仅仅是层的标签，而是该维度的实际值，因此不能采用像"Jan.1980"这样的字符串值。解决这个问题的方法是将netcdf文件保存后，将z维度的值作为数字值添加到其中。不幸的是，这意味着我们不能使用日期/时间变量类型，而必须先将它们转换为数字等价物。在这里，我使用lubridate包来完成这个任务。

# first we write the netcdf file to disk
writeRaster(TempStack, "Temp.nc", overwrite=TRUE, 
            format="CDF",     varname="Temperature", varunit="degC", 
            longname="Temperature -- raster stack to netCDF, monthly average", 
            xname="Longitude",   yname="Latitude", zname='Time', zunit='seconds')

# and open a connection to it to make changes.
# note that we use write=TRUE so that we can change it
nc = nc_open('Temp.nc', write = TRUE)

# now convert the strings to numeric values based on their dates
zvals = lubridate::parse_date_time(names, orders = 'm.y', tz = "UTC")
zvals = as.integer(zvals)

# and we can write these numeric dates to the z dimension
ncdf4::ncvar_put(nc, 'Time', zvals)

将日期写入z维度后，如果想要将数字z值转换回类似于“Jan.1908”等的栅格图层名称，我们还需要反向此过程。同样，lubridate可以帮助实现。

ncb = brick('Temp.nc')
zvals = ncvar_get(nc, 'Time')
zvals =  as.POSIXct(zvals, origin = lubridate::origin, tz = "UTC")
znames = paste0(lubridate::month(zvals, label=T), '.', lubridate::year(zvals))
names(ncb) = znames

让我们来检查一下是否成功：

plot(ncb)