Raster maps with mapsf

mapsf focuses on creating thematic maps based on vector data, maps are mainly generated from sf objects.
However, the latest update of mapsf (v0.9.0) expands its mapping capabilities for raster data. Emulating features from terra (and leveraging its plotting functions internally), mf_raster() now provides three map types for single band rasters: classes, continuous, and interval.
Raster maps benefit from all map layout elements (title, scale bar, north arrow, graticules…) and can be mixed with vector objects.

Download datasets used in examples:

Classes

This type of map is suitable for qualitative data (categories, typologies).
In this example we display a map of Corine Land Cover classes over the Lot region in France.

library(mapsf)

# import raster
clc <- terra::rast("clc_46.tif")

# import vector (Lot region polygon)
dep <- sf::st_read("dep.gpkg", "dep", quiet = TRUE) 

# define a color palette
clc_pal <- c("#E6004D", "#FF0000", "#CC4DF2", "#CC0000", "#E6CCE6", "#A600CC", 
             "#FFE6FF", "#FFFFA8", "#E68000", "#F2A64D", "#E6E64D", "#FFE64D", 
             "#E6CC4D", "#80FF00", "#00A600", "#4DFF00", "#CCF24D", "#A6FF80", 
             "#A6E64D", "#A6F200", "#00CCF2", "#80F2E6")

# map the raster
mf_raster(clc, pal = clc_pal, expandBB = c(0, 0.5, 0.05, 0),
          leg_title = "", leg_size = .7, leg_box_cex = c(.6, .85),
          leg_pos = "topleft")

# add region borders 
mf_map(dep, col = NA, lwd = 2, add = TRUE)

# add graticules
mf_graticule(clc, lty = 2, lwd = .3, cex = 0.5, col = "grey20")

# add map layout elements
mf_title("CORINE Land Cover - Lot region - 2018", tab = FALSE, pos = "center") 
mf_credits("T. Giraud, 2024 | EEA 2018")
mf_arrow(pos = "topright", adjust = clc)
mf_scale(10)

Continuous

This type of map is appropriate for quantitative data such as DEM, surface temperature or vegetation index.
In the following examples we display elevation maps of the Lot region.

# import raster
elev <- terra::rast("elev_46.tif")

# map the raster
mf_raster(elev, leg_pos = "topleft")

# add layout elements
mf_title("Elevation - Lot region") 
mf_credits("SRTM 90m DEM Digital Elevation Database, 2018")
mf_arrow(pos = "topright", adjust = elev)
mf_scale(5)

The colors are set using a palette name or a vector of colors to map the continuous values.
See ?hcl.colors() examples to check palette names and colors.

# select a palette with a palette name
mf_raster(elev, pal = "Rocket", expandBB = c(0,0,0,.3))
mf_title('Continuous palette based on "Rocket"')
# define a palette with colors (from red to black)
mf_raster(elev, pal = c("#fe5234", "#000000"), expandBB = c(0,0,0,.3))
mf_title('Continuous palette from "#fe5234" to "#000000"')

It is also possible to give break values and corresponding colors. Breaks are inflection points of the color scale. In this example the color palette varies from white to khaki between 67 and 150, then from khaki to darkgreen between 150 and 550, and finally from darkgreen to black between 550 and 776.

# set a cartographic theme
mf_theme("green", bg = "black")

# map the raster
mf_raster(elev, breaks = c(67, 150, 550, 776), 
          pal = c("white", "khaki", "darkgreen", "black"),
          leg_val_rnd = 0, leg_pos = "topleft", leg_box_cex = c(.65, 1),
          leg_title = "Elevation (m)", leg_horiz = TRUE, leg_frame = TRUE, 
          leg_frame_border = NA, leg_bg = "grey20")

# add an inset 
mf_inset_on(x = "worldmap", pos = "bottomright", cex = .2)
mf_worldmap(x = elev)
mf_inset_off()

# add layout elements
mf_title("Elevation - Lot region", tab = FALSE, pos = "center") 
mf_credits("T. Giraud, 2024 | SRTM 90m DEM Digital Elevation Database, 2018", 
           bg = "#000000CC")
mf_arrow(pos = "topright", adjust = elev)
mf_scale(5)

# add annotations
max_elev_coord <- terra::xyFromCell(elev, terra::where.max(elev)[2])[1,]
mf_annotation(x = max_elev_coord, "Highest Point", 
              col_arrow = "white", col_txt = "white", 
              pos = "bottomleft", halo = TRUE, cex = 1.2)

Interval

This type of map is also appropriate for quantitative data.
Numerical values are classified using one the methods provided by mf_get_breaks().

# set a cartographic theme
mf_theme("default")

# map the raster
mf_raster(elev, type = "interval", 
          breaks = "q6", pal = "Turku", rev = TRUE,
          leg_val_rnd = 1, leg_pos = "topleft", leg_title = "Elevation (m)")

# add layout elements
mf_title("Elevation - Lot region", tab = FALSE, pos = "center") 
mf_credits("SRTM 90m DEM Digital Elevation Database, 2018", 
           pos = "rightbottom")
mf_arrow(pos = "topright")
mf_scale(5, pos = "bottomleft")

Multiband

The last kind of map that mf_raster() can take care of is RGB composite maps using three bands.

In the following example we download an image from Open Street Map and display it.

library(maptiles)
dep <- sf::st_transform(dep, "EPSG:3857")

# get an OpenStreetMap image on the region
(osm <- get_tiles(x = dep, provider = "OpenStreetMap", crop = TRUE, zoom = 9))

class       : SpatRaster 
dimensions  : 433, 448, 3  (nrow, ncol, nlyr)
resolution  : 305.7481, 305.7481  (x, y)
extent      : 109152.1, 246127.2, 5496740, 5629129  (xmin, xmax, ymin, ymax)
coord. ref. : WGS 84 / Pseudo-Mercator (EPSG:3857) 
source(s)   : memory
colors RGB  : 1, 2, 3 
names       : lyr.1, lyr.2, lyr.3 
min values  :     3,    54,     3 
max values  :   251,   249,   244 

# set them margins
mf_theme(mar = c(0,0,0,0))

# display the map
mf_raster(osm)

# add region limits
mf_map(dep, add = TRUE, col = NA, lwd = 3)

# add attribution
mf_credits(get_credit("OpenStreetMap"))