arreglar los polígonos que se intersectan a sí mismos usando ogr2ogr u otra librería de código abierto

Question

Tengo un shapefile que incluye algunos polígonos que se intersectan entre sí. ¿Cuál sería la forma más fácil de arreglarlos usando el código abierto Python/command line? ¿Es posible utilizar ogr2ogr para hacerlo (las respuestas utilizando fiona, shapely, geopandas, etc. también son bienvenidas)? Quiero hacer esto porque el motor de Google Earth no puede manejar polígonos auto-interseccionados.

La respuesta usando spatialite parece una buena opción, sin embargo mi GDAL no soporta actualmente las funciones de spatialite. He construido GDAL usando Docker.

ogrinfo --formats sqlite | grep 'spatialite' -i

Answer 1

Spatiallite tiene un Método MakeValid que debería solucionar los problemas con los polígonos no válidos. Puede ejecutarlo contra su conjunto de datos utilizando la función sql parámetro en ogr2ogr así:

ogr2ogr -dialect SQLITE \
-sql "select PFAF_ID, ST_MakeValid(Geometry) as Geometry from hybas_lev06_v1c_merged_fiona_upstream_downstream_FAO_V01" \
valid_corrected.shp \
hybas_lev06_v1c_merged_fiona_upstream_downstream_FAO_V01.shp

Answer 2

Ya que mencionaste que tenías acceso a docker, puedes utilizar esta imagen postgis para limpiar/exportar sus geometrías en un solo docker run comando:

docker run -it --rm -e POSTGRES_PASSWORD=mysecretpassword \
 --mount type=bind,source=<host_shapefile_location>,target=/tmp/shapedir \
 mdillon/postgis bash -c "
apt-get update && apt-get install gdal-bin && tee /docker-entrypoint-initdb.d/z.sh <<EOF
psql -c $'CREATE DATABASE ogr WITH TEMPLATE = 'template_postgis''
ogr2ogr -f 'PostgreSQL' -a_srs 'EPSG:<YOUR_SRID>' PG:'host=localhost user=postgres dbname=ogr password=mysecretpassword' /tmp/shapedir/myshp.shp -nlt MULTIPOLYGON -lco GEOMETRY_NAME=geom -nln myshp
ogr2ogr -overwrite -f 'ESRI Shapefile' /tmp/shapedir/output.shp PG:'host=localhost user=postgres dbname=ogr password=mysecretpassword' -sql 'SELECT myfield, st_collectionextract(st_makevalid(st_buffer(st_makevalid(geom),0)),3) FROM myshp'
EOF
/docker-entrypoint.sh postgres"

Este comando hace lo siguiente:

1. Bind monta el directorio del host que contiene el shapefile en /tmp/shapedir en el contenedor.
2. Instala GDAL "vainilla" en el contenedor.
3. Crea una base de datos Postgis adecuada para procesar el shapefile.
4. Utiliza ogr2ogr para importar el shapefile en la base de datos recién creada (y coacciona la geometría a MULTIPOLYGON en aras de la simplicidad y la cordura)
5. Limpia las geometrías mediante st_makevalid(geom) y st_buffer(geom,0) y luego aplica st_collectionextract(geom,3) para garantizar que sólo se devuelvan los multipolígonos.
6. Exporta la geometría limpiada al directorio anfitrión especificado.

Una vez que su shapefile limpio se haya exportado con éxito, sólo tiene que ctrl +c en su terminal para detener y eliminar el contenedor temporal.

Nota: Asegúrese de proporcionar explícitamente el código EPSG correcto a través de a_srs -OGR no es bueno para autodetectarlo. Puede buscar la proyección correcta pegando el contenido del archivo .prj aquí