2

I'd like to see the result after an hour of working while substracting two shapefiles. This is the code, but don't ask me what in those files. Whatever it is, I'd like to know how to skip those errors and see the result even if it's fake. Those files are huge.

shp1 <- rgdal::readOGR(dsn = paste0("LPDF/SHP/",NAME), layer = layer1, encoding = "utf-8", stringsAsFactors = F)
shp2 <- rgdal::readOGR(dsn = paste0("LPDF/SHP/",NAME), layer = layer2, encoding = "utf-8", stringsAsFactors = F)

Here, I'm trying to use some actions to perform validation as found on some different QAs:

p_ss = sf::st_as_sf(shp1)
p_ss <- sf::st_buffer(p_ss, 0)
p_pol <- sf::st_cast(p_ss, "MULTIPOLYGON")
sb1 <- as(p_pol, "Spatial")

p_ss = sf::st_as_sf(shp2) #shp7
p_ss <- sf::st_buffer(p_ss, 0)
p_pol <- sf::st_cast(p_ss, "MULTIPOLYGON")
sb2 <- as(p_pol, "Spatial")

sb1 <- maptools::unionSpatialPolygons(sb1, IDs = 1:length(shp1), threshold = 1, avoidUnaryUnion = T)
sb2 <- maptools::unionSpatialPolygons(sb2, IDs = 1:length(shp2), threshold = 1, avoidUnaryUnion = T)

So up to this point it worked; without previous buffer the 'maptools' function fails. But after an hour this fails:

p_dif <- rgeos::gDifference(spgeom1 = sb1, spgeom2 = sb2, byid = F, unaryUnion_if_byid_false = F, drop_lower_td = T, checkValidity = T)

Error in RGEOSBinTopoFunc(spgeom1, spgeom2, byid, id, drop_lower_td, unaryUnion_if_byid_false, :

TopologyException: side location conflict at 505082.5 427876.5

Not sure what else could I use in R? byid=T crashes due to memory lack. Maybe I should try to use simple loop and just do it one by one thru the whole data (and day possibly)?

Can you give me some suggestions?


Using suggestions from comments I've tried to make use of automated topology corrections:

Just installed lwgeom using R only and together with other libraries:

p_s1 = sf::st_as_sf(shp1)
sv1 <- lwgeom::st_make_valid(p_s1)
p_pol <- sf::st_cast(sv1, "MULTIPOLYGON")
sb1 <- as(p_pol, "Spatial")

p_s2 = sf::st_as_sf(shp2)
sv2 <- lwgeom::st_make_valid(p_s2)
p_pol <- sf::st_cast(sv2, "MULTIPOLYGON")
sb2 <- as(p_pol, "Spatial")

made these shapefiles valid. But in the substraction:

p_dif <- rgeos::gDifference(spgeom1 = sb1, spgeom2 = sb2, byid = F, unaryUnion_if_byid_false = T, drop_lower_td = T, checkValidity = T)

it only worked if with unary union (TRUE) and byid=FALSE, not sure about dropping. The unary union make it very fast as well. This is enough for me.

2

One option might be to try and use the st_make_valid function from sf. For instance, you might identify the features that have a topology issue as follows:

st_is_valid(p, reason = TRUE)

and then try to resolve it as follows:

st_make_valid(p)

These steps are described in more detail here along with an additional option of buffering the features by a zero distance.

3
  • Thank you for those suggestions but none of that removed my Self-intersection errors. – Peter.k Dec 17 '19 at 18:05
  • it may require manual editing to fully resolve. alternatively, you could try GRASS GIS (and if it has to be within R, you could potentially use spgrass7, but I haven't tested the latter for repairing grass geometries from within in R). – dmci Dec 17 '19 at 21:11
  • Thanks again. I've managed to have full use from this post and now it's finally my successful solution. – Peter.k Dec 18 '19 at 10:32

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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