ST_Union falla con TopologyException con polígonos válidos y ST_SnapToGrid

Question

Estoy trabajando con Postgres 9.6 y PostGIS 2.2.2. Estas son mis tablas:

                                     Table "public.allparcels"
  Column  |            Type             |                        Modifiers
----------+-----------------------------+----------------------------------------------------------
 gid      | integer                     | not null default nextval('allparcels_gid_seq'::regclass)
 the_geog | geometry(MultiPolygon,4326) |
 county   | character varying           |

                                       Table "public.parcels"
    Column    |          Type           |                            Modifiers
--------------+-------------------------+------------------------------------------------------------------
 ogc_fid      | integer                 | not null default nextval('parcels_ogc_fid_seq'::regclass)
 wkb_geometry | geometry(Polygon,4326)  |
 county       | character varying       |

Cuando ejecuto esto (a través de psycopg2):

INSERT INTO allparcels(county, the_geog) 
SELECT 'Leeds', ST_MakeValid(ST_Union(wkb_geometry)) FROM parcels WHERE county='Leeds'

Me sale este error:

Traceback (most recent call last):
  File "combine.py", line 21, in <module>
    cursor.execute(q, (county, county))
psycopg2.InternalError: GEOSUnaryUnion: TopologyException: Input geom 1 is invalid: Self-intersection at or near point -1.394465365015132 53.741246412906918 at -1.394465365015132 53.741246412906918

Pensé que debía tener geometrías subyacentes no válidas, pero no veo ninguna:

SELECT * from parcels WHERE NOT ST_isvalid(wkb_geometry);

no produce ningún resultado. Como muestra la consulta, ya estoy ejecutando ST_MakeValid para asegurarse de que el polígono unido es válido antes de insertarlo.

¿Qué estoy haciendo mal?

ACTUALIZACIÓN:

Esta es la última consulta que he probado. La depuración de la consulta indica que es la ST_Union la que falla. Disminuir la precisión (incluso a 0,01) no ayuda:

INSERT INTO allparcels(county, the_geog) 
SELECT 'Leeds', ST_MakeValid(ST_Union(ST_SnapToGrid(wkb_geometry, .000001))) 
FROM parcels WHERE county='Leeds'

Answer 1

Esto ocurre a menudo con ST_Intersection, independientemente de si se utiliza ST_Snap (que es más útil para asegurar una cierta precisión que para corregir errores de geometría) y ST_MakeValid. El problema tiene que ver con el hecho de que cuando se intersectan polígonos, a menudo se encontrarán en un solo punto o a lo largo de una línea, además de producir una (Multi)intersección de polígonos. Es decir, cuando se intersectan dos polígonos cualesquiera, pueden tener múltiples intersecciones, no todas las cuales serán poligonales. En consecuencia, el tipo de datos para esa intersección particular será una GeometryCollection. Al intentar insertar esto en un campo MultiPolígono o Polígono, verá errores sobre la intersección no codificada o las auto-intersecciones. La forma más sencilla que he encontrado para solucionar este problema es utilizar ST_CollectionExtract que permite extraer de la GeometryCollection sólo los tipos Puntos, Cadenas de Líneas o Polígonos. En su caso, el uso del parámetro 3 (para polígonos) y la eliminación de ST_MakeValid debería solucionarlo:

INSERT INTO allparcels(county, the_geog) 
SELECT 'Leeds', ST_CollectionExtract(ST_Union(wkb_geometry), 3) 
 FROM parcels 
WHERE county='Leeds';

El error que estás viendo es consistente con que todos tus polígonos sean válidos en primer lugar, como afirmas. Por lo general, creo que es mejor ejecutar ST_IsValid para comprobar y en el caso de los errores, ST_MakeValid, en las geometrías antes de hacer cualquier intersección - como usted ha hecho. ST_MakeValid y ST_SnapToGrid tienen otros usos, pero generalmente no son la herramienta adecuada para arreglar los problemas de recogida de geometrías dentro de una consulta.