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?
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")
# 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")