0

I am trying to do something fairly simple: I want to clip a vector with point features down to only features that are contained within a raster extent. For raster I am using the package terra, for vector the package sf. I am surprised how interoperable these two common packages are in this regard, seeing that extract(terra::SpatRast, sf::sf_vector_df) works just fine.

My approaches so far have led to only errors. What is a good way to make a terra::raster extent and a sf-vector-dataset interoperable?

vhm = rast("xy.tif")
centers = st_read("centers.shp")

# trying straightforward
centers_vhm = st_intersection(centers, vhm)
>> Error in UseMethod("st_geometry") : no applicable method for 'st_geometry' applied to an object of class "SpatRaster"

# trying with extent
centers_vhm = st_intersection(centers, ext(vhm))
>> Error in geos_op2_geom("intersection", x, y, ...) : 
  st_crs(x) == st_crs(y) ist nicht TRUE

# trying to convert extent to sf object to assign CRS...
vhm_extent = st_as_sf(ext(vhm))
>> Error in UseMethod("st_as_sf") : 
  no applicable method for 'st_as_sf' applied to an object of class "SpatExtent"

1 Answer 1

1

I suggest you start with the bounding box of your raster via sf::st_bbox() and cast it to sfc object via sf::st_as_sfc(). sfc objects can be used to clip a sf vector shapefile.

library(sf)
library(terra)

# a raster, any raster...
raster <- rast(system.file("ex/elev.tif", package="terra"))
# a vector, any vector...
vector <- giscoR::gisco_get_countries()

# this is the action!
bbox <- st_bbox(raster) %>% 
  st_as_sfc()

isection <- st_intersection(vector, bbox)

# validate results
mapview::mapview(isection, zcol = "NAME_ENGL")

image of world countries clipped to the extent of elevation of Luxembourg raster

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.