aio_overpass.element

  1"""Typed result set members."""
  2
  3import math
  4import re
  5from collections import defaultdict
  6from collections.abc import Iterable, Iterator
  7from dataclasses import dataclass
  8from typing import Any, Final, Generic, TypeAlias, TypeVar
  9
 10from aio_overpass import Query
 11from aio_overpass.spatial import GeoJsonDict, Spatial
 12
 13import shapely.geometry
 14import shapely.ops
 15from shapely.geometry import LinearRing, LineString, MultiPolygon, Point, Polygon
 16from shapely.geometry.base import BaseGeometry, BaseMultipartGeometry
 17
 18
 19__docformat__ = "google"
 20__all__ = (
 21    "collect_elements",
 22    "Element",
 23    "Node",
 24    "Way",
 25    "Relation",
 26    "Relationship",
 27    "Bbox",
 28    "GeometryDetails",
 29    "Metadata",
 30)
 31
 32
 33Bbox: TypeAlias = tuple[float, float, float, float]
 34"""
 35The bounding box of a spatial object.
 36
 37This tuple can be understood as any of
 38    - ``(s, w, n, e)``
 39    - ``(minlat, minlon, maxlat, maxlon)``
 40    - ``(minx, miny, maxx, maxy)``
 41"""
 42
 43
 44@dataclass(kw_only=True, slots=True)
 45class Metadata:
 46    """
 47    Metadata concerning the most recent edit of an OSM element.
 48
 49    Attributes:
 50        version: The version number of the element
 51        timestamp: Timestamp (ISO 8601) of the most recent change of this element
 52        changeset: The changeset in which the element was most recently changed
 53        user_name: Name of the user that made the most recent change to the element
 54        user_id: ID of the user that made the most recent change to the element
 55    """
 56
 57    version: int
 58    timestamp: str
 59    changeset: int
 60    user_name: str
 61    user_id: int
 62
 63
 64G = TypeVar("G", bound=BaseGeometry)
 65
 66
 67@dataclass(kw_only=True, slots=True)
 68class GeometryDetails(Generic[G]):
 69    """
 70    Element geometry with more info on its validity.
 71
 72    Shapely validity is based on an [OGC standard](https://www.ogc.org/standard/sfa/).
 73
 74    For MultiPolygons, one assertion is that its elements may only touch at a finite number
 75    of Points, which means they may not share an edge on their exteriors. In terms of
 76    OSM multipolygons, it makes sense to lift this requirement, and such geometries
 77    end up in the ``accepted`` field.
 78
 79    For invalid MultiPolygon and Polygons, we use Shapely's ``make_valid()``. If and only if
 80    the amount of polygons stays the same before and after making them valid,
 81    they will end up in the ``valid`` field.
 82
 83    Attributes:
 84        valid: if set, this is the original valid geometry
 85        accepted: if set, this is the original geometry that is invalid by Shapely standards,
 86                  but accepted by us
 87        fixed: if set, this is the geometry fixed by ``make_valid()``
 88        invalid: if set, this is the original invalid geometry
 89        invalid_reason: if the original geometry is invalid by Shapely standards,
 90                        this message states why
 91    """
 92
 93    valid: G | None = None
 94    accepted: G | None = None
 95    fixed: G | None = None
 96    invalid: G | None = None
 97    invalid_reason: str | None = None
 98
 99    @property
