Obtendo TopologyException: o geom de entrada 1 é inválido devido à interseção automática em R?

24

O erro de auto-interseção 'TopologyException: Input geom 1 is invalid' que surge de geometrias de polígonos inválidas foi amplamente discutido. No entanto, não encontrei uma solução conveniente na Web que dependa apenas da funcionalidade R.

Por exemplo, eu consegui criar um objeto 'SpatialPolygons' a partir da saída da map("state", ...)boa resposta de Josh O'Brien aqui .

library(maps)
library(maptools)

map_states = map("state", fill = TRUE, plot = FALSE)

IDs = sapply(strsplit(map_states$names, ":"), "[[", 1)
spydf_states = map2SpatialPolygons(map_states, IDs = IDs, proj4string = CRS("+init=epsg:4326"))

plot(spydf_states)

estados

O problema com esse conjunto de dados amplamente aplicado é agora que a auto-interseção ocorre no ponto fornecido abaixo.

rgeos::gIsValid(spydf_states)
[1] FALSE
Warning message:
In RGEOSUnaryPredFunc(spgeom, byid, "rgeos_isvalid") :
  Self-intersection at or near point -122.22023214285259 38.060546477866055

Infelizmente, esse problema impede qualquer uso adicional de 'spydf_states', por exemplo, ao chamar rgeos::gIntersection. Como posso resolver esse problema no R?

fdetsch
fonte
11
Se você aumentar o zoom em torno desse ponto: plot(spydf_states, xlim=c(-122.1,-122.3),ylim=c(38,38.1))verá que não há "aparentemente" a respeito - há uma auto-interseção.
Spacedman

Respostas:

39

Usar um buffer de largura zero limpa muitos problemas de topologia em R.

spydf_states <- gBuffer(spydf_states, byid=TRUE, width=0)

No entanto, trabalhar com coordenadas de longo prazo não projetadas pode causar rgeosavisos.

Aqui está um exemplo estendido que reprojeta primeiro uma projeção de Albers:

library(sp)
library(rgeos)

load("~/Dropbox/spydf_states.RData")

# many geos functions require projections and you're probably going to end
# up plotting this eventually so we convert it to albers before cleaning up
# the polygons since you should use that if you are plotting the US
spydf_states <- spTransform(spydf_states, 
                            CRS("+proj=aea +lat_1=29.5 +lat_2=45.5 +lat_0=37.5 +lon_0=-96"))

# simplify the polgons a tad (tweak 0.00001 to your liking)
spydf_states <- gSimplify(spydf_states, tol = 0.00001)

# this is a well known R / GEOS hack (usually combined with the above) to 
# deal with "bad" polygons
spydf_states <- gBuffer(spydf_states, byid=TRUE, width=0)

# any bad polys?
sum(gIsValid(spydf_states, byid=TRUE)==FALSE)

## [1] 0

plot(spydf_states)

insira a descrição da imagem aqui

hrbrmstr
fonte
4
algum comentário / leitura extra sobre por que o gBuffer"hack" funciona?
MichaelChirico
você deseja usar o gSimplify, pois destrói o data.frame e converte o SPDF em objeto de polígono espacial?
wnursal
5
Se você estiver usando, sfvocê também pode usarsf::st_buffer(x, dist = 0)
Phil
também funciona em alguns casos ao usarPostGIS
natsuapo