0

When I save and reopen a raster in terra, I am losing the assigned units. Perhaps there are some filetypes that can save them while others can't, but I haven't found documentation of that.

library(terra)

A <- rast( nrows=10, ncols=10, xmin=0, xmax=10 )
values(A) <- 1:100
units(A) <- "Percent"
writeRaster( x=A, filename="test.tif" )

B <- rast( "test.tif" )
B
units(B)  # empty
Jim Worrall
  • 466
  • 3
  • 11

1 Answers1

2

As far as I know only the netCDF standard supports writing units to files (among raster file formats). The GeoTIFF standard and the data model of the GDAL library used for reading and writing raster data do not have the concept of units as far as I can see.

I have added a work-around for this problem in terra version 1.5.0 (currently the development version available from github)

Both units and timestamps are now saved to file, and with your example I get:

library(terra)
#terra version 1.5.0

A <- rast( nrows=10, ncols=10, xmin=0, xmax=10 )
values(A) <- 1:100
units(A) <- "%"
writeRaster(A, "test.tif", overwrite=TRUE)

B <- rast( "test.tif" )
B
#class       : SpatRaster 
#dimensions  : 10, 10, 1  (nrow, ncol, nlyr)
#resolution  : 1, 18  (x, y)
#extent      : 0, 10, -90, 90  (xmin, xmax, ymin, ymax)
#coord. ref. : lon/lat WGS 84 (EPSG:4326) 
#source      : test.tif 
#name        : lyr.1 
#min value   :     1 
#max value   :   100 
#unit        :     % 

units(B) 
#[1] "%"
Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63
  • Thanks, that's good to know. By the way, I noticed that by just assigning a value to a cell, unit also disappears: – Jim Worrall Nov 25 '21 at 21:15
  • [from previous comment] A[5,5]<-3 erases the unit. Also, assigning an attribute table causes the name of the attribute column to become the raster name. Not sure if that is intentional. With raster values 0:4: RAT <- data.frame( ID=c( 0,1,2,3,4 ), Class=c( "A", "B", "C", "D", "E" )); levels(A) <- RAT Raster name becomes "Class" – Jim Worrall Nov 25 '21 at 21:22
  • 2
    Thanks, the assigning is fixed now. The change in name is intended because you can have multiple attributes and this makes it clear which one is "active". Feel free to file more "issues" at https://github.com/rspatial/terra/issues – Robert Hijmans Nov 26 '21 at 23:16