100    def best(self) -> G | None:
101        """The "best" geometry, prioritizing ``fixed`` over ``invalid``."""
102        return self.valid or self.accepted or self.fixed or self.invalid
103
104
105@dataclass(kw_only=True, repr=False, eq=False)
106class Element(Spatial):
107    """
108    Elements are the basic components of OpenStreetMap's data.
109
110    A query's result set is made up of these elements.
111
112    Objects of this class do not necessarily describe OSM elements in their entirety.
113    The degrees of detail are decided by the ``out`` statements used in an Overpass query:
114    using ``out ids`` f.e. would only include an element's ID, but not its tags, geometry, etc.
115
116    Element geometries have coordinates in the EPSG:4326 coordinate reference system,
117    meaning that the coordinates are (latitude, longitude) tuples on the WGS 84 reference ellipsoid.
118    The geometries are Shapely objects, where the x/y coordinates refer to lat/lon.
119    Since Shapely works on the Cartesian plane, not all operations are useful: distances between
120    Shapely objects are Euclidean distances for instance - not geodetic distances.
121
122    Tags provide the semantics of elements. There are *classifying* tags, where for a certain key
123    there is a limited number of agreed upon values: the ``highway`` tag f.e. is used to identify
124    any kind of road, street or path, and to classify its importance within the road network.
125    Deviating values are perceived as erroneous. Other tags are *describing*, where any value
126    is acceptable for a certain key: the ``name`` tag is a prominent example for this.
127
128    Objects of this class are not meant to represent "derived elements" of the Overpass API,
129    that can have an entirely different data structure compared to traditional elements.
130    An example of this is the structure produced by ``out count``, but also special statements like
131    ``make`` and ``convert``.
132
133    Attributes:
134        id: A number that uniquely identifies an element of a certain type
135            (nodes, ways and relations each have their own ID space).
136            Note that this ID cannot safely be used to link to a specific object on OSM.
137            OSM IDs may change at any time, e.g. if an object is deleted and re-added.
138        tags: A list of key-value pairs that describe the element, or ``None`` if the element's
139              tags not included in the query's result set.
140        meta: Metadata of this element, or ``None`` when not using ``out meta``
141        relations: All relations that are **also in the query's result set**, and that
142                   **are known** to contain this element as a member.
143
144    References:
145        - https://wiki.openstreetmap.org/wiki/Elements
146        - https://wiki.openstreetmap.org/wiki/Map_features
147        - https://wiki.openstreetmap.org/wiki/Overpass_API/Overpass_QL#out
148        - https://wiki.openstreetmap.org/wiki/Overpass_API/Permanent_ID
149    """
150
151    __slots__ = ()
152
153    id: int
154    tags: dict[str, str] | None
155    meta: Metadata | None
156    relations: list["Relationship"]
157
158    @property
159    def base_geometry(self) -> BaseGeometry | None:
160        """
161        The element's geometry, if available.
162
163        For the specifics, refer to the documentation of the ``geometry`` property in each subclass.
164        """
165        match self:
166            case Node() | Way() | Relation():
167                return self.geometry
168            case _:
169                raise NotImplementedError
170
171    @property
172    def base_geometry_details(self) -> GeometryDetails | None:
173        """More info on the validity of ``base_geometry``."""
174        match self:
175            case Way() | Relation():
176                return self.geometry_details
177            case Node():
178                return GeometryDetails(valid=self.geometry)
179            case _:
180                raise NotImplementedError
181
182    def tag(self, key: str, default: str | None = None) -> str | None:
183        """
184        Get the tag value for the given key.
185
186        Returns ``default`` if there is no ``key`` tag.
187
188        References:
189            - https://wiki.openstreetmap.org/wiki/Tags
190        """
191        if not self.tags:
192            return default
193        return self.tags.get(key, default)
194
195    @property
196    def type(self) -> str:
197        """The element's type: "node", "way", or "relation"."""
198        match self:
199            case Node():
200                return "node"
201            case Way():
202                return "way"
203            case Relation():
204                return "relation"
205            case _:
206                raise AssertionError
207
208    @property
209    def link(self) -> str:
210        """This element on openstreetmap.org."""
211        return f"https://www.openstreetmap.org/{self.type}/{self.id}"
212
213    @property
214    def wikidata_id(self) -> str | None:
215        """
216        [Wikidata](https://www.wikidata.org) item ID of this element.
217
218        This is "perhaps, the most stable and reliable manner to obtain Permanent ID
219        of relevant spatial features".
220
221        References:
222            - https://wiki.openstreetmap.org/wiki/Permanent_ID
223            - https://wiki.openstreetmap.org/wiki/Overpass_API/Permanent_ID
224            - https://wiki.openstreetmap.org/wiki/Key:wikidata
225            - https://www.wikidata.org/wiki/Wikidata:Notability
226            - Nodes on Wikidata: https://www.wikidata.org/wiki/Property:P11693
227            - Ways on Wikidata: https://www.wikidata.org/wiki/Property:P10689
228            - Relations on Wikidata: https://www.wikidata.org/wiki/Property:P402
229        """
230        # since tag values are not enforced, use a regex to filter out bad IDs
231        if self.tags and "wikidata" in self.tags and _WIKIDATA_Q_ID.match(self.tags["wikidata"]):
232            return self.tags["wikidata"]
233        return None
234
235    @property
236    def wikidata_link(self) -> str | None:
237        """This element on wikidata.org."""
238        if self.wikidata_id:
239            return f"https://www.wikidata.org/wiki/{self.wikidata_id}"
240        return None
241
242    @property
243    def geojson(self) -> GeoJsonDict:
244        """
245        A mapping of this object, using the GeoJSON format.
246
247        Objects are mapped as the following:
248         - ``Node`` -> ``Feature`` with optional ``Point`` geometry
249         - ``Way`` -> ``Feature`` with optional ``LineString`` or ``Polygon`` geometry
250         - ``Relation`` with geometry -> ``Feature`` with ``Polygon`` or ``MultiPolygon`` geometry
251         - ``Relation`` -> ``FeatureCollection`` (nested ``Relations`` are mapped to unlocated
252           ``Features``)
253
254        ``Feature`` properties contain all the following keys if they are present for the element:
255        ``id``, ``type``, ``role``, ``tags``, ``nodes``, ``bounds``, ``center``, ``timestamp``,
256        ``version``, ``changeset``, ``user``, ``uid``.
257        The JSON object in ``properties`` is therefore very similar to the original JSON object
258        returned by Overpass, skipping only ``geometry``, ``lat`` and ``lon``.
259
260        Additionally, features inside ``FeatureCollections`` receive the special ``__rel__``
261        property. Its value is an object containing all properties of the relation the collection
262        represents. This works around the fact that ``FeatureCollections`` have no ``properties``
263        member. The benefit is that you can take relation tags into consideration when styling
264        and rendering their members on a map (f.e. with Leaflet). The downside is that these
265        properties are duplicated for every feature.
266        """
267        if isinstance(self, Relation) and not self.geometry:
268            return {
269                "type": "FeatureCollection",
270                "features": [_geojson_feature(relship) for relship in self.members],
271            }
272
273        return _geojson_feature(self)
274
275    def __repr__(self) -> str:
276        return f"{type(self).__name__}({self.id})"
277
278
279@dataclass(kw_only=True, slots=True, repr=False, eq=False)
280class Node(Element):
281    """
282    A point in space, at a specific coordinate.
283
284    Nodes are used to define standalone point features (e.g. a bench),
285    or to define the shape or "path" of a way.
286
287    Attributes:
288        geometry: A Point or ``None`` if the coordinate is not included in the query's result set.
289
290    References:
291        - https://wiki.openstreetmap.org/wiki/Node
292    """
293
294    geometry: Point | None
295
296
297@dataclass(kw_only=True, slots=True, repr=False, eq=False)
298class Way(Element):
299    """
300    A way is an ordered list of nodes.
301
302    An open way is a way whose first node is not its last node (e.g. a railway line).
303    A closed way is a way whose first node is also its last node, and may be interpreted either
304    as a closed polyline (e.g. a roundabout), an area (e.g. a patch of grass), or both
305    (e.g. a roundabout surrounding a grassy area).
306
307    Attributes:
308        node_ids: The IDs of the nodes that make up this way, or ``None`` if they are not included
309                  in the query's result set.
310        bounds: The enclosing bounding box of all nodes, or ``None`` when not using ``out bb``.
311                The ``bounds`` property of Shapely geometries can be used as an alternative.
312        center: The center of ``bounds``, or ``None`` when not using ``out center``.
313                If you need a coordinate that is inside the element's geometry, consider Shapely's
314                ``representative_point()`` and ``centroid``.
315        geometry: A Linestring if the way is open, a LinearRing if the way is closed,
316                  a Polygon if the way is closed and its tags indicate that it represents an area,
317                  or ``None`` if the geometry is not included in the query's result set.
318        geometry_details: More info on the validity of ``geometry``.
319
320    References:
321        - https://wiki.openstreetmap.org/wiki/Way
322    """
323
324    node_ids: list[int] | None
325    bounds: Bbox | None
326    center: Point | None
327    geometry: LineString | LinearRing | Polygon | None
328    geometry_details: GeometryDetails[LineString | LinearRing | Polygon] | None
329
330
331@dataclass(kw_only=True, slots=True, repr=False, eq=False)
332class Relation(Element):
333    """
334    A relation is a group of nodes and ways that have a logical or geographic relationship.
335
336    This relationship is described through its tags.
337
338    A relation may define an area geometry, which may have boundaries made up of several
339    unclosed ways. Relations of ``type=multipolygon`` may have boundaries ("outer" role) and
340    holes ("inner" role) made up of several unclosed ways.
341
342    Tags describing the multipolygon always go on the relation. The inner and outer ways are tagged
343    if they describe something in their own right. For example,
344     - a multipolygon relation may be tagged as landuse=forest if it outlines a forest,
345     - its outer ways may be tagged as barrier=fence if the forest is fenced,
346     - and its inner ways may be tagged as natural=water if there is a lake within the forest
347       boundaries.
348
349    Attributes:
350        bounds: The bounding box of this element, or ``None`` when not using ``out bb``.
351                It encloses all node and way members; relations as members have no effect.
352                The ``bounds`` property of Shapely geometries can be used as an alternative.
353        center: The center of ``bounds``, or ``None`` when not using ``out center``.
354                If you need a coordinate that is inside the element's geometry, consider Shapely's
355                ``representative_point()`` and ``centroid``.
356        members: Ordered member elements of this relation, with an optional role
357        geometry: If this relation is deemed to represent an area, these are the complex polygons
358                  whose boundaries and holes are made up of the ways inside the relation. Members
359                  that are not ways, or are not part of any polygon boundary, are not part of the
360                  result geometry. This is ``None`` if the geometry of the relation members is not
361                  included in the query's result set, or if the relation is not deemed to represent
362                  an area.
363        geometry_details: More info on the validity of ``geometry``.
364
365    References:
366        - https://wiki.openstreetmap.org/wiki/Relation
367        - https://wiki.openstreetmap.org/wiki/Relation:multipolygon
368        - https://wiki.openstreetmap.org/wiki/Relation:boundary
369    """
370
371    bounds: Bbox | None
372    center: Point | None
373    members: list["Relationship"]
374    geometry: Polygon | MultiPolygon | None
375    geometry_details: GeometryDetails[Polygon | MultiPolygon] | None
376
377    def __iter__(self) -> Iterator[tuple[str | None, Element]]:
378        """Iterates over all members in the form of ``(role, element)``."""
379        for relship in self.members:
380            yield relship.role, relship.member
381
382
383@dataclass(kw_only=True, slots=True, repr=False)
384class Relationship(Spatial):
385    """
386    The relationship of an element that is part of a relation, with an optional role.
387
388    Attributes:
389        member:     any element
390        relation:   a relation that the member is a part of
391        role:       describes the function of the member in the context of the relation
392
393    References:
394        - https://wiki.openstreetmap.org/wiki/Relation#Roles
395    """
396
397    member: Element
398    relation: Relation
399    role: str | None
400
401    @property
402    def geojson(self) -> GeoJsonDict:
403        """
404        A mapping of ``member``.
405
406        This is ``member.geojson``, with the added properties ``role`` and ``__rel__``.
407        """
408        return _geojson_feature(self)
409
410    def __repr__(self) -> str:
411        role = f" as '{self.role}'" if self.role else " "
412        return f"{type(self).__name__}({self.member}{role} in {self.relation})"
413
414
415_KNOWN_ELEMENTS: Final[set[str]] = {"node", "way", "relation"}
416
417
418_ElementKey: TypeAlias = tuple[str, int]
419"""Elements are uniquely identified by the tuple (type, id)."""
420
421_MemberKey: TypeAlias = tuple[_ElementKey, str]
422"""Relation members are identified by their element key and role."""
423
424_OverpassDict: TypeAlias = dict[str, Any]
425"""A dictionary representing a JSON object returned by the Overpass API."""
426
427
428class _ElementCollector:
429    __slots__ = (
430        "member_dict",
431        "result_set",
432        "typed_dict",
433        "untyped_dict",
434    )
435
436    def __init__(self) -> None:
437        self.result_set: list[_ElementKey] = []
438        self.typed_dict: dict[_ElementKey, Element] = {}
439        self.untyped_dict: dict[_ElementKey, _OverpassDict] = defaultdict(dict)
440        self.member_dict: dict[int, list[_MemberKey]] = defaultdict(list)
441
442
443def collect_elements(query: Query) -> list[Element]:
444    """
445    Produce typed elements from the result set of a query.
446
447    This function exclusively collects elements that are of type "node", "way", or "relation".
448
449    Element data is "conflated", which means that if elements appear more than once in a
450    result set, their data is merged. This is useful f.e. when querying tags for relation members:
451    using ``rel(...); out tags;`` will only print tags for relation itself, not its members.
452    For those you will have to recurse down from the relation, which means members will show
453    up twice in the result set: once untagged as a member of the relation, and once tagged at
454    the top level. This function will have these two occurrences point to the same, single object.
455
456    The order of elements and relation members is retained.
457
458    If you want to query Overpass' "area" elements, you should use the `way(pivot)` or `rel(pivot)`
459    filter to select the element of the chosen type that defines the outline of the given area.
460    Otherwise, they will be ignored.
461
462    Derived elements with other types - f.e. produced by ``make`` and ``convert`` statements
463    or when using ``out count`` - are ignored. If you need to work with other derived elements,
464    you should not use this function.
465
466    Args:
467        query: a finished query
468
469    Returns:
470        the result set as a list of typed elements
471
472    Raises:
473        ValueError: If the input query is unfinished/has no result set.
474        KeyError: The only times there should be missing keys is when either using ``out noids``,
475                  or when building derived elements that are missing common keys.
476    """
477    if not query.done:
478        msg = "query has no result set"
479        raise ValueError(msg)
480
481    collector = _ElementCollector()
482    _collect_untyped(query, collector)
483    _collect_typed(collector)
484    _collect_relationships(collector)
485    return [collector.typed_dict[elem_key] for elem_key in collector.result_set]
486
487
488def _collect_untyped(query: Query, collector: _ElementCollector) -> None:
489    if query.result_set is None:
490        raise AssertionError
491
492    # Here we populate 'untyped_dict' with both top level elements, and
493    # relation members, while conflating their data if they appear as both.
494    # We also populate 'member_dict'.
495    for elem_dict in query.result_set:
496        if elem_dict.get("type") not in _KNOWN_ELEMENTS:
497            continue
498
499        key: _ElementKey = (elem_dict["type"], elem_dict["id"])
500
501        collector.result_set.append(key)
502        collector.untyped_dict[key].update(elem_dict)
503
504        if elem_dict["type"] != "relation":
505            continue
506
507        for mem in elem_dict["members"]:
508            key = (mem["type"], mem["ref"])
509            collector.untyped_dict[key].update(mem)
510            collector.member_dict[elem_dict["id"]].append((key, mem.get("role")))
511
512
513def _collect_typed(collector: _ElementCollector) -> None:
514    for elem_key, elem_dict in collector.untyped_dict.items():
515        elem_type, elem_id = elem_key
516
517        geometry = _geometry(elem_dict)
518
519        center = Point(elem_dict["center"].values()) if "center" in elem_dict else None
520        bounds = tuple(elem_dict["bounds"].values()) if "bounds" in elem_dict else None
521        assert bounds is None or len(bounds) == 4
522
523        meta: Metadata | None = None
524        if "timestamp" in elem_dict:
525            meta = Metadata(
526                timestamp=elem_dict["timestamp"],
527                version=elem_dict["version"],
528                changeset=elem_dict["changeset"],
529                user_name=elem_dict["user"],
530                user_id=elem_dict["uid"],
531            )
532
533        elem: Element
534
535        match elem_type:
536            case "node":
537                assert isinstance(geometry, Point | None)
538                assert geometry is None or geometry.is_valid
539                elem = Node(
540                    id=elem_id,
541                    tags=elem_dict.get("tags"),
542                    geometry=geometry,
543                    meta=meta,
544                    relations=[],  # add later
545                )
546            case "way":
547                assert isinstance(geometry, LineString | LinearRing | Polygon | None)
548                geometry_details = _try_validate_geometry(geometry) if geometry else None
549                elem = Way(
550                    id=elem_id,
551                    tags=elem_dict.get("tags"),
552                    geometry=geometry if geometry_details is None else geometry_details.best,
553                    geometry_details=geometry_details,
554                    meta=meta,
555                    node_ids=elem_dict.get("nodes"),
556                    bounds=bounds,
557                    center=center,
558                    relations=[],  # add later
559                )
560            case "relation":
561                assert isinstance(geometry, Polygon | MultiPolygon | None)
562                geometry_details = _try_validate_geometry(geometry) if geometry else None
563                elem = Relation(
564                    id=elem_id,
565                    tags=elem_dict.get("tags"),
566                    geometry=geometry if geometry_details is None else geometry_details.best,
567                    geometry_details=geometry_details,
568                    meta=meta,
569                    bounds=bounds,
570                    center=center,
571                    relations=[],  # add later
572                    members=[],  # add later
573                )
574            case _:
575                raise AssertionError
576
577        collector.typed_dict[elem_key] = elem
578
579
580def _collect_relationships(collector: _ElementCollector) -> None:
581    for rel_id, mem_roles in collector.member_dict.items():
582        rel = collector.typed_dict[("relation", rel_id)]
583        assert isinstance(rel, Relation)
584
585        for mem_key, mem_role in mem_roles:
586            mem = collector.typed_dict[mem_key]
587            relship = Relationship(member=mem, relation=rel, role=mem_role or None)
588            mem.relations.append(relship)
589            rel.members.append(relship)
590
591
592def _try_validate_geometry(geom: G) -> GeometryDetails[G]:
593    if geom.is_valid:
594        return GeometryDetails(valid=geom)
595
596    invalid_reason = shapely.is_valid_reason(geom)
597
598    if invalid_reason.startswith("Self-intersection") and isinstance(geom, MultiPolygon):
599        # we allow self-intersecting multi-polygons, if
600        # (1) the intersection is just lines or points, and
601        # (2) all the polygons inside are valid
602        intersection = shapely.intersection_all(geom.geoms)
603        accept = not isinstance(intersection, Polygon | MultiPolygon) and all(
604            poly.is_valid for poly in geom.geoms
605        )
606
607        if accept:
608            return GeometryDetails(accepted=geom, invalid_reason=invalid_reason)
609
610        return GeometryDetails(invalid=geom, invalid_reason=invalid_reason)
611
612    if isinstance(geom, Polygon):
613        valid_polygons = [g for g in _flatten(shapely.make_valid(geom)) if isinstance(g, Polygon)]
614        if len(valid_polygons) == 1:
615            return GeometryDetails(
616                fixed=valid_polygons[0],
617                invalid=geom,
618                invalid_reason=invalid_reason,
619            )
620
621    if isinstance(geom, MultiPolygon):
622        valid_polygons = [g for g in _flatten(shapely.make_valid(geom)) if isinstance(g, Polygon)]
623        if len(valid_polygons) == len(geom.geoms):
624            return GeometryDetails(
625                fixed=MultiPolygon(valid_polygons),
626                invalid=geom,
627                invalid_reason=invalid_reason,
628            )
629
630    return GeometryDetails(invalid=geom, invalid_reason=invalid_reason)
631
632
633def _geometry(raw_elem: _OverpassDict) -> BaseGeometry | None:
634    """
635    Construct the geometry a given OSM element makes up.
636
637    Args:
638        raw_elem: an element from a query's result set
639
640    Returns:
641        - None if there is no geometry available for this element.
642        - Point when given a node.
643        - LineString when given an open way.
644        - LinearRing when given a closed way, that supposedly is *not* an area.
645        - Polygon when given a closed way, that supposedly is an area.
646        - (Multi-)Polygons containing given a (multipolygon) relation.
647          Relation members that are not ways, or are not part of any polygon boundary, are
648          not part of the result geometry.
649
650    Raises:
651        ValueError: if element is not of type 'node', 'way', 'relation', or 'area'
652    """
653    if raw_elem.get("type") not in _KNOWN_ELEMENTS:
654        msg = "expected element of type 'node', 'way', 'relation', or 'area'"
655        raise ValueError(msg)
656
657    if raw_elem["type"] == "node":
658        lat, lon = raw_elem.get("lat"), raw_elem.get("lon")
659        if lat and lon:
660            return Point(lat, lon)
661
662    if raw_elem["type"] == "way":
663        ls = _line(raw_elem)
664        if ls and ls.is_ring and _is_area_element(raw_elem):
665            return Polygon(ls)
666        return ls
667
668    if _is_area_element(raw_elem):
669        outers = (
670            ls
671            for ls in (
672                _line(mem) for mem in raw_elem.get("members", ()) if mem.get("role") == "outer"
673            )
674            if ls
675        )
676        inners = (
677            ls
678            for ls in (
679                _line(mem) for mem in raw_elem.get("members", ()) if mem.get("role") == "inner"
680            )
681            if ls
682        )
683
684        shells = [ls for ls in _flatten(shapely.ops.linemerge(outers)) if ls.is_closed]
685        holes = [ls for ls in _flatten(shapely.ops.linemerge(inners)) if ls.is_closed]
686
687        polys = [
688            Polygon(shell=shell, holes=[hole for hole in holes if shell.contains(hole)])
689            for shell in shells
690        ]
691
692        if len(polys) == 1:
693            return polys[0]
694
695        return MultiPolygon(polys)
696
697    return None
698
699
700def _line(way: _OverpassDict) -> LineString | LinearRing | None:
701    """Returns the geometry of a way in the result set."""
702    if "geometry" not in way or len(way["geometry"]) < 2:
703        return None
704    is_ring = way["geometry"][0] == way["geometry"][-1]
705    cls = LinearRing if is_ring else LineString
706    return cls((c["lat"], c["lon"]) for c in way["geometry"])
707
708
709def _flatten(obj: BaseGeometry) -> Iterable[BaseGeometry]:
710    """Recursively flattens multipart geometries."""
711    if isinstance(obj, BaseMultipartGeometry):
712        return (nested for contained in obj.geoms for nested in _flatten(contained))
713    return (obj,)
714
715
716def _is_area_element(el: _OverpassDict) -> bool:
717    """
718    Decide if ``el`` likely represents an area, and should be viewed as a (multi-)polygon.
719
720    Args:
721        el: a way or relation from a query's result set
722
723    Returns:
724        ``False`` if the input is not a relation or closed way.
725        ``False``, unless there are tags which indicate that the way represents an area.
726
727    References:
728        - https://wiki.openstreetmap.org/wiki/Overpass_API/Areas
729        - https://github.com/drolbr/Overpass-API/blob/master/src/rules/areas.osm3s
730          (from 2018-04-09)
731        - https://wiki.openstreetmap.org/wiki/Overpass_turbo/Polygon_Features
732        - https://github.com/tyrasd/osm-polygon-features/blob/master/polygon-features.json
733          (from 2016-11-03)
734    """
735    # Check if the element is explicitly an area
736    if el["type"] == "area":
737        return True
738
739    # Check if a given way is open
740    if el["type"] == "way" and ("geometry" not in el or el["geometry"][0] != el["geometry"][-1]):
741        return False
742
743    # Assume not an area if there are no tags available
744    if "tags" not in el:
745        return False
746
747    tags = el["tags"]
748
749    # Check if the element is explicitly tagged as not area
750    if tags.get("area") == "no":
751        return False
752
753    # Check if there is a tag where any value other than 'no' suggests area
754    # (note: Overpass may require the "name" tag as well)
755    if any(tags.get(name, "no") != "no" for name in _AREA_TAG_NAMES):
756        return True
757
758    # Check if there are tag values that suggest area
759    # (note: Overpass may require the "name" tag as well)
760    return any(
761        (
762            v in _AREA_TAG_VALUES_ONE_OF.get(k, ())
763            or v not in _AREA_TAG_VALUES_NONE_OF.get(k, (v,))
764            for k, v in tags.items()
765        )
766    )
767
768
769_AREA_TAG_NAMES: Final[set[str]] = {
770    "area",
771    "area:highway",
772    "amenity",
773    "boundary",
774    "building",
775    "building:part",
776    "craft",
777    "golf",
778    "historic",
779    "indoor",
780    "landuse",
781    "leisure",
782    "military",
783    "office",
784    "place",
785    "public_transport",
786    "ruins",
787    "shop",
788    "tourism",
789    # for relations
790    "admin_level",
791    "postal_code",
792    "addr:postcode",
793}
794
795_AREA_TAG_VALUES_ONE_OF: Final[dict[str, set[str]]] = {
796    "barrier": {"city_wall", "ditch", "hedge", "retaining_wall", "wall", "spikes"},
797    "highway": {"services", "rest_area", "escape", "elevator"},
798    "power": {"plant", "substation", "generator", "transformer"},
799    "railway": {"station", "turntable", "roundhouse", "platform"},
800    "waterway": {"riverbank", "dock", "boatyard", "dam"},
801    # for relations
802    "type": {"multipolygon"},
803}
804
805_AREA_TAG_VALUES_NONE_OF: Final[dict[str, set[str]]] = {
806    "aeroway": {"no", "taxiway"},
807    "man_made": {"no", "cutline", "embankment", "pipeline"},
808    "natural": {"no", "coastline", "cliff", "ridge", "arete", "tree_row"},
809}
810
811_WIKIDATA_Q_ID: Final[re.Pattern[str]] = re.compile(r"^Q\d+$")
812
813
814def _geojson_properties(obj: Element | Relationship) -> GeoJsonDict:
815    elem = obj if isinstance(obj, Element) else obj.member
816
817    properties = {
818        "id": elem.id,
819        "type": elem.type,
820        "tags": elem.tags,
821        "timestamp": elem.meta.timestamp if elem.meta else None,
822        "version": elem.meta.version if elem.meta else None,
823        "changeset": elem.meta.changeset if elem.meta else None,
824        "user": elem.meta.user_name if elem.meta else None,
825        "uid": elem.meta.user_id if elem.meta else None,
826        "nodes": getattr(elem, "nodes", None),
827    }
828
829    if isinstance(elem, Way | Relation):
830        # TODO: these are lat/lon order, as opposed to the geometry in GeoJSON - problem?
831        properties["bounds"] = elem.bounds
832        properties["center"] = elem.center.coords[0] if elem.center else None
833
834    properties = {k: v for k, v in properties.items() if v is not None}
835
836    if isinstance(obj, Relationship):
837        properties["role"] = obj.role or ""
838        properties["__rel__"] = _geojson_properties(obj.relation)
839
840    return properties
841
842
843def _geojson_geometry(obj: Element | Relationship) -> GeoJsonDict | None:
844    elem = obj if isinstance(obj, Element) else obj.member
845
846    if not elem.base_geometry:
847        return None
848
849    # Flip coordinates for GeoJSON compliance.
850    geom = shapely.ops.transform(lambda lat, lon: (lon, lat), elem.base_geometry)
851
852    # GeoJSON-like mapping that implements __geo_interface__.
853    mapping = shapely.geometry.mapping(geom)
854
855    # This geometry does not exist in GeoJSON.
856    if mapping["type"] == "LinearRing":
857        mapping["type"] = "LineString"
858
859    return mapping
860
861
862def _geojson_bbox(obj: Element | Relationship) -> Bbox | None:
863    elem = obj if isinstance(obj, Element) else obj.member
864
865    geom = elem.base_geometry
866    if not geom:
867        return elem.bounds if isinstance(elem, Way | Relation) else None
868
869    bounds = geom.bounds  # can be (nan, nan, nan, nan)
870    if not any(math.isnan(c) for c in bounds):
871        (minlat, minlon, maxlat, maxlon) = bounds
872        return minlon, minlat, maxlon, maxlat
873
874    return None
875
876
877def _geojson_feature(obj: Element | Relationship) -> GeoJsonDict:
878    feature: GeoJsonDict = {
879        "type": "Feature",
880        "geometry": _geojson_geometry(obj),
881        "properties": _geojson_properties(obj),
882    }
883
884    bbox = _geojson_bbox(obj)
885    if bbox:
886        feature["bbox"] = bbox
887
888    return feature