Problem solving overlapping vectors in terra using terra::union

user_3808018 0
Viewed 58

I have several polygons of habitat types, and depending on the ones I select I need to get the percentage of a country. The issue is that many polygons intersect, so I need to solve the intersections so that I don’t add the area twice and assume that the percentage of the country is higher than I should calculate.

library(sf)
library(terra)
hab1 <- sf::read_sf("/vsicurl/https://github.com/Sustainscapes/BiodiversityCounsil/raw/master/ForQuestion/Overdrev.shp") |>
  terra::vect()

hab2 <- sf::read_sf("/vsicurl/https://github.com/Sustainscapes/BiodiversityCounsil/raw/master/ForQuestion/Strandeng.shp") |>
  terra::vect()

Test1 <- terra::union(hab1, hab2)

This is not really working as it should, since the number of features in test1 is 21528, however the length of each variable is only 3

as seen in the following code, and that wont allow for subsetting:

length(Test1)
# [1] 21528
length(Test1$Overdrev)
# [1] 3
length(Test1$Strandeng)
# [1] 3

(a more elaborate version of this question is available here)

r-raster shapefile gis r
1 Answers

I think this was fixed a while ago, but not in the CRAN version yet. Can you try this with the development version please? You can install it with install.packages('terra', repos='https://rspatial.r-universe.dev')

With that version, I get:

library(terra)
hab1 <- vect("/vsicurl/https://github.com/Sustainscapes/BiodiversityCounsil/raw/master/ForQuestion/Overdrev.shp") 
hab2 <- vect("/vsicurl/https://github.com/Sustainscapes/BiodiversityCounsil/raw/master/ForQuestion/Strandeng.shp")
Test1 <- terra::union(hab1, hab2)
length(Test1)
#[1] 3

#data.frame(Test1)
#  Overdrev Strandeng
#1 Overdrev Strandeng
#2 Overdrev      <NA>
#3     <NA> Strandeng

Test1
# class       : SpatVector 
# geometry    : polygons 
# dimensions  : 3, 2  (geometries, attributes)
# extent      : 445333, 892786.7, 6050590, 6400331  (xmin, xmax, ymin, ymax)
# coord. ref. : ETRS89 / UTM zone 32N (EPSG:25832) 
# names       : Overdrev Strandeng
# type        :    <chr>     <chr>
# values      : Overdrev Strandeng
#               Overdrev        NA
#                     NA Strandeng

And the overlapping areas:

ovp <- Test1[!is.na(Test1$Overdrev) & !is.na(Test1$Strandeng), ]

More generally, it seems that you are using raster data that was transformed to polygons. That is asking for punishment. Can't you use the original raster data?

user_6352450
Related