Oddziel wielokąty na podstawie przecięcia za pomocą PostGIS

Mam tabelę wielokątów PostGIS, w której niektóre przecinają się ze sobą. Oto co próbuję zrobić:

• Dla danego wielokąta wybranego przez id, podaj mi wszystkie wielokąty, które się przecinają. Gruntownie,select the_geom from the_table where ST_Intersects(the_geom, (select the_geom from the_table where source_id = '123'))
• Z tych wielokątów muszę utworzyć nowe wielokąty, aby przecięcie stało się nowym wielokątem. Więc jeśli wielokąt A przecina się z wielokątem B, otrzymam 3 nowe wielokąty: A minus AB, AB i B minus AB.

Jakieś pomysły?

Ponieważ powiedziałeś, że masz grupę przecinających się wielokątów dla każdego interesującego cię wielokąta, możesz chcieć stworzyć coś, co nazywa się „nakładką wielokątów”.

To nie jest dokładnie to, co robi rozwiązanie Adama. Aby zobaczyć różnicę, spójrz na to zdjęcie skrzyżowania ABC:

Wierzę, że rozwiązanie Adama stworzy wielokąt „AB” obejmujący zarówno obszar „AB! C” i „ABC”, jak i wielokąt „AC” obejmujący „AC! B” i „ABC” oraz „ Wielokąt BC, czyli „BC! A” i „ABC”. Tak więc wielokąty wyjściowe „AB”, „AC” i „BC” nakładałyby się na obszar „ABC”.

Nakładka wielokąta wytwarza wielokątów, więc AB! C byłby jednym wielokątem, a ABC byłby jednym wielokątem.

Tworzenie nakładki wielokątów w PostGIS jest w rzeczywistości dość proste.

Istnieją w zasadzie trzy kroki.

Krok 1 to wyodrębnienie szkicu [Zauważ, że używam zewnętrznego pierścienia wielokąta, robi się to trochę bardziej skomplikowane, jeśli chcesz poprawnie obsługiwać otwory]:

SELECT ST_ExteriorRing(polygon_col) AS the_geom FROM my_table) AS lines

Krok 2 polega na „węzłowaniu” linii (tworzenie węzła na każdym skrzyżowaniu). Niektóre biblioteki, takie jak JTS, mają klasy „Noder”, których możesz użyć, ale w PostGIS funkcja ST_Union robi to za Ciebie:

SELECT ST_Union(the_geom) AS the_geom FROM (...your lines...) AS noded_lines

Krok 3 polega na utworzeniu wszystkich możliwych nienakładających się wielokątów, które mogą pochodzić ze wszystkich tych linii, wykonanych przez funkcję ST_Polygonize :

SELECT ST_Polygonize(the_geom) AS the_geom FROM (...your noded lines...)

Możesz zapisać dane wyjściowe każdego z tych kroków w tabeli tymczasowej lub połączyć je wszystkie w jedną instrukcję:

CREATE TABLE my_poly_overlay AS
SELECT geom FROM ST_Dump((
    SELECT ST_Polygonize(the_geom) AS the_geom FROM (
        SELECT ST_Union(the_geom) AS the_geom FROM (
            SELECT ST_ExteriorRing(polygon_col) AS the_geom FROM my_table) AS lines
        ) AS noded_lines
    )
)

Używam ST_Dump, ponieważ dane wyjściowe ST_Polygonize są kolekcją geometrii i (zazwyczaj) wygodniej jest mieć tabelę, w której każdy wiersz jest jednym z wielokątów tworzących nakładkę wieloboku.

Jeśli dobrze rozumiem, chcesz wziąć (A jest lewą geometrią, B jest prawą):

Obraz A∪B http://img838.imageshack.us/img838/3996/intersectab1.png

I wyodrębnij:

A ∖ AB

Zdjęcie A ∖ AB http://img830.imageshack.us/img830/273/intersectab2.png

AB

Obraz AB http://img828.imageshack.us/img828/7413/intersectab3.png

i B ∖ AB

Obraz B ∖ AB http://img839.imageshack.us/img839/5458/intersectab4.png

To znaczy - trzy różne geometrie dla każdej przecinającej się pary.

Najpierw utwórzmy widok wszystkich przecinających się geometrii. Zakładając, że twoja nazwa tabeli to polygons_table, użyjemy:

CREATE OR REPLACE VIEW p_intersections AS    -- Create a view with the 
SELECT t1.the_geom as t1_geom,               -- intersecting geoms. Each pair
       t2.the_geom as t2_geom                -- appears once (t2.id<t2.id)
    FROM polygons_table t1, polygons_table t2  
         WHERE t1.id<t2.id AND t1.the_geom && t2.the_geom 
                           AND intersects t1.the_geom, t2.the_geom;

Teraz mamy widok (praktycznie tabelę tylko do odczytu), która przechowuje pary przecinających się geom, gdzie każda para pojawia się tylko raz ze względu na t1.id<t2.idwarunek.

Teraz zebrać skrzyżowań - A∖AB, ABi B∖AB, przy użyciu SQL na UNIONwszystkie trzy pytania:

--AB:
SELECT ST_intersection(t1.the_geom, t2.the_geom) 
    AS geom 
    FROM p_intersections

UNION 

--A∖AB:
SELECT ST_Difference(t1.the_geom, t2.the_geom) 
    AS geom 
    FROM p_intersections

UNION

--B∖AB:
SELECT ST_Difference(t2.the_geom, t1.the_geom) 
    AS geom 
    FROM p_intersections;

Uwagi:

1. &&Operator jest używany jako filtr, przed intersectsoperatorem w celu poprawy wydajności.
2. Zdecydowałem się utworzyć VIEWzamiast jednego gigantycznego zapytania; Jest to wyłącznie dla wygody.
3. Jeśli chodziło Ci ABo połączenie, a nie przecięcie, Ai B- Użyj ST_Union zamiast st_intersection w UNIONzapytaniu w odpowiednich miejscach.
4. ∖Znak jest znak Unicode dla różnicy zbiorów; usuń go z kodu, jeśli spowoduje to pomyłkę w bazie danych.
5. Zdjęcia dzięki uprzejmości ładnej kategorii Teoria zbiorów .

To, co opisujesz, to sposób, w jaki operator unijny działa w ArcGIS, ale jest trochę inny niż Union lub skrzyżowanie w świecie GEOS. Podręcznik Shapely zawiera przykłady działania zestawów w GEOS . Jednak wiki PostGIS ma dobry przykład użycia linii, która powinna załatwić sprawę.

Alternatywnie możesz obliczyć trzy rzeczy:

1. ST_Intersection (geom_a, geom_b)
2. ST_Difference (geom_a, geom_b)
3. ST_Difference (geom_b, geom_a)

To powinny być trzy wielokąty, o których wspomniałeś w drugim punkcie kuli.

Coś jak:

INSERT INTO new_table VALUES ((wybierz id, the_geom ze old_table gdzie st_intersects (the_geom, (wybierz the_geom ze old_table gdzie id = '123')) = true

EDYCJA: potrzebujesz rzeczywistego przecięcia wielokąta.

INSERT INTO new_table values ​​((wybierz id, ST_Intersection (the_geom, (wybierz the_geom ze starego, gdzie id = 123))

zobacz, czy to się uda.