4

I am generating landscape classification output rasters and adapted the code from WAF’s stackoverflow answer to the this question to generate GeoTIFFs with embedded standardized color and attribute information.

library(rgdal)
library(raster)

addColorTable <- function(inRstName, outRstName, rat.df){
  r<- readGDAL(inRstName)
  rat.df$color<- as.character(rat.df$color)
  rat.df$attribute<- as.character(rat.df$attribute)
  outRst <- writeGDAL(r, outRstName, type="Byte", 
                      colorTable=list(rat.df$color), 
                      catNames=list(rat.df$attribute), mvFlag=11L)
  return(raster(outRst))
}

# create dummy data set
r <- raster(nrow=10, ncol=10)
r[] <- 0
r[51:100] <- 1
r[3:6, 1:5] <- 2
r[1, 1] <- 3
writeRaster(r,'dummy_raster.tif',overwrite=T)

#create color and attribute table
# This defines the values, the color and the attribute
valT <- c(0,1,2,3)
colT <-  c("#FF0000", "#FF9900" ,"#99FF00","#0000FF")
attT <- c('Forest','Water body','City','Cropland')
rat.df <- data.frame(value=valT,color=colT,attribute=attT)

# apply the magic function
rnew <- addColorTable('dummy_raster.tif', 'dummy_raster_with_symbology.tif', rat.df)

plot(rnew)

enter image description here

The code works perfectly within R to generate the categorical land cover rasters with the proper colors (see figure above). However, no matter what I do, when I visualize the categorical raster outside of R (e.g., QGIS or ArcGIS) the raster legend includes legend elements for 0-255 raster values, even if my raster only has values between 0 and 3 (see figure below). Most importantly, element 11 in the legend is always displayed with a 0% opacity symbology despite nothing in my code specifying opacity.

Does anyone have any solutions to correctly save raster color and attribute info in R?

enter image description here In the case above, the transparency of class 11 is not an issue because there are only 4 categorical classes, but in my actual classification rasters, my 11th class gets 'erased' as well.

5
  • 2
    Unfortunately, this is the way it has to be as the raster is encoded to 8 bit so, has 255 values. If you look at common landcover data, distributed with embedded color, (eg., NLCD *.img) you will see the same thing. Apr 29 at 17:29
  • Thanks Jeff. But what explains the weird 0% opacity of class 11? Every categorical raster with attribute table I export from r has this 'bug' in it... Apr 29 at 17:33
  • Not sure what is going on with value 11. My first guess is that it is being assigned the NoData value but, this seems like an odd NA value. It is also difficult to identify this as an R, raster or GDAL issue. Have to tried to assign an explicit NA value to the 8-bit range (eg., 255=NA in the rat.df)? Apr 29 at 18:27
  • Hmm... just realized that the opacity issue is only present when visualizing the data in qgis not arcgis. Apr 29 at 19:41
  • 1
    in QGIS properties > symbology if you press classify only the four classes remain
    – Elio Diaz
    Apr 29 at 23:59
4

This should work better with terra.

library(terra)
#terra version 1.2.4
r <- rast(nrow=10, ncol=10)
r[] <- 0
r[51:100] <- 1
r[3:6, 1:5] <- 2
r[1, 1] <- 3

coltab(r) <- c("#FF0000", "#FF9900" ,"#99FF00","#0000FF")
levels(r) <- c('Forest','Water body','City','Cropland')
x <- writeRaster(r, "test.tif", overwrite=TRUE)
plot(x)

enter image description here

1
  • Thanks! I have been slowly transitioning to terra and it is incredibly fast and powerful. I am trying to get my colleagues to switch as well. May 5 at 21:56

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.