Identyfikacja „długich i wąskich” wielokątów w PostGIS

Mam zestaw wielokątów reprezentujących duże obszary, powiedzmy dzielnice miast. Chcę zidentyfikować duże nakładające się obszary między nimi.

Ale jest problem: czasami wielokąty zachodzą na siebie na obwodzie (ponieważ zostały narysowane z małą precyzją). Spowoduje to generowanie długich i wąskich nakładek, na których mnie nie obchodzi.

Ale innym razem będą istnieć duże nakładki solidnych wielokątów, co oznacza duże obszary, w których wielokąt sąsiedztwa nakłada się na inny. Chcę wybrać tylko te.

Zobacz zdjęcie poniżej nakładających się na siebie. Wyobraź sobie, że chciałem wybrać tylko niebieski wielokąt w lewym dolnym rogu.

Mogłem patrzeć na obszary, ale czasami wąskie są tak długie, że ostatecznie mają obszary tak duże jak niebieski wielokąt. Próbowałem zrobić stosunek powierzchni do obwodu, ale to również dało mieszane wyniki.

Próbowałem nawet użyć ST_MinimumClearance, ale czasami duże obszary będą miały przymocowaną wąską część lub dwa bardzo bliskie wierzchołki.

Jakieś pomysły na inne podejścia?

Ostatecznie najlepiej działało mi użycie bufora ujemnego, jak sugerują poniższe @Cyril i @FGreg.

Użyłem czegoś takiego jak:

ST_Area(ST_Buffer(geom, -10)) as neg_buffer_area

W moim przypadku jednostkami były metry, więc 10 m bufor ujemny.

W przypadku wąskich wielokątów obszar ten zwrócił zero (geometria byłaby pusta). Następnie użyłem tej kolumny do odfiltrowania wąskich wielokątów.

Spróbowałbym utworzyć bufor ujemny, jeśli zjada cienkie wielokąty, to dobrze, jeśli nie je wielokąta, to jest mój ... :-)

uruchom ten skrypt, wcześniej ustawiając 2/3 szerokości wielokątów liniowych ...

create table name_table as
SELECT ST_Buffer(
(ST_Dump(
(ST_Union(
ST_Buffer(
(geom),-0.0001))))).geom,
0.0001)) as geom from source_table

System operacyjny: -) ...

Zamiast obszaru / obwodu lepiej jest użyć obszaru podzielonego przez kwadrat obwodu (lub jego odwrotność).

Nazywa się to również „wskaźnikiem kształtu”. Kwadrat obwodu podzielonego przez obszar ma minimalną wartość 4 * Pi () (w przypadku dysku, który jest najbardziej kompaktową geometrią 2D), więc można go znormalizować za pomocą 4 * Pi () interpretacja (znormalizowane wartości bliskie 1 oznaczają wtedy, że masz bardzo zwarte obiekty, a kwadraty mają wartości około 1,27).

EDYCJA: Próg na tym obszarze byłby przydatny do usunięcia bardzo małych artefaktów, które mogłyby być zwarte. Wtedy wskaźnik kształtu pokazałby lepszy kontrast. EDYCJA: oprócz tej odpowiedzi użycie ST_Snap może pomóc w rozwiązaniu problemu przed jego wystąpieniem.

Jedną z opcji byłoby użycie stosunku powierzchni wielokąta do najdłuższej linii, którą można narysować za pomocą jego kończyn. Identyfikacja długich wąskich wielokątów.

select * from polygons where ST_Length(ST_LongestLine(geom, geom)) < ST_Area(geom) * 4

Działa to całkiem dobrze w przypadku wielokątów płaskich. Możesz dopasować współczynnik (przez pomnożenie obszaru) do własnych potrzeb i projekcji.

Wygląda na to, że może to pasować do Twojego przypadku użycia: Wyeliminuj wybrane wielokąty

Łączy wybrane wielokąty warstwy wejściowej z pewnymi sąsiadującymi wielokątami, usuwając ich wspólną granicę. Przylegający wielokąt może być albo tym, który ma największy lub najmniejszy obszar, albo tym, który dzieli największą wspólną granicę z wielokątem, który ma zostać wyeliminowany.

Eliminacja jest zwykle stosowana do pozbycia się wielokątów płaskich, tj. Małych wielokątów, które są wynikiem procesów przecinania się wielokątów, w których granice danych wejściowych są podobne, ale nie identyczne.

Wygląda na to, że chcesz wypróbować opcję „Największa wspólna granica”.

Wygląda mi to na idealny przypadek użycia dla rozszerzenia topologii PostGIS . Parametr tolerancji topologii określa, jak daleko pozwalasz wierzchołkom przyciągać się do innych istniejących wielokątów, radzić sobie z niską precyzją danych źródłowych i czyścić je.

Krótko mówiąc, strategia jest następująca:

1. Włącz rozszerzenie topologii

CREATE EXTENSION postgis_topology;

2. Utwórz nową pustą topologię

SELECT topology.CreateTopology('neighborhoods_topo', 4326, 1e-7);

Trzecim parametrem jest tolerancja w jednostkach CRS; wybierz mądrze. Idealnie, potrzebujesz CRS, gdzie jednostką jest metry. Jeśli jednostka CRS nie jest licznikiem, tak jak w przypadku WGS 84 aka 4326, użyj ST_Transformponownie rzutować swoje wielokąty.

3. Dodaj kolumnę TopoGeometry do tabeli wielokątów

SELECT topology.AddTopoGeometryColumn('neighborhoods_topo', 'public', 'neighborhoods', 'topogeom', 'POLYGON');

To zwraca nowy layer_id. Zapisz to, będzie potrzebne później. Będzie to warstwa, 1jeśli zaczniesz od zera, i będzie rosła przy każdym nowym połączeniu.

4. Dodaj wszystkie wielokąty do topologii

UPDATE public.neighborhoods
SET topogeom = topology.toTopoGeom(geom, 'neighborhoods_topo', 1, 1e-7);

Duży zestaw danych może potrwać kilka godzin, bądź cierpliwy. 1jest zwrócony wcześniej identyfikator_warstwy.

5. Znajdź twarze pojawiające się w kilku dzielnicach

Znajdź wszystkie twarze z topologii, które są obecne w 2 lub więcej topogeometrii. Zostawię zapytanie jako ćwiczenie. Najłatwiej jest prawdopodobnie z GetTopoGeomElementsfunkcją, a następnie pogrupuj według identyfikatora twarzy i spójrz na te z liczbą 2 lub więcej. Alternatywnie, możesz utworzyć nową tabelę z oczyszczoną geometrią z kolumny topogeom, po prostu rzutuj ją na standardową geometrię topogeom::geometryi powtórz to, co już masz, ale teraz z czystym zestawem danych bez nakładania się taśmy.