3

I am writing a rasterbrick with multiple layers to a netcdf file using the raster package in R:

writeRaster(gridfile, "Data", "CDF", overwrite=TRUE,
     varname="Data",varunit="mm.d-1",longname="Monthly data",
     xname="lon",yname="lat", zname="time", zunit="month",bylayer=FALSE,NAflag=-9999)

Is there a way to have the data with time coordinates in the written ncdf file? E.g. in the NetCDF Viewer Panoply I can plot the data in [lon][lat] but not in [time][lat].

Also how can specify the timesteps as specific months in a year instead of the timesteps starting from 1 and increasing?

  • I have found a solution using the ncdf package: stat.ethz.ch/pipermail/r-sig-geo/2013-January/017265.html However, is there a way to directly do this through the raster package? – Dom Apr 19 '13 at 12:36
  • Can you add some data, or at least include str(gridfile) and sessionInfo() so we know your R, raster, and ncdf/ncdf4 versions? – mdsumner Apr 22 '13 at 21:20
  • I wrote that answer and unfortunately I am still looking for a way to add time slices into netcdf files using the raster package. – thiagoveloso Oct 22 '15 at 17:57
3

You should use the R-Package netcdf, it has much more possibilities to taylor your resulting data fields. To hold the GIS compatibility, use the Lon/Lat attribute schema beside the time attribute (Spatial temporal advises ).

If you want to create and use spatio temporal datasets like in climatology you should read the best practice article of unidata.ucar.edu. There are several tools which work with netcdf spatio temporal standard structures.

Mapping netcdf - NCL | NCO netCDF-Operator | CDO - Climate Data Operators


NCL Example

enter image description here

2

I know this is an old question, but if anyone is still looking for an answer on how to do this in R the following code that will make a netCDF with 2 variables (prec, tmax), and 3 dimensions (long, lat, time). Pretty straight forward to change this suit your number of variables or time dimension etc.

library(raster)
library(ncdf4)
prec <- getData("worldclim", res = 10, var = "prec")
tmax <- getData("worldclim", res = 10, var = "tmax")
tmax[] <- tmax[]/10 # tmax is scaled by 10
# Change NA values to -999
prec[is.na(prec)] <- tmax[is.na(tmax)] <- -999

# make sure the rasters are identical
compareRaster(prec, tmax, res = TRUE, orig = TRUE) ## TRUE

# output filename
filename <- "testMultiDim.nc"

# Longitude and Latitude data
xvals <- unique(values(init(prec, "x")))
yvals <- unique(values(init(prec, "y")))
nx <- length(xvals)
ny <- length(yvals)
lon <- ncdim_def("longitude", "degrees_east", xvals)
lat <- ncdim_def("latitude", "degrees_north", yvals)

# Missing value to use
mv <- -999

# Time component
time <- ncdim_def(name = "Time", 
                  units = "months", 
                  vals = 1:12, 
                  unlim = TRUE,
                  longname = "Month_of_year")

# Define the precipitation variables
var_prec <- ncvar_def(name = "precipitation",
                      units = "mm/month",
                      dim = list(lon, lat, time),
                      longname = "Monthly_Total_Precipitation",
                      missval = mv,
                      compression = 9)

# Define the temperature variables
var_temp <- ncvar_def(name = "temperature",
                      units = "degC",
                      dim = list(lon, lat, time),
                      longname = "Monthly_Avg_Temperature_degC",
                      missval = mv,
                      compression = 9)

# Add the variables to the file
ncout <- nc_create(filename, list(var_prec, var_temp), force_v4 = TRUE)
print(paste("The file has", ncout$nvars,"variables"))
print(paste("The file has", ncout$ndim,"dimensions"))

# add some global attributes
ncatt_put(ncout, 0, "Title", "MultiDimesionsalNCDFTest")
ncatt_put(ncout, 0, "Source", "Some example data from the raster package")
ncatt_put(ncout, 0, "References", "See the raster package")
ncatt_put(ncout, 0, "Created on", date())

# Place the precip and tmax values in the file
# need to loop through the layers to get them 
# to match to correct time index
for (i in 1:nlayers(prec)) { 
  #message("Processing layer ", i, " of ", nlayers(prec))
  ncvar_put(nc = ncout, 
            varid = var_prec, 
            vals = values(prec[[i]]), 
            start = c(1, 1, i), 
            count = c(-1, -1, 1))
  ncvar_put(ncout, var_temp, values(tmax[[i]]), 
            start = c(1, 1, i), 
            count = c(-1, -1, 1))
}
# Close the netcdf file when finished adding variables
nc_close(ncout)

You can open it and check the dimensions are correct

# Open the netcdf file
nc <- nc_open(filename)
nc
File testMultiDim.nc (NC_FORMAT_NETCDF4):

     2 variables (excluding dimension variables):
        float precipitation[longitude,latitude,Time]   (Chunking: [1080,450,1])  (Compression: level 9)
            units: mm/month
            _FillValue: -999
            long_name: Monthly_Total_Precipitation
        float temperature[longitude,latitude,Time]   (Chunking: [1080,450,1])  (Compression: level 9)
            units: degC
            _FillValue: -999
            long_name: Monthly_Avg_Temperature_degC

     3 dimensions:
        longitude  Size:2160
            units: degrees_east
            long_name: longitude
        latitude  Size:900
            units: degrees_north
            long_name: latitude
        Time  Size:12   *** is unlimited ***
            units: months
            long_name: Month_of_year

    4 global attributes:
        Title: MultiDimesionsalNCDFTest
        Source: Some example data from the raster package
        References: See the raster package
        Created on: Wed Nov 28 10:35:53 2018
1

Valid units strings for CF-compliant netcdf files look like "seconds since 1970-01-01 00:00:00" or "days since 2013-04-19 00:00:00"

I think "months since ...." isn't supported by udunits since months are not equally sized.

You might do better thinking of the z dimension as a record dimension and then add a time variable with ncdf that locates each record in time.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.