Saturday, 17 March 2018

raster - How to correct link a set of colors to a rasterLayer


I have downloaded a raster with land cover data from http://www.diva-gis.org/gdata.



library(tmap)
library(raster)

grid_cov <- raster('VEN_cov/VEN_cov.grd')
grid_cov

class : RasterLayer
dimensions : 1416, 1656, 2344896 (nrow, ncol, ncell)
resolution : 0.008333333, 0.008333333 (x, y)
extent : -73.5, -59.7, 0.5, 12.3 (xmin, xmax, ymin, ymax)

coord. ref. : +proj=longlat +ellps=WGS84
data source : /Users/altons/Dropbox/datacamp/geospatial course/VEN_cov/VEN_cov.grd
names : VEN_cov
values : 1, 22 (min, max)

Values in raster are:


unique(values(grid_cov))
20 15 14 2 13 16 12 1 19 22 7 17 18 8

which according to this document(pages 19 & 20)represent the type of land cover.



Now I have tried to map specific colors to land cover codes but my attempts have failed so far:


Attempt 1: pass color in ascending order as per land cover codes


cover_palette <- c("#008B00", "#00CD00", "#7CCD7C", "#98FB98", "#FFEFDB", "#FAEBD7", "#A4D3EE", "#CDC673", "#FF7256", "#FF7F50", "#EEE8CD", "#1874CD", "#F0F8FF", "#CD2626")

tm_shape(grid_alt) + tm_raster(legend.show = F,palette=cover_palette,style='cat' )

Bits in red/blue are mountains... which is wrong based on my definition of colors.


enter image description here


Attempt 2: Pass colors in the same order as per values(grid_cov)


cover_palette <- c('#F0F8FF',

'#CDC673',
'#A4D3EE',
'#00CD00',
'#FAEBD7',
'#FF7256',
'#FFEFDB',
'#008B00',
'#1874CD',
'#CD2626',
'#7CCD7C',

'#FF7F50',
'#EEE8CD',
'#98FB98')

tm_shape(grid_alt) + tm_raster(legend.show = F,palette=cover_palette,style='order' )
save_tmap(filename="~/Desktop/attempt2.png",dpi=75)

enter image description here


The last plot - the red bits are mountains which I thought I have colored in white....


Anyway - I think my question is: If I created a data.frame with all the colors then how can I merge it with the RasterLayer object? and therefore use colors$color instead of guessing!



my data frame looks like this:


lc_id = c(20,15,14,2,13,16,12,1,19,22,7,17,18,8)
cover_palette <- c('#F0F8FF','#CDC673',
'#A4D3EE','#00CD00',
'#FAEBD7','#FF7256',
'#FFEFDB','#008B00',
'#1874CD','#CD2626',
'#7CCD7C','#FF7F50',
'#EEE8CD','#98FB98')


colors <- data.frame(id = lc_id,color=cover_palette)

Answer



Install the rasterVis package, for we are going to use its levelplot functionality...


Prep:


library(raster)
library(rasterVis)

Read the raster in:


gc = raster("./VEN_cov.grd")


Convert to categorical:


gc = ratify(gc)
levels(gc)

Now I don't like numeric categorical levels, because its easy to mess them up with indices, so I'll always make character ones. If you have short names for the categories (like "Urban", "Forest", "Woodland") use them here. Long names will mess up the key. Anyway:


levels(gc)[[1]]$code = paste0("GLC-",levels(gc)[[1]]$ID)
levels(gc)

Your colors data frame is good, but use character colours instead of factors so they don't get converted to numbers inadvertently, and sort it so its in the same order as levels(gc):


colors <- data.frame(id = lc_id,color=cover_palette,stringsAsFactors=FALSE)

csort = colors[order(colors$id),]

Now the magic:


levelplot(gc,col.regions=csort$color)

enter image description here


which I think is right. There's a few bits of red urban, and a lot of forest classes.


To do in tmap, use:


> tm_shape(gc)+tm_raster(palette=csort$color)

No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...