Uzyskiwanie TopologyException: Dane wejściowe geom 1 są nieprawidłowe, co wynika z przecięcia się w R?


24

Powszechnie dyskutowano na temat błędu własnego przecięcia „TopologyException: Input geom 1 is invalid”, który wynika z nieprawidłowych geometrii wielokątów. Jednak nie znalazłem w Internecie wygodnego rozwiązania, które opierałoby się wyłącznie na funkcjonalności R.

Na przykład, udało mi się utworzyć obiekt 'SpatialPolygons spacerem od wyjścia map("state", ...)po miłą odpowiedź Josh O'Briena tutaj .

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)

stany

Problem z tym szeroko stosowanym zestawem danych polega na tym, że samo-przecięcie występuje w punkcie podanym poniżej.

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

Niestety ten problem uniemożliwia dalsze użycie „spydf_states”, np. Podczas dzwonienia rgeos::gIntersection. Jak mogę rozwiązać ten problem z poziomu R?


1
Jeśli zbliżysz się do tego punktu: plot(spydf_states, xlim=c(-122.1,-122.3),ylim=c(38,38.1))zobaczysz, że nie ma w tym „pozornie” - istnieje skrzyżowanie.
Spacedman

Odpowiedzi:


39

Użycie bufora o zerowej szerokości usuwa wiele problemów topologicznych w R.

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

Jednak praca z nieprojektowanymi współrzędnymi długości może powodować generowanie rgeosostrzeżeń.

Oto rozszerzony przykład, który najpierw przestawia się na projekcję Albersa:

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)

wprowadź opis zdjęcia tutaj


4
jakiś dodatkowy komentarz / czytanie na temat tego, dlaczego działa gBuffer„hack”?
MichaelChirico,

czy chcesz użyć gSimplify, który zrywa dane.frame i przekształca SPDF w przestrzenny obiekt wielokąta?
wnursal

5
Jeśli używasz sf, możesz także użyćsf::st_buffer(x, dist = 0)
Phil

działa również w niektórych przypadkach przy użyciuPostGIS
natsuapo
Korzystając z naszej strony potwierdzasz, że przeczytałeś(-aś) i rozumiesz nasze zasady używania plików cookie i zasady ochrony prywatności.
Licensed under cc by-sa 3.0 with attribution required.