rspatial / terra

R package for spatial data handling https://rspatial.github.io/terra/reference/terra-package.html
GNU General Public License v3.0
540 stars 89 forks source link

`simplifyGeom` with topology preserve #587

Closed AMBarbosa closed 2 years ago

AMBarbosa commented 2 years ago

Hi, terra::simplifyGeom is better than rgeos::gSimplify in that it preserves the attribute table, but it's worse in that it doesn't have an option for preserving topology. Any chance this could be included? Cheers!

rhijmans commented 2 years ago

In fact it is the other way around. That is, the option is hard-coded, because I could not think of a situation where I would not want that. Am I wrong? See:

> getMethod("simplifyGeom", "SpatVector")
Method Definition:

function (x, ...) 
{
    .local <- function (x, tolerance = 0.1) 
    {
        preserveTopology <- TRUE
        x@ptr <- x@ptr$simplify(tolerance, preserveTopology)
        messages(x, "simplifyGeom")
    }
    .local(x, ...)
}
AMBarbosa commented 2 years ago

Excellent point ;) I actually meant that rgeos::gSimplify (with topologyPreserve = TRUE) apparently prevents adjacent polygons from overlapping. But I notice now that it also does some strange things to some of the polygons (e.g. bottom right in the image below):

library(terra)
library(raster)
v <- vect(system.file("ex/lux.shp", package="terra"))

v1 <- terra::simplifyGeom(v, tolerance = 0.2)
v2 <- rgeos::gSimplify(as(v, "Spatial"), tol = 0.2, topologyPreserve = TRUE)

par(mfrow = c(1, 2))
plot(v1, main = "simplifyGeom")
plot(v, border = "darkgrey", add = TRUE)
plot(v2, main = "gSimplify")
plot(v, border = "darkgrey", add = TRUE)

image

And a closer look reveals there are still some overlaps there as well:

par(mfrow = c(1, 1))
plot(v2, main = "gSimplify", xlim = c(6, 6.2), ylim = c(49.5, 49.8), col = 1:length(v2), alpha = 0.5)
plot(v, border = "darkgrey", add = TRUE)

image

So it seems we're actually better off with terra::simplifyGeom overall :)

rhijmans commented 2 years ago

Thanks. I have also added argument preserveTopology=TRUE so that you can set it to FALSE for when that is useful (e.g. remove small polygons)