cartography 0.11.1

Cartography is a map rendering library for Geographic features expressed using [georust](https://georust.org/) libraries.
Documentation
//! OSM PBF file reader

use std::collections::HashMap;
use std::path::Path;

use geo::Contains as _;
use osmpbf::{Element, ElementReader};

use crate::{FeaturesVecLayer, Result};

use super::OsmFeature;
use super::stitching;

fn assign_inner_rings(
  outer_rings: &[geo::LineString],
  inner_rings: &[geo::LineString],
) -> Vec<(geo::LineString, Vec<geo::LineString>)>
{
  let mut assigned = vec![Vec::new(); outer_rings.len()];

  for inner in inner_rings
  {
    let Some(inner_anchor) = inner.points().next()
    else
    {
      continue;
    };

    if let Some((idx, _)) = outer_rings
      .iter()
      .enumerate()
      .find(|(_, outer)| geo::Polygon::new((*outer).clone(), vec![]).contains(&inner_anchor))
    {
      assigned[idx].push(inner.clone());
    }
  }

  outer_rings.iter().cloned().zip(assigned).collect()
}

/// Internal representation for collecting ways by ID during passes
#[derive(Clone)]
struct WayGeometry
{
  geometry: geo::LineString,
}

/// Parse an OSM PBF file and return a layer with all nodes, ways, and multipolygon relations.
pub fn read_pbf(path: impl AsRef<Path>) -> Result<FeaturesVecLayer<OsmFeature>>
{
  let path = path.as_ref();

  let reader = ElementReader::from_path(path)?;

  let node_coords = reader.par_map_reduce(
    |element| match element
    {
      Element::Node(node) =>
      {
        let coord = geo::coord! {
          x: node.lon(),
          y: node.lat(),
        };
        vec![(node.id(), coord)]
      }
      Element::DenseNode(node) =>
      {
        let coord = geo::coord! {
          x: node.lon(),
          y: node.lat(),
        };
        vec![(node.id(), coord)]
      }
      _ => vec![],
    },
    Vec::new,
    |mut a, b| {
      a.extend(b);
      a
    },
  )?;

  let node_coords: HashMap<i64, geo::Coord> = node_coords.into_iter().collect();

  let reader = ElementReader::from_path(path)?;
  let mut features = Vec::new();
  let mut ways_by_id: HashMap<i64, WayGeometry> = HashMap::new();

  reader.for_each(|element| match element
  {
    Element::Node(node) =>
    {
      let tags = node
        .tags()
        .map(|(k, v)| (k.to_owned(), v.to_owned()))
        .collect();
      let geometry = geo::Geometry::Point(geo::Point::new(node.lon(), node.lat()));
      features.push(OsmFeature {
        id: node.id(),
        geometry,
        tags,
      });
    }
    Element::Way(way) =>
    {
      let tags: Vec<(String, String)> = way
        .tags()
        .map(|(k, v)| (k.to_owned(), v.to_owned()))
        .collect();

      let node_refs: Vec<i64> = way.refs().collect();
      if node_refs.len() < 2
      {
        return;
      }

      let coords: Vec<geo::Coord> = node_refs
        .iter()
        .filter_map(|&node_id| node_coords.get(&node_id).copied())
        .collect();

      if coords.len() != node_refs.len()
      {
        return;
      }

      let linestring = geo::LineString::new(coords);

      ways_by_id.insert(
        way.id(),
        WayGeometry {
          geometry: linestring.clone(),
        },
      );

      let geometry = if node_refs.first() == node_refs.last() && linestring.0.len() >= 3
      {
        geo::Geometry::Polygon(geo::Polygon::new(linestring, vec![]))
      }
      else if let Some(area_tag) = tags.iter().find(|(k, _)| k == "area")
      {
        if area_tag.1 == "yes" && linestring.0.len() >= 3
        {
          geo::Geometry::Polygon(geo::Polygon::new(linestring, vec![]))
        }
        else
        {
          geo::Geometry::LineString(linestring)
        }
      }
      else
      {
        geo::Geometry::LineString(linestring)
      };

      features.push(OsmFeature {
        id: way.id(),
        geometry,
        tags,
      });
    }
    _ =>
    {}
  })?;

  let reader = ElementReader::from_path(path)?;
  reader.for_each(|element| {
    if let Element::Relation(relation) = element
    {
      let is_multipolygon = relation
        .tags()
        .any(|(k, v)| k == "type" && v == "multipolygon");
      if !is_multipolygon
      {
        return;
      }

      let tags: Vec<(String, String)> = relation
        .tags()
        .map(|(k, v)| (k.to_owned(), v.to_owned()))
        .collect();

      let mut outer_ways = Vec::new();
      let mut inner_ways = Vec::new();

      for member in relation.members()
      {
        let member_id = member.member_id;
        let role = member.role().unwrap_or("");

        if member.member_type != osmpbf::RelMemberType::Way
        {
          continue;
        }

        if let Some(way_geom) = ways_by_id.get(&member_id)
        {
          if role == "outer"
          {
            outer_ways.push(way_geom.geometry.clone());
          }
          else if role == "inner"
          {
            inner_ways.push(way_geom.geometry.clone());
          }
        }
      }

      if outer_ways.is_empty()
      {
        return;
      }

      let outer_rings = match stitching::stitch_rings(outer_ways)
      {
        Ok(rings) => rings,
        Err(_) => return,
      };

      let inner_rings = if !inner_ways.is_empty()
      {
        stitching::stitch_rings(inner_ways).unwrap_or_default()
      }
      else
      {
        vec![]
      };

      let geometry = if outer_rings.len() == 1
      {
        let exterior = outer_rings[0].clone();
        let interiors = assign_inner_rings(&outer_rings, &inner_rings)
          .into_iter()
          .next()
          .map(|(_, holes)| holes)
          .unwrap_or_default();
        let polygon = geo::Polygon::new(exterior, interiors);
        geo::Geometry::Polygon(polygon)
      }
      else
      {
        let polygons: Vec<geo::Polygon> = assign_inner_rings(&outer_rings, &inner_rings)
          .into_iter()
          .map(|(outer, interiors)| geo::Polygon::new(outer, interiors))
          .collect();
        geo::Geometry::MultiPolygon(geo::MultiPolygon::new(polygons))
      };

      features.push(OsmFeature {
        id: relation.id(),
        geometry,
        tags,
      });
    }
  })?;

  Ok(FeaturesVecLayer::from(features))
}

#[cfg(test)]
mod tests
{
  use super::*;

  #[test]
  fn assign_inner_rings_matches_containing_outer()
  {
    let outer_a = geo::LineString::from(vec![
      (0.0, 0.0),
      (4.0, 0.0),
      (4.0, 4.0),
      (0.0, 4.0),
      (0.0, 0.0),
    ]);
    let outer_b = geo::LineString::from(vec![
      (10.0, 10.0),
      (14.0, 10.0),
      (14.0, 14.0),
      (10.0, 14.0),
      (10.0, 10.0),
    ]);
    let inner = geo::LineString::from(vec![
      (1.0, 1.0),
      (2.0, 1.0),
      (2.0, 2.0),
      (1.0, 2.0),
      (1.0, 1.0),
    ]);

    let assigned = assign_inner_rings(
      &[outer_a.clone(), outer_b.clone()],
      std::slice::from_ref(&inner),
    );
    assert_eq!(assigned.len(), 2);
    assert_eq!(assigned[0].1, vec![inner]);
    assert!(assigned[1].1.is_empty());
  }
}