Non-multipolygon relations get bounding box as geometry instead of actual geometry

[laurentpellegrino]

Description

When processing a planet PBF with ohsome-planet contributions, all non-multipolygon relation types receive their bounding box as geometry instead of their actual assembled geometry, even with --multipolygon-member-limit -1.

For example, relation 2195878 (Río Aragón,type=waterway) outputs a rectangular polygon matching its bbox coordinates exactly:

POLYGON ((-1.7876927 42.2204245, -1.7876927 42.8105252, -0.5049919 42.8105252, -0.5049919 42.2204245, -1.7876927 42.2204245))

This affects a large number of relations. On a recent planet extract (~1.45 billion features), approximately 5.4 million relations have exactly 5 points (bounding box rectangle). Affected relation types include waterway, route, boundary (non-multipolygon), and others.

Expected behavior

Relations should have their actual geometry assembled from member ways/nodes, similar to how multipolygon relations are handled via --multipolygon-member-limit.

Steps to reproduce

ohsome-planet contributions --pbf planet-latest.osm.pbf --data output --multipolygon-member-limit -1

Then query the output:

-- Show bbox-only relations by type
SELECT tags['type'] as relation_type, COUNT(*) as cnt
FROM 'output/contributions/latest/*.parquet'
WHERE osm_type = 'relation' AND ST_NPoints(geometry) = 5
GROUP BY relation_type ORDER BY cnt DESC;

For context, setting --multipolygon-member-limit to -1 (unlimited) only enables geometry assembly for multipolygon relations. It seems there is currently no option to enable geometry assembly for other relation types (waterway, route, boundary, etc.).

[Hagellach37]

Hey @laurentpellegrino ,
thanks for reporting this. This points to an interesting aspect of ohsome-planet's geometry building and let's discuss the details here and what behaviour you would expect from the tool.

The --multipolygon-member-limit is only coming into effect for relations for which the geometry represents (or should represent a multipolygon).

The example (https://www.openstreetmap.org/relation/2195878) that you have posted indeed returns a bounding box geometry, because this is a non-multipolygon relation. In your example this is maybe a multilinestring.

For relations in OSM all sorts of geometry types can be combined and these will usually result in a GeometryCollection.

Currently, ohsome-planet provides a bounding box geometry for these GeometryCollections in the geometry attribute of the parquet files. But you can also extract the individual relation member geometries from the members column(s). We've put together a simple SQL example, how to extract these geometries here: https://github.com/GIScience/ohsome-planet/blob/main/docs/useful_queries.md#extract-relation-member-geometries

The main question for us at the moment: Is it useful to populate the geometry attribute of the parquet files with a GeometryCollection? We find it often difficult to make sense of these geometry collections without the extra information that the members have.

[rtroilo]

Hey @laurentpellegrino ,

Let's discuss this together. We are very open to further adjustments to the schema.
OSM relations are always problematic and cause us a lot of headaches.
For relations that are not multipolygons, we currently replace the actual geometry with the bounding box of that relation. However, you can determine the actual type in the geometry_type column.
This means that geometry_type = ‘GeometryCollection’ results in a bounding box geometry (polygon) instead of a GeometryCollection.
Why do we do this? GeometryCollection are difficult to handle in themselves, qgis has problems with them and you lose all the role information of the members as well. We often faced the problem, for example, that we wanted to analyze a bus line, but were only interested in the lines and not the stops and platforms.
With duckdb and its mighty list functions we could easily do for example something like that, as @Hagellach37 already mentioned:

select 
  osm_type, 
  osm_id, 
  tags, 
  ST_Collect( -- creates a MULTILINESTRING
    list_transform(
      list_filter(members,  m -> m.geometry_type = 'LineString'),  -- filter only members which are LineStrings
      m -> ST_GeomFromWKB(m.geometry)) -- convert and select the geometry of those members
  ) as geometry from read_parquet('ohsome-planet.parquet') 
where geometry_type = 'GeometryCollection'  and map_contains_entry(tags, 'route', 'bicycle')

Of course, you need quite a bit of knowledge to do this. And it's not exactly beginner-friendly. And it seems strange why the geometry is suddenly just a bounding box instead of a collection. But OSM relations, and especially non-polygon relations, are complicated in themselves.
So I'm curious to hear what you think about it.

[laurentpellegrino]

Thanks for your answers and help.

I created this issue because we received reports from a few users who thought there was a bug when a bounding box geometry appears instead of a GeometryCollection in the GeoParquet files we publish:

https://openplanetdata.com/datasets/openstreetmap/

I now better understand the context and, since there is a workaround, I will extract the individual relation member geometries from the members column when needed for now.