Próbuję wykonać połączenie na wspólnym polu po scaleniu dwóch sąsiednich plików kształtów. Pliki kształtów kończą się co najmniej jednym cienkim kawałkiem przestrzeni między nimi. Podczas próby połączenia pojawia się następujący błąd osieroconej dziury:
Błąd w createPolygonsComment (p): rgeos_PolyCreateComment: osierocona dziura, nie można znaleźć zawierającego wielokąt dla dziury o indeksie 17
Pod tym linkiem przesłałem do Dropbox odtwarzalny przykład .
Oto kod, aby odtworzyć problem:
#loading required packages
require(sp)
require(rgdal)
require(maptools)
require(rgeos)
#load example data, set "dsn=" to your working directory or specify the path
example <- readOGR(dsn=".",layer="ReproducibleExample")
#Attempting a UnionSpatialPolygons based on the COUNTY field
example.df <- as(example, "data.frame")
countycol <- example.df$COUNTY
example.diss <- unionSpatialPolygons(example, countycol)
Zwroty:
Błąd w createPolygonsComment (p): rgeos_PolyCreateComment: osierocona dziura, nie można znaleźć zawierającego wielokąt dla dziury o indeksie 17
Wypróbowanie poprawki zaproponowanej tutaj i tutaj :
slot(example, "polygons") <- lapply(slot(example, "polygons"), checkPolygonsHoles)
Zwraca ten sam błąd, który pochodzi z próby połączenia, ale z innym numerem indeksu:
rgeos_PolyCreateComment: osierocona dziura, nie można znaleźć zawierającego wielokąta dla dziury o indeksie 30
Wypróbowanie poprawki zaproponowanej w pomocnym samouczku Rogera Bivanda
fix <- slot(example, "polygons")
fixa <- lapply(fix, checkPolygonsHoles)
Zwraca ten sam błąd przy indeksie 30 jak powyżej.
Inni podnieśli ten problem tu i tutaj , i chociaż powyższe rozwiązania wydają się działać w niektórych przypadkach, inne nie zostały rozwiązane. Jeden użytkownik użył QGIS, aby rozwiązać problem, a drugi miał naprawione 2 z 3 elementów, ale brak ostatecznego rozwiązania.
Wygląda na to, że ludzie nadal mają problemy pomimo tego, że kod działa od czasu do czasu. Czy ktoś znalazł rozwiązanie w R?
W ArcGIS wykonałem narzędzie „napraw geometrię” i to rozwiązało problem, ale wygląda na to, że powinna być poprawka w R.
Odpowiedzi:
Przeanalizowałem problemy z geometrią w dołączonych danych i wydaje się, że nie tylko,
orphaned holes
ale takżegeometry validity issues
. To prawda, że anorphaned hole
jest w pewnym sensie kwestią poprawności geometrii, ale rgeos nie traktuje tego w taki sam sposób, jak w przypadku osieroconych otworów, zamiast zwykłego ostrzeżenia pojawia się błąd. Jak wskazujesz, są one wskazówkami, aby sprawdzić otwory wielokąta, ale nie zawsze udaje się, gdy są stosowane w celu naprawy osieroconych otworów.Więc:
wyczyść swoje dane (co jest wymagane, jeśli chcesz wykonywać geoprzetwarzanie jak związek)
użyj oczyszczonych danych w procesie związkowym
1. Czyszczenie geometrii Naprawianie geometrii w R może być czasem trudne, dlatego próbowałem zbudować eksperymentalny pakiet R (patrz https://github.com/eblondel/cleangeo ), który ma na celu ułatwienie czyszczenia
sp
obiektów (obecnie ograniczone do kształty wielokątne). Możesz zainstalować pakiet za pomocą:Na początek dobrze jest zobaczyć, jakie są problemy geometrii z danymi źródłowymi. W tym celu możesz wykonać następujące czynności (twoje dane są duże, więc może to zająć trochę czasu):
Dzięki temu zobaczysz, że twoje dane mają 2 rodzaje problemów:
orphaned holes
igeometry validity issues
. Obie (i nie tylko osierocone dziury) prawdopodobnie spowodują, żeunion
proces się nie powiedzie, więc dane powinny być wcześniej wyczyszczone, w miarę możliwości w sposób zautomatyzowany. W celu szybkiego odtworzenia pierwszy przykładowy kod poniżej bierze tylko podzbiór funkcji, które są oznaczone jako podejrzane (z wyjątkiem najnowszej, z indeksem = 9002 w oryginalnych danych - zobacz moją notatkę poniżej na ten temat)Jeśli
clgeo_Clean
dobrze spełnia swoje zadanie, wszystkie geometrie powinny być teraz prawidłowe. Możesz zastosować to do pełnego zestawu danych (z wyjątkiem indeksu funkcji = 9002)2. Proces Unii Teraz zobaczmy, czy
union
działa na tym zestawie danych:Uwaga: jak powiedziano wcześniej, usunąłem jedną funkcję (indeks = 9002) .Planując ją:
plot(sp[9002,])
zobaczysz, że ta funkcja jest bardzo (bardzo) złożona. Wyłączyłem go z próbki tylko dlatego, że sprawdzanie otworów zajmowało zbyt dużo czasu. Zobaczmy teraz, czy ten sam problem występuje przy użyciureadShapePoly
(zmaptools
) do odczytu danych ...3. Przejdź do readShapePoly vs. readOGR do odczytu danych (UPDATE)
readOGR
nie jest jedyną funkcją dostępną do odczytu plików kształtów. Możesz także użyćreadShapePoly
zmaptools
pakietu, ogólnie bardziej wydajnego niż pierwszy:Oprócz szybszego działania:
jeśli użyjesz powyższego kodu opartego na
clgeo_CollectionReport
, nie ma problemu z osieroconymi otworami, ale nadal problemy z geometrią.Czyszczenie geometrii za pomocą
clgeo_Clean
również działa dobrze, a teraz nie zacina się z indeksem funkcji 9002I ... proces związkowy działa.
Zobacz poniżej wynik wykresu:
Wniosek : wolą maptools, aby odczytać dane pliku kształtu , i rozważ użycie cleangeo do czyszczenia danych przed jakimkolwiek przetwarzaniem geoprzestrzennym .
źródło
Wygodnym rozwiązaniem, które nadal działa dla mnie w R, jest zastosowanie bufora o zerowej szerokości :
unionSpatialPolygons zajmuje trochę czasu z tym zestawem danych, ale wydaje się, że działa dobrze.
źródło