¿Usa ST_Difference para eliminar características superpuestas?

11

Estoy tratando de usar ST_Difference para crear un conjunto de polígonos (processing.trimmedparcelsnew) que no contienen ninguna parte del área cubierta por otro conjunto de polígonos (test.single_geometry_1) usando PostGis 2.1 (y Postgres SQL 9.3). Aquí está mi consulta:

CREATE TABLE processing.trimmedparcelsnew AS
SELECT
    orig.id, ST_Difference(orig.geom, cont.geom) AS difference
FROM 
    test.single_geometry_1 cont,
    test.multi_geometry_1 orig;

Pero los polígonos resultantes no se han recortado, sino que parecen haberse dividido donde se cruzan con la otra capa. He intentado ejecutar el select sin poner el resultado en una tabla y todo lo demás que se me ocurre, pero parece que no puedo hacer que esta función funcione.

He adjuntado una foto del resultado

ingrese la descripción de la imagen aquí


Después de los comentarios, he intentado agregar una cláusula WHERE. Quiero las parcelas que no tienen intersecciones, y las áreas de intersección de las otras parcelas eliminadas (la prueba de capa. Single_geometry representa la contaminación que quiero eliminar de mis parcelas). Intenté una intersección pero, por supuesto, realmente quiero las no intersecciones, así que ahora estoy intentando una disyunción. También intenté agregar el origen a mi tabla, pero la documentación de ST_Difference ( http://postgis.net/docs/ST_Difference.html ) dice que devuelve la geometría exacta que necesito (una geometría que representa esa parte de la geometría A que no se cruza con la geometría B), por lo que estoy confundido sobre por qué querría el polígono original en mi tabla. De todos modos, aquí está mi código modificado:

CREATE TABLE processing.trimmedparcelsnew AS
SELECT
    orig.id, ST_Difference(orig.geom, cont.geom) AS difference, orig.geom AS geom
FROM 
    test.single_geometry_1 cont,
    test.multi_geometry_1 orig
WHERE ST_Disjoint(orig.geom, cont.geom);

Siguiendo la respuesta de Dbaston, ahora he intentado:

CREATE TABLE processing.parcels_trimmed AS
SELECT id, COALESCE(ST_Difference(geom, (SELECT ST_Union(b.geom) 
                                         FROM test.single_geometry_1 b
                                         WHERE ST_Intersects(a.geom, b.geom)
                                         AND a.id != b.id)), a.geom)
FROM test.multi_geometry_1 a;

El resultado de esto es solo una copia de test.multi_geometry_1. Aunque ahora la división ya no se produce.

Probé la versión anterior, pero nuevamente recibí una copia de test.multi_geometry_1:

CREATE TABLE processing.parcels_trimmed_no_coalesce AS
SELECT id, COALESCE(ST_Difference(geom, (SELECT ST_Union(b.geom) 
                                         FROM test.single_geometry_1 b
                                         WHERE ST_Intersects(a.geom, b.geom)
                                         AND a.id != b.id)), a.geom)
FROM test.multi_geometry_1 a;

Estoy empezando a preguntarme si hay algo más que estoy haciendo mal. La declaración que sigue es:

DROP TABLE IF EXISTS processing.parcels_trimmed_no_coalesce;

Y estoy ejecutando las consultas desde la ventana de consulta SQL de PostgreSQL y Openjump.

La declaración que uso para ver la tabla es:

SELECT * FROM processing.parcels_trimmed_no_coalesce;

En aras de la simplificación, ahora he reducido esta consulta a solo:

SELECT id, COALESCE(ST_Difference(geom, (SELECT ST_Union(b.geom) 
                                         FROM test.geometriestocutagainst b
                                         WHERE ST_Intersects(a.geom, b.geom)
                                         AND a.id != b.id)), a.geom)
FROM test.geometriestocut a;

Esto todavía da como resultado los polígonos originales (test.geometriestocut) cuando el resultado deseado es el original recortado contra test.geometriestocuta contra.

Mercado
fuente
No especificó una WHEREcláusula, por lo que puede tener una expansión polinómica en la tabla resultante. ¿Cuántas filas hay trimmedparcelsnew?
Vince
Si solo desea la diferencia donde se cruzan, puede intentar agregar WHERE ST_Intersects (orig.geom, cont.geom). De lo contrario, la diferencia de dos polígonos que no se cruzan es el polígono original.
John Powell
Hay 24 filas en la parcela recortada nueva, quiero la diferencia incluso cuando no se cruzan, entonces, ¿corregiría que necesito usar orig.geom en la tabla en lugar de la diferencia?
Martes
Una prueba disjunta debería producir una expansión polinómica: cada entidad en cont aparecerá una vez para cada entidad en origen que no se superponga, y la diferencia nunca cambiará una geometría de entrada
Vince
Ok, gracias por la aclaración, pero todavía no estoy seguro de por qué el código original no funciona. Si ST_Difference (orig.geom, cont.geom) devuelve las geometrías en a que no se cruzan con b, entonces ¿por qué la tabla contiene las geometrías divididas en lugar de las geometrías en a que no intersecan b?
Mart

Respuestas:

14

Una autounión le permite operar en la relación entre pares de dos características. Pero no creo que le interesen los pares: para cada característica, desea operar en la relación entre esa característica y todas las demás características en su conjunto de datos. Puede lograr esto con una expresión de subconsulta:

CREATE TABLE parcels_trimmed AS
SELECT id, ST_Difference(geom, (SELECT ST_Union(b.geom) 
                                FROM parcels b
                                WHERE ST_Intersects(a.geom, b.geom)
                                  AND a.id != b.id))
FROM parcels a;

Sin embargo, es posible que vea algo extraño en los resultados. ¡Las parcelas que no tienen superposiciones se están eliminando por completo! Esto se debe a que el ST_Unionagregado en un conjunto de registros vacío será NULL, y ST_Difference(geom, NULL)es NULL. Para suavizar esto, debe completar su ST_Differencellamada en COALESCE:

CREATE TABLE parcels_trimmed AS
SELECT id, COALESCE(ST_Difference(geom, (SELECT ST_Union(b.geom) 
                                         FROM parcels b
                                         WHERE ST_Intersects(a.geom, b.geom)
                                         AND a.id != b.id)), a.geom)
FROM parcels a;

Esto significa que si el resultado de ST_Differencees NULL, la expresión fusionada se evaluará a la geometría original.

La consulta anterior eliminará las áreas superpuestas de su dominio por completo. Si en su lugar desea elegir un ganador, puede hacerlo a.id < b.id, o algún otro criterio, en lugar de a.id != b.id.

dbaston
fuente
Gracias por la respuesta, desafortunadamente estoy teniendo problemas para hacer que esto funcione para mí, en cambio solo termino con el polígono original (a). Editaré mi pregunta con más información. Gracias de nuevo.
Martes
2

Yo tenía el mismo problema que tú. No sé si ya encontró una solución a su problema, pero modifiqué la respuesta aceptada arriba y obtuve lo que quería.

CREATE TABLE parcels_trimmed AS
SELECT id, COALESCE(ST_Difference(geom, (SELECT ST_Collect(b.geom) 
                                         FROM parcels b
                                         WHERE ST_Intersects(a.geom, b.geom)
                                         )), a.geom)
FROM parcels a;
Trozo de cuero
fuente
1

Yo uso ST_DifferenceAgg () de los complementos PostGIS . Debe fusionar las dos tablas juntas, tener un identificador único y un índice en la columna de geometría. Aquí hay un breve ejemplo:

WITH overlappingtable AS (
  SELECT 1 id, ST_GeomFromText('POLYGON((0 1, 3 2, 3 0, 0 1), (1.5 1.333, 2 1.333, 2 0.666, 1.5 0.666, 1.5 1.333))') geom
  UNION ALL
  SELECT 2 id, ST_GeomFromText('POLYGON((1 1, 3.8 2, 4 0, 1 1))')
  UNION ALL
  SELECT 3 id, ST_GeomFromText('POLYGON((2 1, 4.6 2, 5 0, 2 1))')
  UNION ALL
  SELECT 4 id, ST_GeomFromText('POLYGON((3 1, 5.4 2, 6 0, 3 1))')
  UNION ALL
  SELECT 5 id, ST_GeomFromText('POLYGON((3 1, 5.4 2, 6 0, 3 1))')
)
SELECT a.id, ST_DifferenceAgg(a.geom, b.geom) geom
FROM overlappingtable a,
     overlappingtable b
WHERE a.id = b.id OR -- Make sure to pass at least once the polygon with itself
      ((ST_Contains(a.geom, b.geom) OR -- Select all the containing, contained and overlapping polygons
        ST_Contains(b.geom, a.geom) OR
        ST_Overlaps(a.geom, b.geom)) AND
       (ST_Area(a.geom) < ST_Area(b.geom) OR -- Make sure bigger polygons are removed from smaller ones
        (ST_Area(a.geom) = ST_Area(b.geom) AND -- If areas are equal, arbitrarily remove one from the other but in a determined order so it's not done twice.
         a.id < b.id)))
GROUP BY a.id
HAVING ST_Area(ST_DifferenceAgg(a.geom, b.geom)) > 0 AND NOT ST_IsEmpty(ST_DifferenceAgg(a.geom, b.geom));

Esto fusionará las partes superpuestas con el polígono superpuesto más grande. Si desea mantener la parte superpuesta separada, mire el ejemplo ST_splitAgg ().

Pierre Racine
fuente