Home > Net >  Dealing with raster data in ggplot
Dealing with raster data in ggplot

Time:10-05

I have a vector spatial data for the boundary of a county and topographical data for the same county in raster format.

I need to create a map with ggplot so that it only shows those data in raster format that are within the county boundaries (which in turn are in a vector spatial file).

In other words, I need to remove raster data that is outside the county outline. Is it possible to do this with ggplot?

enter image description here

Reproducible example:

# load packages

library(elevatr) 
library(terra) 
library(geobr)

# get the municipality shapefile (vectorized spatial data)

municipality_shape <- read_municipality(code_muni = 3305802)

plot(municipality_shape$geom) 



# get the raster topographical data

prj_dd <- "EPSG:4674" 

t <- elevatr::get_elev_raster(locations = municipality_shape, 
                              z = 10,
                              prj = prj_dd) 

obj_raster <- rast(t) 

plot(obj_raster) 



# create the ggplot map

df_tere_topo <- obj_raster %>%
  as.data.frame(xy = TRUE) %>%
  rename(altitude = file40ac737835de)

ggplot() 
  geom_raster(data = df_tere_topo, aes(x = x, y = y, fill = `altitude`)) 
  geom_sf(municipality_shape, mapping = aes(), color = 'red', fill = NA)

CodePudding user response:

I think the most efficient way is to crop your SpatRaster to the extent of your vector data. With this approach the plotting is more efficient since you are not using data that you don't wnat to plot.

Another option is to set limits in coord_sf().

On this reprex I am using the package tidyterra as well, that has some functions to work with ggplot2 terra (I am the developer of tidyterra):

# load packages
library(elevatr)
library(terra)
library(geobr)
library(dplyr)
library(ggplot2)
# Use tidyterra
library(tidyterra)

# get the municipality shapefile (vectorized spatial data)

municipality_shape <- read_municipality(code_muni = 3305802)



# get the raster topographical data

prj_dd <- "EPSG:4674"

t <- elevatr::get_elev_raster(
  locations = municipality_shape,
  z = 10,
  prj = prj_dd
)

obj_raster <- rast(t)

# Option1: Crop first

obj_raster_crop <- crop(obj_raster, vect(municipality_shape))



# create the ggplot map
# using tidyterra
ggplot()  
  geom_spatraster(data = obj_raster_crop)  
  geom_sf(municipality_shape, mapping = aes(), color = "white", fill = NA)  
  coord_sf(expand = FALSE)  
  labs(fill="altitude")

enter image description here

 # Option 2: Use limits and no crop
lims <- sf::st_bbox(municipality_shape)

ggplot()  
  geom_spatraster(
    data = obj_raster,
    # Avoid resampling
    maxcell = ncell(obj_raster)
  )  
  geom_sf(municipality_shape, mapping = aes(), color = "white", fill = NA)  
  coord_sf(
    expand = FALSE,
    xlim = lims[c(1, 3)],
    ylim = lims[c(2, 4)]
  )  
  labs(fill="altitude")

enter image description here

  • Related