La intersección entre dos multipolígonos produce un objeto GeometryCollection anómalo lleno de LineString's y Polygon's (intentando obtener el área de intersección)

Question

Normalmente cuando llamo al método #intersección en un multipolígono para encontrar la forma de intersección en otro multipolígono obtengo un multipolígono de vuelta. Eso o algo más que puedo llamar a la función #área.

Pues bien, 150 uniones espaciales en una lista de shapefiles de zonificación (es decir, este es un caso raro) estoy escudriñando mi código arroja un error porque no puede llamar al método #area en un CAPIGeometryCollection . Estoy teniendo problemas para vincular qgis con mi db postgresql y no sé cómo crear un kml para un objeto "GeometryCollection"

Este es el aspecto de la salida de la cadena WKT para la forma de intersección:

"GEOMETRYCOLLECTION (LINESTRING (1480424.245856002 596629.5658560097, 1480362.875856012 596566.8158560097), LINESTRING (1480362.875856012 596566.8158560097, 1480293.995856002 596496.4358560145), LINESTRING (1480293.995856002 596496.4358560145, 1480235.625856012 596435.0658560097), LINESTRING (1480235.625856012 596435.0658560097, 1480186.495856002 596385.8158560097), LINESTRING (1480186.495856002 596385.8158560097, 1480146.125856012 596345.3758560121), LINESTRING (1480146.125856012 596345.3758560121, 1480141.125856012 596341.3158560097), LINESTRING (1480141.125856012 596341.3158560097, 1480090.995856002 596301.3158560097), LINESTRING (1480090.995856002 596301.3158560097, 1480044.875856012 596263.8758560121), LINESTRING (1480044.875856012 596263.8758560121, 1480005.745856002 596231.995856002), POLYGON ((1480512.6034600246 596678.2322904326, 1480512.598981008 596678.1839810014, 1480511.7765895594 596677.7314431852, 1480512.6034600246 596678.2322904326)), POLYGON ((1480005.745856002 596231.995856002, 1479874.7139810026 596136.0839810073, 1479816.8752310127 596088.7502310127, 1479758.3752310127 596040.8127310127, 1479707.547731012 595999.19689875, 1479758.375856012 596040.8158560097, 1479816.875856012 596088.8158560097, 1479875.375856012 596136.6258560121, 1480005.745856002 596231.995856002)))"

Podría omitir este error con un rescue bloque o un next identificador y seguir adelante por caso pero eso no es muy profesional. La documentación sobre este tema es extremadamente escasa, así que pensé en contribuir con una pregunta de stackexchange

Además (acabo de comprobarlo), no puedo cargar GEOMETRYCOLLECTION wkts en qgis (actualmente pasando a mano los polígonos y las cadenas de líneas)

(nota: alguien con >300 puntos o lo que sea que marque esta pregunta con estas nuevas etiquetas de pregunta: 'geometry-collection' & 'multi-polygon', son bastante básicos para gis creo, esp los multipolígonos, gracias)

EDITAR: Podría sacar los polígonos de estas colecciones, pero me pregunto qué significan estas cadenas de líneas (es decir, tal vez son los anillos exteriores de los polígonos que forman parte de la forma de intersección)

Answer 1

Solución por ahora: ignorar los LineStrings, estos son probablemente producidos por una unión espacial de tipo "st_touches" donde sólo se tocan los perímetros de los multipolígonos

RGeo te permite separar las piezas de una colección de geometría, así que escribí un código para atrapar la forma de intersección si es del tipo GeometryCollection, enviándola a un área diferente para extraer los polígonos, sumar su tamaño, marcarlo como el tamaño de la intersección, y seguir adelante

    x = shape.intersection(p.proj_shape_2264)
    if x.geometry_type.to_s == "GeometryCollection"
      area = 0
      (0..(x.count - 1)).to_a.each do |k| 
        collection_shape = x.geometry_n(k)
        next if collection_shape.geometry_type.to_s == "LineString"
        area  += collection_shape.area            
      end
      z.intersect_size = area
    else

Answer 2

Recientemente tuve dos geometrías cuya intersección produjo un GeometryCollection de Polygon y diminutos LineString . Mi "arreglo" fue erosionar ligeramente la geometría de origen antes de la intersección:

shape.buffer(1e-5).intersection(p.proj_shape_2264)

Ahora la intersección produce sólo una Polygon .

Actualización : Me encontré con esto de nuevo en una situación en la que la aplicación requiere que la geometría no sea erosionada. Escribí algo similar a Respuesta de boulder_ruby , pero sólo paso por Polígono y Multipolígono. Esta es mi implementación:

  # Returns geometry of type Polygon or MultiPolygon and never a
  # GeometryCollection for the given geometry.
  def as_polygonal(geom)
    # non-intersecting geometries will return an empty
    # GeometryCollection... we only want to "fix" things which have
    # linestrings inside the intersection.
    if geom.geometry_type == RGeo::Feature::GeometryCollection and geom.count > 0
      acceptable_types = [RGeo::Feature::Polygon, RGeo::Feature::MultiPolygon]
      polygonal = geom.select {|g| acceptable_types.include?(g.geometry_type) }

      geom = RGeo::Geos.factory(srid: 4326).parse_wkt('POLYGON EMPTY')
      polygonal.each {|p| geom = geom.union(p)}
    end

    geom
  end