diff --git a/Cargo.lock b/Cargo.lock index a0dcb4c4d5..cd0ca2d5be 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -1902,9 +1902,9 @@ dependencies = [ [[package]] name = "geom" version = "0.1.0" +source = "git+https://github.com/a-b-street/geom#34194e73db3fd6343ad223e8a4b1075201094f1a" dependencies = [ "anyhow", - "bincode", "earcutr", "fs-err", "geo", @@ -1913,8 +1913,6 @@ dependencies = [ "instant", "ordered-float", "polylabel", - "rand", - "rand_xorshift", "rstar", "serde", "serde_json", @@ -3362,10 +3360,22 @@ version = "6.4.1" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "9b7820b9daea5457c9f21c69448905d723fbd21136ccf521748f23fd49e723ee" +[[package]] +name = "osm2lanes" +version = "0.1.0" +source = "git+https://github.com/a-b-street/osm2streets#086a347e6cad7c556d7f7d5bf57d4b270c8a030f" +dependencies = [ + "abstutil", + "anyhow", + "enumset", + "geom", + "serde", +] + [[package]] name = "osm2streets" version = "0.1.0" -source = "git+https://github.com/a-b-street/osm2streets#c098792f790ce9cefce074f9e4a9c8b1f27900bc" +source = "git+https://github.com/a-b-street/osm2streets#086a347e6cad7c556d7f7d5bf57d4b270c8a030f" dependencies = [ "abstutil", "anyhow", @@ -3373,8 +3383,9 @@ dependencies = [ "geo", "geojson", "geom", - "itertools 0.10.5", + "itertools 0.11.0", "log", + "osm2lanes", "petgraph", "serde", "serde_json", @@ -4555,7 +4566,7 @@ checksum = "a8f112729512f8e442d81f95a8a7ddf2b7c6b8a1a6f509a95864142b30cab2d3" [[package]] name = "streets_reader" version = "0.1.0" -source = "git+https://github.com/a-b-street/osm2streets#c098792f790ce9cefce074f9e4a9c8b1f27900bc" +source = "git+https://github.com/a-b-street/osm2streets#086a347e6cad7c556d7f7d5bf57d4b270c8a030f" dependencies = [ "abstutil", "anyhow", diff --git a/Cargo.toml b/Cargo.toml index c30201638e..6180c40d2c 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -10,7 +10,6 @@ members = [ "cli", "collisions", "convert_osm", - "geom", "headless", "importer", "kml", @@ -49,6 +48,7 @@ futures = { version = "0.3.27"} futures-channel = { version = "0.3.29"} geo = "0.26.0" geojson = { version = "0.24.1", features = ["geo-types"] } +geom = { git = "https://github.com/a-b-street/geom" } getrandom = "0.2.11" instant = "0.1.7" log = "0.4.20" @@ -65,7 +65,6 @@ web-sys = "0.3.65" # due to the 2 core dependency crates listed below. This patch is required to # avoid Cargo from getting confused. [patch."https://github.com/a-b-street/abstreet/"] -geom = { path = "geom" } abstutil = { path = "abstutil" } [patch.crates-io] diff --git a/apps/fifteen_min/Cargo.toml b/apps/fifteen_min/Cargo.toml index f4db617f83..2a8bb1f5f9 100644 --- a/apps/fifteen_min/Cargo.toml +++ b/apps/fifteen_min/Cargo.toml @@ -16,7 +16,7 @@ abstio = { path = "../../abstio" } abstutil = { path = "../../abstutil" } contour = { workspace = true } geojson = { workspace = true } -geom = { path = "../../geom" } +geom = { workspace = true } getrandom = { workspace = true, optional = true } log = { workspace = true } map_gui = { path = "../../map_gui" } diff --git a/apps/game/Cargo.toml b/apps/game/Cargo.toml index a60d3cb021..779980d510 100644 --- a/apps/game/Cargo.toml +++ b/apps/game/Cargo.toml @@ -29,7 +29,7 @@ fs-err = { workspace = true } futures-channel = { workspace = true } geo = { workspace = true } geojson = { workspace = true } -geom = { path = "../../geom" } +geom = { workspace = true } getrandom = { workspace = true, optional = true } instant = { workspace = true } kml = { path = "../../kml" } diff --git a/apps/ltn/Cargo.toml b/apps/ltn/Cargo.toml index 7cd02e96d3..52b17e7bef 100644 --- a/apps/ltn/Cargo.toml +++ b/apps/ltn/Cargo.toml @@ -22,7 +22,7 @@ flate2 = { workspace = true } futures-channel = { workspace = true } geo = { workspace = true } geojson = { workspace = true } -geom = { path = "../../geom" } +geom = { workspace = true } getrandom = { workspace = true, optional = true } home = { version = "0.5.5", optional = true } lazy_static = "1.4.0" diff --git a/apps/map_editor/Cargo.toml b/apps/map_editor/Cargo.toml index 33f911d06c..26423cc74e 100644 --- a/apps/map_editor/Cargo.toml +++ b/apps/map_editor/Cargo.toml @@ -15,7 +15,7 @@ wasm = ["getrandom/js", "wasm-bindgen", "widgetry/wasm-backend"] abstio = { path = "../../abstio" } abstutil = { path = "../../abstutil" } fs-err = { workspace = true } -geom = { path = "../../geom" } +geom = { workspace = true } getrandom = { workspace = true, optional = true } log = { workspace = true } raw_map = { path = "../../raw_map" } diff --git a/apps/osm_viewer/Cargo.toml b/apps/osm_viewer/Cargo.toml index d20892ea8e..27bf461cfc 100644 --- a/apps/osm_viewer/Cargo.toml +++ b/apps/osm_viewer/Cargo.toml @@ -14,7 +14,7 @@ wasm = ["getrandom/js", "map_gui/wasm", "wasm-bindgen", "widgetry/wasm-backend"] [dependencies] abstio = { path = "../../abstio" } abstutil = { path = "../../abstutil" } -geom = { path = "../../geom" } +geom = { workspace = true } getrandom = { workspace = true, optional = true } map_gui = { path = "../../map_gui" } map_model = { path = "../../map_model" } diff --git a/apps/parking_mapper/Cargo.toml b/apps/parking_mapper/Cargo.toml index e743179bcc..9fa2a689c5 100644 --- a/apps/parking_mapper/Cargo.toml +++ b/apps/parking_mapper/Cargo.toml @@ -12,7 +12,7 @@ abstio = { path = "../../abstio" } abstutil = { path = "../../abstutil" } anyhow = { workspace = true } fs-err = { workspace = true } -geom = { path = "../../geom" } +geom = { workspace = true } log = { workspace = true } map_gui = { path = "../../map_gui" } map_model = { path = "../../map_model" } diff --git a/apps/santa/Cargo.toml b/apps/santa/Cargo.toml index d5bbcfe7e9..206a333b3e 100644 --- a/apps/santa/Cargo.toml +++ b/apps/santa/Cargo.toml @@ -15,7 +15,7 @@ wasm = ["getrandom/js", "map_gui/wasm", "wasm-bindgen", "widgetry/wasm-backend"] abstio = { path = "../../abstio" } abstutil = { path = "../../abstutil" } anyhow = { workspace = true } -geom = { path = "../../geom" } +geom = { workspace = true } getrandom = { workspace = true, optional = true } kml = { path = "../../kml" } log = { workspace = true } diff --git a/blockfinding/Cargo.toml b/blockfinding/Cargo.toml index 08a6e5903e..7834a39bd4 100644 --- a/blockfinding/Cargo.toml +++ b/blockfinding/Cargo.toml @@ -6,7 +6,7 @@ edition = "2021" [dependencies] abstutil = { path = "../abstutil" } anyhow = { workspace = true } -geom = { path = "../geom" } +geom = { workspace = true } log = { workspace = true } map_model = { path = "../map_model" } serde = { workspace = true } diff --git a/cli/Cargo.toml b/cli/Cargo.toml index 4b5dc8d30a..eb4cc9d031 100644 --- a/cli/Cargo.toml +++ b/cli/Cargo.toml @@ -11,7 +11,7 @@ convert_osm = { path = "../convert_osm" } csv = { workspace = true } fs-err = { workspace = true } geo = { workspace = true } -geom = { path = "../geom" } +geom = { workspace = true } importer = { path = "../importer" } log = { workspace = true } map_model = { path = "../map_model" } diff --git a/collisions/Cargo.toml b/collisions/Cargo.toml index 593f1d3cba..6e1beb26a5 100644 --- a/collisions/Cargo.toml +++ b/collisions/Cargo.toml @@ -5,7 +5,7 @@ authors = ["Dustin Carlino "] edition = "2021" [dependencies] -geom = { path = "../geom" } +geom = { workspace = true } kml = { path = "../kml" } log = { workspace = true } serde = { workspace = true } diff --git a/convert_osm/Cargo.toml b/convert_osm/Cargo.toml index 6185c18417..78f6983c99 100644 --- a/convert_osm/Cargo.toml +++ b/convert_osm/Cargo.toml @@ -10,7 +10,7 @@ abstutil = { path = "../abstutil" } anyhow = { workspace = true } csv = { workspace = true } fs-err = { workspace = true } -geom = { path = "../geom" } +geom = { workspace = true } kml = { path = "../kml" } log = { workspace = true } osm2streets = { git = "https://github.com/a-b-street/osm2streets" } diff --git a/geom/Cargo.toml b/geom/Cargo.toml deleted file mode 100644 index 222148656f..0000000000 --- a/geom/Cargo.toml +++ /dev/null @@ -1,24 +0,0 @@ -[package] -name = "geom" -version = "0.1.0" -authors = ["Dustin Carlino "] -edition = "2021" - -[dependencies] -anyhow = { workspace = true } -earcutr = "0.4.3" -fs-err = { workspace = true } -geo = { workspace = true } -geojson = { workspace = true } -histogram = "0.6.9" -instant = { workspace = true } -ordered-float = { version = "3.7.0", features=["serde"] } -polylabel = "2.5.0" -rstar = "0.11.0" -serde = { workspace = true } -serde_json = { workspace = true } - -[dev-dependencies] -bincode = { workspace = true } -rand = { workspace = true } -rand_xorshift = { workspace = true } diff --git a/geom/README.md b/geom/README.md deleted file mode 100644 index be0f4226aa..0000000000 --- a/geom/README.md +++ /dev/null @@ -1,17 +0,0 @@ -# geom - -This crate contains primitive types used by A/B Street. It's unclear if other -apps will have any use for this crate. In some cases, `geom` just wraps much -more polished APIs, like `rust-geo`. In others, it has its own geometric -algorithms, but they likely have many bugs and make use-case-driven assumptions. -So, be warned if you use this. - -## Contents - -Many of the types are geometric: `Pt2D`, `Ring`, `Distance`, `Line`, -`InfiniteLine`, `FindClosest`, `Circle`, `Angle`, `LonLat`, `Bounds`, -`GPSBounds`, `PolyLine`, `Polygon`, `Triangle`. - -Some involve time: `Time`, `Duration`, `Speed`. - -And there's also a `Percent` wrapper and a `Histogram`. diff --git a/geom/src/angle.rs b/geom/src/angle.rs deleted file mode 100644 index dddcc80aca..0000000000 --- a/geom/src/angle.rs +++ /dev/null @@ -1,137 +0,0 @@ -use std::fmt; - -use serde::{Deserialize, Serialize}; - -/// An angle, stored in radians. -#[derive(Clone, Copy, Debug, Serialize, Deserialize, PartialEq, PartialOrd)] -pub struct Angle(f64); - -impl Angle { - pub const ZERO: Angle = Angle(0.0); - - /// Create an angle in radians. - // TODO Normalize here, and be careful about % vs euclid_rem - pub fn new_rads(rads: f64) -> Angle { - // Retain more precision for angles... - Angle((rads * 10_000_000.0).round() / 10_000_000.0) - } - - /// Create an angle in degrees. - pub fn degrees(degs: f64) -> Angle { - Angle::new_rads(degs.to_radians()) - } - - /// Invert the direction of this angle. - pub fn opposite(self) -> Angle { - Angle::new_rads(self.0 + std::f64::consts::PI) - } - - pub(crate) fn invert_y(self) -> Angle { - Angle::new_rads(2.0 * std::f64::consts::PI - self.0) - } - - /// Rotates this angle by some degrees. - pub fn rotate_degs(self, degrees: f64) -> Angle { - Angle::new_rads(self.0 + degrees.to_radians()) - } - - /// Returns [0, 2pi) - pub fn normalized_radians(self) -> f64 { - if self.0 < 0.0 { - // TODO Be more careful about how we store the angle, to make sure this works - self.0 + (2.0 * std::f64::consts::PI) - } else { - self.0 - } - } - - /// Returns [0, 360) - pub fn normalized_degrees(self) -> f64 { - self.normalized_radians().to_degrees() - } - - /// Returns [-180, 180] - pub fn simple_shortest_rotation_towards(self, other: Angle) -> f64 { - // https://math.stackexchange.com/questions/110080/shortest-way-to-achieve-target-angle - ((self.normalized_degrees() - other.normalized_degrees() + 540.0) % 360.0) - 180.0 - } - - /// Logically this returns [-180, 180], but keep in mind when we print this angle, it'll - /// normalize to be [0, 360]. - pub fn shortest_rotation_towards(self, other: Angle) -> Angle { - Angle::degrees(self.simple_shortest_rotation_towards(other)) - } - - /// True if this angle is within some degrees of another, accounting for rotation - pub fn approx_eq(self, other: Angle, within_degrees: f64) -> bool { - self.simple_shortest_rotation_towards(other).abs() < within_degrees - } - - /// True if this angle is within some degrees of another, accounting for rotation and allowing - /// two angles that point in opposite directions - pub fn approx_parallel(self, other: Angle, within_degrees: f64) -> bool { - self.approx_eq(other, within_degrees) || self.opposite().approx_eq(other, within_degrees) - } - - /// I don't know how to describe what this does. Use for rotating labels in map-space and making - /// sure the text is never upside-down. - pub fn reorient(self) -> Angle { - let theta = self.normalized_degrees().rem_euclid(360.0); - let mut result = self; - if theta > 90.0 { - result = result.opposite(); - } - if theta > 270.0 { - result = result.opposite(); - } - result - } - - /// Calculates the average of some angles. - pub fn average(input: Vec) -> Angle { - // Thanks https://rosettacode.org/wiki/Averages/Mean_angle - let num = input.len() as f64; - let mut cos_sum = 0.0; - let mut sin_sum = 0.0; - for x in input { - cos_sum += x.0.cos(); - sin_sum += x.0.sin(); - } - Angle::new_rads((sin_sum / num).atan2(cos_sum / num)) - } -} - -impl fmt::Display for Angle { - fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { - write!(f, "Angle({} degrees)", self.normalized_degrees()) - } -} - -impl std::ops::Add for Angle { - type Output = Angle; - - fn add(self, other: Angle) -> Angle { - Angle::new_rads(self.0 + other.0) - } -} - -impl std::ops::Neg for Angle { - type Output = Angle; - - fn neg(self) -> Angle { - Angle::new_rads(-self.0) - } -} - -#[cfg(test)] -mod tests { - use super::*; - - #[test] - fn test_average() { - assert_eq!( - Angle::degrees(30.0), - Angle::average(vec![Angle::degrees(30.0)]) - ); - } -} diff --git a/geom/src/bounds.rs b/geom/src/bounds.rs deleted file mode 100644 index 56afd4deef..0000000000 --- a/geom/src/bounds.rs +++ /dev/null @@ -1,333 +0,0 @@ -use serde::{Deserialize, Serialize}; - -use rstar::primitives::{GeomWithData, Rectangle}; -use rstar::{RTree, SelectionFunction, AABB}; - -use crate::{Circle, Distance, LonLat, Polygon, Pt2D, Ring}; - -/// Represents a rectangular boundary of `Pt2D` points. -#[derive(Clone, Copy, Debug, PartialEq, Serialize, Deserialize)] -pub struct Bounds { - pub min_x: f64, - pub min_y: f64, - pub max_x: f64, - pub max_y: f64, -} - -impl Bounds { - /// A boundary including no points. - pub fn new() -> Bounds { - Bounds { - min_x: f64::MAX, - min_y: f64::MAX, - max_x: f64::MIN, - max_y: f64::MIN, - } - } - - pub fn zero() -> Self { - Bounds { - min_x: 0.0, - min_y: 0.0, - max_x: 0.0, - max_y: 0.0, - } - } - - /// Create a boundary covering some points. - pub fn from(pts: &[Pt2D]) -> Bounds { - let mut b = Bounds::new(); - for pt in pts { - b.update(*pt); - } - b - } - - /// Create a boundary covering some polygons. - pub fn from_polygons(polygons: &[Polygon]) -> Bounds { - let mut b = Bounds::new(); - for poly in polygons { - for ring in poly.get_rings() { - for pt in ring.points() { - b.update(*pt); - } - } - } - b - } - - /// Update the boundary to include this point. - pub fn update(&mut self, pt: Pt2D) { - self.min_x = self.min_x.min(pt.x()); - self.max_x = self.max_x.max(pt.x()); - self.min_y = self.min_y.min(pt.y()); - self.max_y = self.max_y.max(pt.y()); - } - - /// Unions two boundaries. - pub fn union(&mut self, other: Bounds) { - self.update(Pt2D::new(other.min_x, other.min_y)); - self.update(Pt2D::new(other.max_x, other.max_y)); - } - - /// Expand the existing boundary by some distance evenly on all sides. - pub fn add_buffer(&mut self, sides: Distance) { - self.min_x -= sides.inner_meters(); - self.max_x += sides.inner_meters(); - self.min_y -= sides.inner_meters(); - self.max_y += sides.inner_meters(); - } - - /// Transform the boundary by scaling its corners. - pub fn scale(mut self, factor: f64) -> Self { - self.min_x *= factor; - self.min_y *= factor; - self.max_x *= factor; - self.max_y *= factor; - self - } - - /// True if the point is within the boundary. - pub fn contains(&self, pt: Pt2D) -> bool { - pt.x() >= self.min_x && pt.x() <= self.max_x && pt.y() >= self.min_y && pt.y() <= self.max_y - } - - /// Converts the boundary to the format used by `rstar`. - pub fn as_bbox(&self, data: T) -> GeomWithData, T> { - GeomWithData::new( - Rectangle::from_corners([self.min_x, self.min_y], [self.max_x, self.max_y]), - data, - ) - } - - /// Creates a rectangle covering this boundary. - pub fn get_rectangle(&self) -> Polygon { - Ring::must_new(vec![ - Pt2D::new(self.min_x, self.min_y), - Pt2D::new(self.max_x, self.min_y), - Pt2D::new(self.max_x, self.max_y), - Pt2D::new(self.min_x, self.max_y), - Pt2D::new(self.min_x, self.min_y), - ]) - .into_polygon() - } - - /// Creates a circle centered in the middle of this boundary. Always uses the half width as a - /// radius, so if the width and height don't match, this is pretty meaningless. - pub fn to_circle(&self) -> Circle { - Circle::new(self.center(), Distance::meters(self.width() / 2.0)) - } - - /// The width of this boundary. - // TODO Really should be Distance - pub fn width(&self) -> f64 { - self.max_x - self.min_x - } - - /// The height of this boundary. - pub fn height(&self) -> f64 { - self.max_y - self.min_y - } - - /// The center point of this boundary. - pub fn center(&self) -> Pt2D { - Pt2D::new( - self.min_x + self.width() / 2.0, - self.min_y + self.height() / 2.0, - ) - } -} - -/// Represents a rectangular boundary of `LonLat` points. After building one of these, `LonLat`s -/// can be transformed into `Pt2D`s, treating the top-left of the boundary as (0, 0), and growing -/// to the right and down (screen-drawing order, not Cartesian) in meters. -#[derive(Clone, Debug, PartialEq, Serialize, Deserialize)] -pub struct GPSBounds { - pub min_lon: f64, - pub min_lat: f64, - pub max_lon: f64, - pub max_lat: f64, -} - -impl GPSBounds { - /// A boundary including no points. - pub fn new() -> GPSBounds { - GPSBounds { - min_lon: f64::MAX, - min_lat: f64::MAX, - max_lon: f64::MIN, - max_lat: f64::MIN, - } - } - - /// Create a boundary covering some points. - pub fn from(pts: Vec) -> GPSBounds { - let mut b = GPSBounds::new(); - for pt in pts { - b.update(pt); - } - b - } - - /// Update the boundary to include this point. - pub fn update(&mut self, pt: LonLat) { - self.min_lon = self.min_lon.min(pt.x()); - self.max_lon = self.max_lon.max(pt.x()); - self.min_lat = self.min_lat.min(pt.y()); - self.max_lat = self.max_lat.max(pt.y()); - } - - /// True if the point is within the boundary. - pub fn contains(&self, pt: LonLat) -> bool { - pt.x() >= self.min_lon - && pt.x() <= self.max_lon - && pt.y() >= self.min_lat - && pt.y() <= self.max_lat - } - - /// The bottom-right corner of the boundary, in map-space. - // TODO cache this - pub fn get_max_world_pt(&self) -> Pt2D { - let width = LonLat::new(self.min_lon, self.min_lat) - .gps_dist(LonLat::new(self.max_lon, self.min_lat)); - let height = LonLat::new(self.min_lon, self.min_lat) - .gps_dist(LonLat::new(self.min_lon, self.max_lat)); - Pt2D::new(width.inner_meters(), height.inner_meters()) - } - - /// Converts the boundary to map-space. - pub fn to_bounds(&self) -> Bounds { - let mut b = Bounds::new(); - b.update(Pt2D::new(0.0, 0.0)); - b.update(self.get_max_world_pt()); - b - } - - /// Convert all points to map-space, failing if any points are outside this boundary. - pub fn try_convert(&self, pts: &[LonLat]) -> Option> { - let mut result = Vec::new(); - for gps in pts { - if !self.contains(*gps) { - return None; - } - result.push(gps.to_pt(self)); - } - Some(result) - } - - /// Convert all points to map-space. The points may be outside this boundary. - pub fn convert(&self, pts: &[LonLat]) -> Vec { - pts.iter().map(|gps| gps.to_pt(self)).collect() - } - - /// Convert map-space points back to `LonLat`s. This is only valid if the `GPSBounds` used - /// is the same as the one used to originally produce the `Pt2D`s. - pub fn convert_back(&self, pts: &[Pt2D]) -> Vec { - pts.iter().map(|pt| pt.to_gps(self)).collect() - } - - pub fn convert_back_xy(&self, x: f64, y: f64) -> LonLat { - let (width, height) = { - let pt = self.get_max_world_pt(); - (pt.x(), pt.y()) - }; - - let lon = (x / width * (self.max_lon - self.min_lon)) + self.min_lon; - let lat = self.min_lat + ((self.max_lat - self.min_lat) * (height - y) / height); - LonLat::new(lon, lat) - } - - /// Returns points in order covering this boundary. - pub fn get_rectangle(&self) -> Vec { - vec![ - LonLat::new(self.min_lon, self.min_lat), - LonLat::new(self.max_lon, self.min_lat), - LonLat::new(self.max_lon, self.max_lat), - LonLat::new(self.min_lon, self.max_lat), - LonLat::new(self.min_lon, self.min_lat), - ] - } -} - -// Convenient wrapper around rstar -#[derive(Clone)] -pub struct QuadTree(RTree, T>>); - -impl QuadTree { - pub fn new() -> Self { - Self(RTree::new()) - } - - pub fn builder() -> QuadTreeBuilder { - QuadTreeBuilder::new() - } - - /// Slow, prefer bulk_load or QuadTreeBuilder - pub fn insert(&mut self, entry: GeomWithData, T>) { - self.0.insert(entry); - } - pub fn insert_with_box(&mut self, data: T, bbox: Bounds) { - self.0.insert(bbox.as_bbox(data)); - } - - pub fn bulk_load(entries: Vec, T>>) -> Self { - Self(RTree::bulk_load(entries)) - } - - pub fn query_bbox_borrow(&self, bbox: Bounds) -> impl Iterator + '_ { - let envelope = AABB::from_corners([bbox.min_x, bbox.min_y], [bbox.max_x, bbox.max_y]); - self.0 - .locate_in_envelope_intersecting(&envelope) - .map(|x| &x.data) - } -} - -impl QuadTree { - pub fn query_bbox(&self, bbox: Bounds) -> impl Iterator + '_ { - let envelope = AABB::from_corners([bbox.min_x, bbox.min_y], [bbox.max_x, bbox.max_y]); - self.0 - .locate_in_envelope_intersecting(&envelope) - .map(|x| x.data) - } -} - -impl QuadTree { - // TODO Inefficient - pub fn remove(&mut self, data: T) -> Option { - self.0 - .remove_with_selection_function(Selector(data)) - .map(|item| item.data) - } -} - -struct Selector(T); - -impl SelectionFunction, T>> for Selector { - fn should_unpack_parent(&self, _envelope: &AABB<[f64; 2]>) -> bool { - // TODO Maybe we can ask the caller to estimate where the previous object was? - true - } - - fn should_unpack_leaf(&self, leaf: &GeomWithData, T>) -> bool { - self.0 == leaf.data - } -} - -pub struct QuadTreeBuilder(Vec, T>>); - -impl QuadTreeBuilder { - pub fn new() -> Self { - Self(Vec::new()) - } - - pub fn add(&mut self, entry: GeomWithData, T>) { - self.0.push(entry); - } - pub fn add_with_box(&mut self, data: T, bbox: Bounds) { - self.0.push(bbox.as_bbox(data)); - } - - pub fn build(self) -> QuadTree { - QuadTree::bulk_load(self.0) - } -} diff --git a/geom/src/circle.rs b/geom/src/circle.rs deleted file mode 100644 index c63d77113c..0000000000 --- a/geom/src/circle.rs +++ /dev/null @@ -1,105 +0,0 @@ -use std::fmt; - -use anyhow::Result; -use serde::{Deserialize, Serialize}; - -use crate::{Angle, Bounds, Distance, Polygon, Pt2D, Ring, Tessellation}; - -const TRIANGLES_PER_CIRCLE: usize = 60; - -/// A circle, defined by a center and radius. -#[derive(Serialize, Deserialize, Debug, Clone)] -pub struct Circle { - pub center: Pt2D, - pub radius: Distance, -} - -impl Circle { - /// Creates a circle. - pub fn new(center: Pt2D, radius: Distance) -> Circle { - Circle { center, radius } - } - - /// True if the point is inside the circle. - pub fn contains_pt(&self, pt: Pt2D) -> bool { - // avoid sqrt by squaring radius instead - (pt.x() - self.center.x()).powi(2) + (pt.y() - self.center.y()).powi(2) - < self.radius.inner_meters().powi(2) - } - - /// Get the boundary containing this circle. - pub fn get_bounds(&self) -> Bounds { - Bounds { - min_x: self.center.x() - self.radius.inner_meters(), - max_x: self.center.x() + self.radius.inner_meters(), - min_y: self.center.y() - self.radius.inner_meters(), - max_y: self.center.y() + self.radius.inner_meters(), - } - } - - /// Renders the circle as a polygon. - pub fn to_polygon(&self) -> Polygon { - self.to_ring().into_polygon() - } - - /// Renders some percent, between [0, 1], of the circle. The shape starts from 0 degrees. - pub fn to_partial_tessellation(&self, percent_full: f64) -> Tessellation { - #![allow(clippy::float_cmp)] - assert!((0. ..=1.).contains(&percent_full)); - let mut pts = vec![self.center]; - let mut indices = Vec::new(); - for i in 0..TRIANGLES_PER_CIRCLE { - pts.push(self.center.project_away( - self.radius, - Angle::degrees((i as f64) / (TRIANGLES_PER_CIRCLE as f64) * percent_full * 360.0), - )); - indices.push(0); - indices.push(i + 1); - if i != TRIANGLES_PER_CIRCLE - 1 { - indices.push(i + 2); - } else if percent_full == 1.0 { - indices.push(1); - } else { - indices.pop(); - indices.pop(); - } - } - Tessellation::new(pts, indices) - } - - /// Returns the ring around the circle. - fn to_ring(&self) -> Ring { - let mut pts: Vec = (0..=TRIANGLES_PER_CIRCLE) - .map(|i| { - self.center.project_away( - self.radius, - Angle::degrees((i as f64) / (TRIANGLES_PER_CIRCLE as f64) * 360.0), - ) - }) - .collect(); - // With some radii, we get duplicate adjacent points - pts.dedup(); - Ring::must_new(pts) - } - - /// Creates an outline around the circle, strictly contained with the circle's original radius. - pub fn to_outline(&self, thickness: Distance) -> Result { - if self.radius <= thickness { - bail!( - "Can't make Circle outline with radius {} and thickness {}", - self.radius, - thickness - ); - } - - let bigger = self.to_ring(); - let smaller = Circle::new(self.center, self.radius - thickness).to_ring(); - Ok(Polygon::with_holes(bigger, vec![smaller])) - } -} - -impl fmt::Display for Circle { - fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { - write!(f, "Circle({}, {})", self.center, self.radius) - } -} diff --git a/geom/src/conversions.rs b/geom/src/conversions.rs deleted file mode 100644 index ce80bc421c..0000000000 --- a/geom/src/conversions.rs +++ /dev/null @@ -1,14 +0,0 @@ -//! Conversions between this crate and `geo`. Long-term, we should think about directly using `geo` -//! or wrapping it, but in the meantime... -//! -//! TODO Also, there's no consistency between standalone methods like this and From/Into impls. - -use crate::Pt2D; - -pub fn pts_to_line_string(raw_pts: &[Pt2D]) -> geo::LineString { - let pts: Vec = raw_pts - .iter() - .map(|pt| geo::Point::new(pt.x(), pt.y())) - .collect(); - pts.into() -} diff --git a/geom/src/distance.rs b/geom/src/distance.rs deleted file mode 100644 index a0a6952622..0000000000 --- a/geom/src/distance.rs +++ /dev/null @@ -1,309 +0,0 @@ -use std::{cmp, f64, fmt, ops}; - -use serde::{Deserialize, Serialize}; - -use crate::{deserialize_f64, serialize_f64, trim_f64, Duration, Speed, UnitFmt, EPSILON_DIST}; - -/// A distance, in meters. Can be negative. -#[derive(Clone, Copy, Debug, PartialEq, PartialOrd, Serialize, Deserialize)] -pub struct Distance( - #[serde(serialize_with = "serialize_f64", deserialize_with = "deserialize_f64")] f64, -); - -// By construction, Distance is a finite f64 with trimmed precision. -impl Eq for Distance {} - -#[allow(clippy::derive_ord_xor_partial_ord)] // false positive -impl Ord for Distance { - fn cmp(&self, other: &Distance) -> cmp::Ordering { - self.partial_cmp(other).unwrap() - } -} - -impl Distance { - pub const ZERO: Distance = Distance::const_meters(0.0); - - /// Creates a distance in meters. - pub fn meters(value: f64) -> Distance { - if !value.is_finite() { - panic!("Bad Distance {}", value); - } - - Distance(trim_f64(value)) - } - - // TODO Can't panic inside a const fn, seemingly. Don't pass in anything bad! - pub const fn const_meters(value: f64) -> Distance { - Distance(value) - } - - /// Creates a distance in inches. - pub fn inches(value: f64) -> Distance { - Distance::meters(0.0254 * value) - } - - /// Creates a distance in miles. - pub fn miles(value: f64) -> Distance { - Distance::meters(1609.34 * value) - } - - /// Creates a distance in centimeters. - pub fn centimeters(value: usize) -> Distance { - Distance::meters((value as f64) / 100.0) - } - - /// Creates a distance in feet. - pub fn feet(value: f64) -> Distance { - Distance::meters(value * 0.3048) - } - - /// Returns the absolute value of this distance. - pub fn abs(self) -> Distance { - if self.0 > 0.0 { - self - } else { - Distance(-self.0) - } - } - - /// Returns the square root of this distance. - pub fn sqrt(self) -> Distance { - Distance::meters(self.0.sqrt()) - } - - /// Returns the distance in meters. Prefer to work with type-safe `Distance`s. - // TODO Remove if possible. - pub fn inner_meters(self) -> f64 { - self.0 - } - - /// Returns the distance in feet. - pub fn to_feet(self) -> f64 { - self.0 * 3.28084 - } - - /// Returns the distance in miles. - pub fn to_miles(self) -> f64 { - self.to_feet() / 5280.0 - } - - /// Describes the distance according to formatting rules. Rounds to 1 decimal place for both - /// small (feet and meters) and large (miles and kilometers) units. - pub fn to_string(self, fmt: &UnitFmt) -> String { - if fmt.metric { - if self.0 < 1000.0 { - format!("{}m", (self.0 * 10.0).round() / 10.0) - } else { - let km = self.0 / 1000.0; - format!("{}km", (km * 10.0).round() / 10.0) - } - } else { - let feet = self.to_feet(); - let miles = self.to_miles(); - if miles >= 0.1 { - format!("{} miles", (miles * 10.0).round() / 10.0) - } else { - format!("{} ft", (feet * 10.0).round() / 10.0) - } - } - } - - /// Calculates a percentage, usually in [0.0, 1.0], of self / other. If the denominator is - /// zero, returns 0%. - pub fn safe_percent(self, other: Distance) -> f64 { - if other == Distance::ZERO { - return 0.0; - } - self / other - } - - /// Rounds this distance up to a higher, more "even" value to use for buckets along a plot's - /// axis. Always rounds for imperial units (feet). - pub fn round_up_for_axis(self) -> Distance { - let ft = self.to_feet(); - let miles = ft / 5280.0; - if ft <= 0.0 { - Distance::ZERO - } else if ft <= 10.0 { - Distance::feet(ft.ceil()) - } else if ft <= 100.0 { - Distance::feet(10.0 * (ft / 10.0).ceil()) - } else if miles < 0.1 { - Distance::feet(100.0 * (ft / 100.0).ceil()) - } else if miles <= 1.0 { - Distance::miles((miles * 10.0).ceil() / 10.0) - } else if miles <= 10.0 { - Distance::miles(miles.ceil()) - } else if miles <= 100.0 { - Distance::miles(10.0 * (miles / 10.0).ceil()) - } else { - self - } - } - - pub(crate) fn to_u64(self) -> u64 { - (self.0 / EPSILON_DIST.0) as u64 - } - - pub(crate) fn from_u64(x: u64) -> Distance { - (x as f64) * EPSILON_DIST - } -} - -impl fmt::Display for Distance { - fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { - write!(f, "{}m", self.0) - } -} - -impl ops::Add for Distance { - type Output = Distance; - - fn add(self, other: Distance) -> Distance { - Distance::meters(self.0 + other.0) - } -} - -impl ops::AddAssign for Distance { - fn add_assign(&mut self, other: Distance) { - *self = *self + other; - } -} - -impl ops::Sub for Distance { - type Output = Distance; - - fn sub(self, other: Distance) -> Distance { - Distance::meters(self.0 - other.0) - } -} - -impl ops::Neg for Distance { - type Output = Distance; - - fn neg(self) -> Distance { - Distance::meters(-self.0) - } -} - -impl ops::SubAssign for Distance { - fn sub_assign(&mut self, other: Distance) { - *self = *self - other; - } -} - -impl ops::Mul for Distance { - type Output = Distance; - - fn mul(self, scalar: f64) -> Distance { - Distance::meters(self.0 * scalar) - } -} - -impl ops::Mul for f64 { - type Output = Distance; - - fn mul(self, other: Distance) -> Distance { - Distance::meters(self * other.0) - } -} - -impl ops::MulAssign for Distance { - fn mul_assign(&mut self, other: f64) { - *self = *self * other; - } -} - -impl ops::Div for Distance { - type Output = f64; - - fn div(self, other: Distance) -> f64 { - if other == Distance::ZERO { - panic!("Can't divide {} / {}", self, other); - } - self.0 / other.0 - } -} - -impl ops::Div for Distance { - type Output = Distance; - - fn div(self, scalar: f64) -> Distance { - if scalar == 0.0 { - panic!("Can't divide {} / {}", self, scalar); - } - Distance::meters(self.0 / scalar) - } -} - -impl ops::Div for Distance { - type Output = Duration; - - fn div(self, other: Speed) -> Duration { - if other == Speed::ZERO { - panic!("Can't divide {} / 0 mph", self); - } - Duration::seconds(self.0 / other.inner_meters_per_second()) - } -} - -impl std::iter::Sum for Distance { - fn sum(iter: I) -> Distance - where - I: Iterator, - { - let mut sum = Distance::ZERO; - for x in iter { - sum += x; - } - sum - } -} - -impl Default for Distance { - fn default() -> Distance { - Distance::ZERO - } -} - -#[cfg(test)] -mod tests { - use super::*; - - #[test] - fn round_up_for_axis() { - let fmt = UnitFmt { - metric: false, - round_durations: false, - }; - - for (input, expected) in [ - (-3.0, 0.0), - (0.0, 0.0), - (3.2, 4.0), - (30.2, 40.0), - (300.2, 400.0), - ( - Distance::miles(0.13).to_feet(), - Distance::miles(0.2).to_feet(), - ), - ( - Distance::miles(0.64).to_feet(), - Distance::miles(0.7).to_feet(), - ), - ( - Distance::miles(2.6).to_feet(), - Distance::miles(3.0).to_feet(), - ), - ( - Distance::miles(2.9).to_feet(), - Distance::miles(3.0).to_feet(), - ), - ] { - assert_eq!( - Distance::feet(input).round_up_for_axis().to_string(&fmt), - Distance::feet(expected).to_string(&fmt) - ); - } - } -} diff --git a/geom/src/duration.rs b/geom/src/duration.rs deleted file mode 100644 index b9b78d2991..0000000000 --- a/geom/src/duration.rs +++ /dev/null @@ -1,438 +0,0 @@ -use std::{cmp, ops}; - -use anyhow::Result; -use instant::Instant; -use serde::{Deserialize, Serialize}; - -use crate::{deserialize_f64, serialize_f64, trim_f64, Distance, Speed, UnitFmt}; - -/// A duration, in seconds. Can be negative. -#[derive(Clone, Copy, Debug, PartialEq, PartialOrd, Serialize, Deserialize)] -pub struct Duration( - #[serde(serialize_with = "serialize_f64", deserialize_with = "deserialize_f64")] f64, -); - -// By construction, Duration is a finite f64 with trimmed precision. -impl Eq for Duration {} - -#[allow(clippy::derive_ord_xor_partial_ord)] // false positive -impl Ord for Duration { - fn cmp(&self, other: &Duration) -> cmp::Ordering { - self.partial_cmp(other).unwrap() - } -} - -impl Duration { - pub const ZERO: Duration = Duration::const_seconds(0.0); - pub const EPSILON: Duration = Duration::const_seconds(0.0001); - - /// Creates a duration in seconds. - pub fn seconds(value: f64) -> Duration { - if !value.is_finite() { - panic!("Bad Duration {}", value); - } - - Duration(trim_f64(value)) - } - - /// Creates a duration in minutes. - pub fn minutes(mins: usize) -> Duration { - Duration::seconds((mins as f64) * 60.0) - } - - /// Creates a duration in hours. - pub fn hours(hours: usize) -> Duration { - Duration::seconds((hours as f64) * 3600.0) - } - - /// Creates a duration in minutes. - pub fn f64_minutes(mins: f64) -> Duration { - Duration::seconds(mins * 60.0) - } - - /// Creates a duration in milliseconds. - pub fn milliseconds(value: f64) -> Duration { - Duration::seconds(value / 1000.0) - } - - pub const fn const_seconds(value: f64) -> Duration { - Duration(value) - } - - pub(crate) fn to_u64(self) -> u64 { - (self.0 / Duration::EPSILON.0) as u64 - } - - pub(crate) fn from_u64(x: u64) -> Duration { - (x as f64) * Duration::EPSILON - } - - pub fn abs(&self) -> Self { - Self(self.0.abs()) - } - - /// Returns the duration in seconds. Prefer working in typesafe `Duration`s. - // TODO Remove if possible. - pub fn inner_seconds(self) -> f64 { - self.0 - } - - /// Splits the duration into (hours, minutes, seconds, deciseconds). - // TODO Could share some of this with Time -- the representations are the same - fn get_parts(self) -> (usize, usize, usize, usize) { - // Force positive - let mut remainder = self.inner_seconds().abs(); - let hours = (remainder / 3600.0).floor(); - remainder -= hours * 3600.0; - let minutes = (remainder / 60.0).floor(); - remainder -= minutes * 60.0; - let seconds = remainder.floor(); - remainder -= seconds; - let decis = (remainder / 0.1).round(); - - ( - hours as usize, - minutes as usize, - seconds as usize, - decis as usize, - ) - } - - /// Parses a duration such as "3:00" to `Duration::minutes(3)`. - // TODO This is NOT the inverse of Display! - pub fn parse(string: &str) -> Result { - let parts: Vec<&str> = string.split(':').collect(); - if parts.is_empty() { - bail!("Duration {}: no :'s", string); - } - - let mut seconds: f64 = 0.0; - if parts.last().unwrap().contains('.') { - let last_parts: Vec<&str> = parts.last().unwrap().split('.').collect(); - if last_parts.len() != 2 { - bail!("Duration {}: no . in last part", string); - } - seconds += last_parts[1].parse::()? / 10.0; - seconds += last_parts[0].parse::()?; - } else { - seconds += parts.last().unwrap().parse::()?; - } - - match parts.len() { - 1 => Ok(Duration::seconds(seconds)), - 2 => { - seconds += 60.0 * parts[0].parse::()?; - Ok(Duration::seconds(seconds)) - } - 3 => { - seconds += 60.0 * parts[1].parse::()?; - seconds += 3600.0 * parts[0].parse::()?; - Ok(Duration::seconds(seconds)) - } - _ => bail!("Duration {}: weird number of parts", string), - } - } - - /// If two durations are within this amount, they'll print as if they're the same. - pub fn epsilon_eq(self, other: Duration) -> bool { - let eps = Duration::seconds(0.1); - match self.cmp(&other) { - cmp::Ordering::Greater => self - other < eps, - cmp::Ordering::Less => other - self < eps, - cmp::Ordering::Equal => true, - } - } - - /// Returns the duration elapsed from this moment in real time. - pub fn realtime_elapsed(since: Instant) -> Duration { - let dt = since.elapsed(); - Duration::seconds((dt.as_secs() as f64) + (f64::from(dt.subsec_nanos()) * 1e-9)) - } - - /// Rounds a duration up to the nearest whole number multiple. - pub fn round_up(self, multiple: Duration) -> Duration { - let remainder = self % multiple; - if remainder == Duration::ZERO { - self - } else { - self + multiple - remainder - } - } - - /// Returns the duration as a number of minutes, rounded up. - pub fn num_minutes_rounded_up(self) -> usize { - let (hrs, mins, secs, rem) = self.get_parts(); - let mut result = mins + 60 * hrs; - if secs != 0 || rem != 0 { - result += 1; - } - result - } - - // TODO Do something fancier? http://vis.stanford.edu/papers/tick-labels - // TODO Unit test me - /// Returns (rounded max, the boundaries) - pub fn make_intervals_for_max(self, num_labels: usize) -> (Duration, Vec) { - // Example 1: 43 minutes, max 5 labels... raw_mins_per_interval is 8.6 - let raw_mins_per_interval = (self.num_minutes_rounded_up() as f64) / (num_labels as f64); - // So then this rounded up to 10 minutes - let mut mins_per_interval = Duration::seconds(60.0 * raw_mins_per_interval) - .round_up(Duration::minutes(5)) - .num_minutes_rounded_up(); - - // Example 2: 8 minutes, max 5 labels... raw_mins_per_interval is 1.6 - // If we're under 25 minutes, this is going to be weird. - if self < (num_labels as f64) * Duration::minutes(5) { - // rounded up to 5 mins? 1 min increments - // up to 10? 2 min increments - // up to 15? 3 - // up to 20? 4 - // then after that the normal behavior - mins_per_interval = (self.round_up(Duration::minutes(5)) / (num_labels as f64)) - .num_minutes_rounded_up(); - } - - let max = (num_labels as f64) * Duration::minutes(mins_per_interval); - let labels = (0..=num_labels) - .map(|i| Duration::minutes(i * mins_per_interval)) - .collect(); - - if max < self { - panic!( - "Wait max of {} with {} labels wound up with rounded max of {}", - self, num_labels, max - ); - } - (max, labels) - } - - /// Shows only the largest unit (hours, minute, seconds), rounded to `precision` decimal points. - /// - /// ``` - /// use geom::Duration; - /// assert_eq!(Duration::seconds(3600.0).to_rounded_string(0), "1hr"); - /// assert_eq!(Duration::seconds(3600.0).to_rounded_string(1), "1.0hr"); - /// assert_eq!(Duration::seconds(7800.0).to_rounded_string(0), "2hr"); - /// assert_eq!(Duration::seconds(800.0).to_rounded_string(1), "13.3min"); - /// assert_eq!(Duration::seconds(-800.0).to_rounded_string(1), "-13.3min"); - /// assert_eq!(Duration::seconds(0.0).to_rounded_string(0), "0"); - /// assert_eq!(Duration::seconds(12.5).to_rounded_string(1), "12.5s"); - /// assert_eq!(Duration::seconds(12.5).to_rounded_string(2), "12.50s"); - /// ``` - pub fn to_rounded_string(self, precision: usize) -> String { - let (hours, minutes, seconds, remainder) = self.get_parts(); - if hours == 0 && minutes == 0 && seconds == 0 && remainder == 0 { - return "0".to_string(); - } - - let sign = if self < Duration::ZERO { "-" } else { "" }; - - let (whole, part, unit) = { - if hours != 0 { - let whole = hours as f64; - let part = minutes as f64 / 60.0; - let unit = "hr"; - (whole, part, unit) - } else if minutes != 0 { - let whole = minutes as f64; - let part = seconds as f64 / 60.0; - let unit = "min"; - (whole, part, unit) - } else { - let whole = seconds as f64; - let part = remainder as f64 / 10.0; - let unit = "s"; - (whole, part, unit) - } - }; - - let number = format!("{:.1$}", whole + part, precision); - return format!("{}{}{}", sign, number, unit); - } - - /// Describes the duration according to formatting rules. - pub fn to_string(self, fmt: &UnitFmt) -> String { - let mut s = String::new(); - if self < Duration::ZERO { - s = "-".to_string(); - } - let (hours, minutes, seconds, remainder) = self.get_parts(); - if hours == 0 && minutes == 0 && seconds == 0 && remainder == 0 { - return "0s".to_string(); - } - - if hours != 0 { - s = format!("{}{}hr ", s, hours); - } - if minutes != 0 { - s = format!("{}{}min ", s, minutes); - } - if remainder != 0 { - if fmt.round_durations { - s = format!("{}{}s", s, seconds); - } else { - s = format!("{}{}.{:1}s", s, seconds, remainder); - } - } else if seconds != 0 { - s = format!("{}{}s", s, seconds); - } - // Trim trailing whitespace, in case we have non-zero hours/minutes, but zero seconds - s.trim_end().to_string() - } -} - -impl std::fmt::Display for Duration { - fn fmt(&self, f: &mut std::fmt::Formatter) -> std::fmt::Result { - write!( - f, - "{}", - (*self).to_string(&UnitFmt { - metric: false, - round_durations: false - }) - ) - } -} - -impl ops::Add for Duration { - type Output = Duration; - - fn add(self, other: Duration) -> Duration { - Duration::seconds(self.0 + other.0) - } -} - -impl ops::AddAssign for Duration { - fn add_assign(&mut self, other: Duration) { - *self = *self + other; - } -} - -impl ops::SubAssign for Duration { - fn sub_assign(&mut self, other: Duration) { - *self = *self - other; - } -} - -impl ops::Sub for Duration { - type Output = Duration; - - fn sub(self, other: Duration) -> Duration { - Duration::seconds(self.0 - other.0) - } -} - -impl ops::Neg for Duration { - type Output = Duration; - - fn neg(self) -> Duration { - Duration::seconds(-self.0) - } -} - -impl ops::Mul for Duration { - type Output = Duration; - - fn mul(self, other: f64) -> Duration { - Duration::seconds(self.0 * other) - } -} - -// TODO Both of these work. Use a macro or crate to define both, so we don't have to worry about -// order for commutative things like multiplication. :P -impl ops::Mul for f64 { - type Output = Duration; - - fn mul(self, other: Duration) -> Duration { - Duration::seconds(self * other.0) - } -} - -impl ops::Mul for Duration { - type Output = Distance; - - fn mul(self, other: Speed) -> Distance { - Distance::meters(self.0 * other.inner_meters_per_second()) - } -} - -impl ops::Div for Duration { - type Output = f64; - - fn div(self, other: Duration) -> f64 { - if other.0 == 0.0 { - panic!("Can't divide {} / {}", self, other); - } - self.0 / other.0 - } -} - -impl ops::Div for Duration { - type Output = Duration; - - fn div(self, other: f64) -> Duration { - if other == 0.0 { - panic!("Can't divide {} / {}", self, other); - } - Duration::seconds(self.0 / other) - } -} - -impl ops::Rem for Duration { - type Output = Duration; - - fn rem(self, other: Duration) -> Duration { - Duration::seconds(self.0 % other.0) - } -} - -impl std::iter::Sum for Duration { - fn sum(iter: I) -> Duration - where - I: Iterator, - { - let mut sum = Duration::ZERO; - for x in iter { - sum += x; - } - sum - } -} - -impl Default for Duration { - fn default() -> Duration { - Duration::ZERO - } -} - -#[cfg(test)] -mod tests { - use super::*; - - #[test] - fn print_durations() { - let dont_round = UnitFmt { - metric: false, - round_durations: false, - }; - let round = UnitFmt { - metric: false, - round_durations: true, - }; - - assert_eq!("0s", Duration::ZERO.to_string(&dont_round)); - assert_eq!("0s", Duration::seconds(0.001).to_string(&dont_round)); - assert_eq!( - "1min 30.1s", - Duration::seconds(90.123).to_string(&dont_round) - ); - assert_eq!("1min 30s", Duration::seconds(90.123).to_string(&round)); - assert_eq!( - "2hr 33min 5s", - (Duration::hours(2) + Duration::minutes(33) + Duration::seconds(5.0)) - .to_string(&dont_round) - ); - assert_eq!("3hr", Duration::hours(3).to_string(&dont_round)); - assert_eq!("42min", Duration::minutes(42).to_string(&dont_round)); - } -} diff --git a/geom/src/find_closest.rs b/geom/src/find_closest.rs deleted file mode 100644 index dd04473e3c..0000000000 --- a/geom/src/find_closest.rs +++ /dev/null @@ -1,98 +0,0 @@ -use std::collections::BTreeMap; - -use geo::{ClosestPoint, Contains, EuclideanDistance, Intersects}; - -use crate::conversions::pts_to_line_string; -use crate::{Bounds, Distance, Polygon, Pt2D, QuadTree}; - -// TODO Maybe use https://crates.io/crates/spatial-join proximity maps - -/// A quad-tree to quickly find the closest points to some polylines. -#[derive(Clone)] -pub struct FindClosest { - // TODO maybe any type of geo:: thing - geometries: BTreeMap, - quadtree: QuadTree, -} - -impl FindClosest -where - K: Clone + Ord + std::fmt::Debug, -{ - pub fn new() -> FindClosest { - FindClosest { - geometries: BTreeMap::new(), - quadtree: QuadTree::new(), - } - } - - /// Add an object to the quadtree, remembering some key associated with the points. - /// TODO This doesn't properly handle single points, and will silently fail by never returning - /// any matches. - pub fn add(&mut self, key: K, pts: &[Pt2D]) { - self.geometries.insert(key.clone(), pts_to_line_string(pts)); - self.quadtree.insert_with_box(key, Bounds::from(pts)); - } - - /// Adds the outer ring of a polygon to the quadtree. - pub fn add_polygon(&mut self, key: K, polygon: &Polygon) { - self.add(key, polygon.get_outer_ring().points()); - } - - /// For every object within some distance of a query point, return the (object's key, point on - /// the object's polyline, distance away). - pub fn all_close_pts( - &self, - query_pt: Pt2D, - max_dist_away: Distance, - ) -> Vec<(K, Pt2D, Distance)> { - let query_geom = geo::Point::new(query_pt.x(), query_pt.y()); - - self.quadtree - .query_bbox_borrow( - Polygon::rectangle_centered(query_pt, max_dist_away, max_dist_away).get_bounds(), - ) - .filter_map(|key| { - if let geo::Closest::SinglePoint(pt) = - self.geometries[key].closest_point(&query_geom) - { - let dist = Distance::meters(pt.euclidean_distance(&query_geom)); - if dist <= max_dist_away { - Some((key.clone(), Pt2D::new(pt.x(), pt.y()), dist)) - } else { - None - } - } else if self.geometries[key].contains(&query_geom) { - // TODO Yay, FindClosest has a bug. :P - Some((key.clone(), query_pt, Distance::ZERO)) - } else { - None - } - }) - .collect() - } - - /// Finds the closest point on the existing geometry to the query pt. - pub fn closest_pt(&self, query_pt: Pt2D, max_dist_away: Distance) -> Option<(K, Pt2D)> { - self.all_close_pts(query_pt, max_dist_away) - .into_iter() - .min_by_key(|(_, _, dist)| *dist) - .map(|(k, pt, _)| (k, pt)) - } - - /// Find all objects with a point inside the query polygon - pub fn all_points_inside(&self, query: &Polygon) -> Vec { - let query_geo: geo::Polygon = query.clone().into(); - - self.quadtree - .query_bbox_borrow(query.get_bounds()) - .filter_map(|key| { - if self.geometries[key].intersects(&query_geo) { - Some(key.clone()) - } else { - None - } - }) - .collect() - } -} diff --git a/geom/src/gps.rs b/geom/src/gps.rs deleted file mode 100644 index e10c977841..0000000000 --- a/geom/src/gps.rs +++ /dev/null @@ -1,169 +0,0 @@ -use std::fmt; - -use anyhow::Result; -use geojson::{GeoJson, Value}; -use ordered_float::NotNan; -use serde::{Deserialize, Serialize}; - -use crate::{Distance, GPSBounds, Pt2D}; - -/// Represents a (longitude, latitude) point. -#[derive(Copy, Clone, PartialEq, Eq, PartialOrd, Ord, Debug, Serialize, Deserialize)] -pub struct LonLat { - longitude: NotNan, - latitude: NotNan, -} - -impl LonLat { - /// Note the order of arguments! - pub fn new(lon: f64, lat: f64) -> LonLat { - LonLat { - longitude: NotNan::new(lon).unwrap(), - latitude: NotNan::new(lat).unwrap(), - } - } - - /// Returns the longitude of this point. - pub fn x(self) -> f64 { - self.longitude.into_inner() - } - - /// Returns the latitude of this point. - pub fn y(self) -> f64 { - self.latitude.into_inner() - } - - /// Transform this to a world-space point. Can go out of bounds. - pub fn to_pt(self, b: &GPSBounds) -> Pt2D { - let (width, height) = { - let pt = b.get_max_world_pt(); - (pt.x(), pt.y()) - }; - - let x = (self.x() - b.min_lon) / (b.max_lon - b.min_lon) * width; - // Invert y, so that the northernmost latitude is 0. Screen drawing order, not Cartesian - // grid. - let y = height - ((self.y() - b.min_lat) / (b.max_lat - b.min_lat) * height); - Pt2D::new(x, y) - } - - /// Returns the Haversine distance to another point. - pub(crate) fn gps_dist(self, other: LonLat) -> Distance { - let earth_radius_m = 6_371_000.0; - let lon1 = self.x().to_radians(); - let lon2 = other.x().to_radians(); - let lat1 = self.y().to_radians(); - let lat2 = other.y().to_radians(); - - let delta_lat = lat2 - lat1; - let delta_lon = lon2 - lon1; - - let a = (delta_lat / 2.0).sin().powi(2) - + (delta_lon / 2.0).sin().powi(2) * lat1.cos() * lat2.cos(); - let c = 2.0 * a.sqrt().atan2((1.0 - a).sqrt()); - Distance::meters(earth_radius_m * c) - } - - /// Pretty meaningless units, for comparing distances very roughly - pub fn fast_dist(self, other: LonLat) -> NotNan { - NotNan::new((self.x() - other.x()).powi(2) + (self.y() - other.y()).powi(2)).unwrap() - } - - /// Finds the average of a set of coordinates. - pub fn center(pts: &[LonLat]) -> LonLat { - if pts.is_empty() { - panic!("Can't find center of 0 points"); - } - let mut x = 0.0; - let mut y = 0.0; - for pt in pts { - x += pt.x(); - y += pt.y(); - } - let len = pts.len() as f64; - LonLat::new(x / len, y / len) - } - - /// Parses a WKT-style line-string into a list of coordinates. - pub fn parse_wkt_linestring(raw: &str) -> Option> { - // Input is something like LINESTRING (-111.9263026 33.4245036, -111.9275146 33.4245016, - // -111.9278751 33.4233106) - let mut pts = Vec::new(); - // -111.9446 33.425474, -111.9442814 33.4254737, -111.9442762 33.426894 - for pair in raw - .strip_prefix("LINESTRING (")? - .strip_suffix(')')? - .split(", ") - { - let mut nums = Vec::new(); - for x in pair.split(' ') { - nums.push(x.parse::().ok()?); - } - if nums.len() != 2 { - return None; - } - pts.push(LonLat::new(nums[0], nums[1])); - } - if pts.len() < 2 { - return None; - } - Some(pts) - } - - /// Extract polygons from a raw GeoJSON string. For multipolygons, only returns the first - /// member. If the GeoJSON feature has a property called `name`, this will also be returned. - pub fn parse_geojson_polygons(raw: String) -> Result, Option)>> { - let geojson = raw.parse::()?; - let features = match geojson { - GeoJson::Feature(feature) => vec![feature], - GeoJson::FeatureCollection(feature_collection) => feature_collection.features, - _ => bail!("Unexpected geojson: {:?}", geojson), - }; - let mut polygons = Vec::new(); - for mut feature in features { - let points = match feature.geometry.take().map(|g| g.value) { - Some(Value::MultiPolygon(multi_polygon)) => multi_polygon[0][0].clone(), - Some(Value::Polygon(polygon)) => polygon[0].clone(), - _ => bail!("Unexpected feature: {:?}", feature), - }; - let name = feature - .property("name") - .and_then(|value| value.as_str()) - .map(|x| x.to_string()); - polygons.push(( - points - .into_iter() - .map(|pt| LonLat::new(pt[0], pt[1])) - .collect(), - name, - )); - } - Ok(polygons) - } - - /// Reads a GeoJSON file and returns coordinates from the one polygon contained. - pub fn read_geojson_polygon(path: &str) -> Result> { - let raw = fs_err::read_to_string(path)?; - let mut list = Self::parse_geojson_polygons(raw)?; - if list.len() != 1 { - bail!("{path} doesn't contain exactly one polygon"); - } - Ok(list.pop().unwrap().0) - } - - pub fn to_geojson(self) -> geojson::Geometry { - geojson::Geometry::new(geojson::Value::Point(vec![self.x(), self.y()])) - } -} - -impl fmt::Display for LonLat { - fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { - write!(f, "LonLat({0}, {1})", self.x(), self.y()) - } -} - -impl From for geo::Point { - fn from(pt: LonLat) -> Self { - geo::Point::new(pt.x(), pt.y()) - } -} diff --git a/geom/src/lib.rs b/geom/src/lib.rs deleted file mode 100644 index e120dbcf2b..0000000000 --- a/geom/src/lib.rs +++ /dev/null @@ -1,208 +0,0 @@ -#![allow(clippy::new_without_default)] - -#[macro_use] -extern crate anyhow; - -use serde::{Deserialize, Deserializer, Serialize, Serializer}; - -pub use crate::angle::Angle; -pub use crate::bounds::{Bounds, GPSBounds, QuadTree, QuadTreeBuilder}; -pub use crate::circle::Circle; -pub use crate::distance::Distance; -pub use crate::duration::Duration; -pub use crate::find_closest::FindClosest; -pub use crate::gps::LonLat; -pub use crate::line::{InfiniteLine, Line}; -pub use crate::percent::Percent; -pub use crate::polygon::Polygon; -pub use crate::polyline::{ArrowCap, PolyLine}; -pub use crate::pt::{HashablePt2D, Pt2D}; -pub use crate::ring::Ring; -pub use crate::speed::Speed; -pub use crate::stats::{HgramValue, Histogram, Statistic}; -pub use crate::tessellation::{Tessellation, Triangle}; -pub use crate::time::Time; - -mod angle; -mod bounds; -mod circle; -mod conversions; -mod distance; -mod duration; -mod find_closest; -mod gps; -mod line; -mod percent; -mod polygon; -mod polyline; -mod pt; -mod ring; -mod speed; -mod stats; -mod tessellation; -mod time; -mod utils; - -// About 0.4 inches... which is quite tiny on the scale of things. :) -pub const EPSILON_DIST: Distance = Distance::const_meters(0.01); - -/// Reduce the precision of an f64. This helps ensure serialization is idempotent (everything is -/// exactly the same before and after saving/loading). Ideally we'd use some kind of proper -/// fixed-precision type instead of f64. -pub fn trim_f64(x: f64) -> f64 { - (x * 10_000.0).round() / 10_000.0 -} - -/// Serializes a trimmed `f64` as an `i32` to save space. -fn serialize_f64(x: &f64, s: S) -> Result { - // So a trimmed f64's range becomes 2**31 / 10,000 =~ 214,000, which is plenty - // We MUST round here, the same as trim_f64. The unit test demonstrates why. - let int = (x * 10_000.0).round() as i32; - int.serialize(s) -} - -/// Deserializes a trimmed `f64` from an `i32`. -fn deserialize_f64<'de, D: Deserializer<'de>>(d: D) -> Result { - let x = ::deserialize(d)?; - Ok(x as f64 / 10_000.0) -} - -/// Specifies how to stringify different geom objects. -#[derive(Clone, Serialize, Deserialize, Copy)] -pub struct UnitFmt { - /// Round `Duration`s to a whole number of seconds. - pub round_durations: bool, - /// Display in metric; US imperial otherwise. - pub metric: bool, -} - -impl UnitFmt { - /// Default settings using metric. - pub fn metric() -> Self { - Self { - round_durations: true, - metric: true, - } - } - - /// Default settings using imperial. - pub fn imperial() -> Self { - Self { - round_durations: true, - metric: false, - } - } -} - -#[derive(Clone, Copy, Debug)] -pub struct CornerRadii { - pub top_left: f64, - pub top_right: f64, - pub bottom_right: f64, - pub bottom_left: f64, -} - -impl CornerRadii { - pub fn uniform(radius: f64) -> Self { - Self { - top_left: radius, - top_right: radius, - bottom_right: radius, - bottom_left: radius, - } - } - - pub fn zero() -> Self { - Self::uniform(0.0) - } -} - -impl std::convert::From for CornerRadii { - fn from(uniform: f64) -> Self { - Self::uniform(uniform) - } -} - -impl Default for CornerRadii { - fn default() -> Self { - Self::zero() - } -} - -/// Create a GeoJson with one feature per geometry, with the specified properties. -// TODO Rethink after https://github.com/georust/geojson/issues/170 -pub fn geometries_with_properties_to_geojson( - input: Vec<( - geojson::Geometry, - serde_json::Map, - )>, -) -> geojson::GeoJson { - let mut features = Vec::new(); - for (geom, properties) in input { - let mut f = geojson::Feature::from(geom); - f.properties = Some(properties); - features.push(f); - } - geojson::GeoJson::from(geojson::FeatureCollection { - bbox: None, - features, - foreign_members: None, - }) -} - -/// Create a GeoJson with one feature per geometry, and no properties. -pub fn geometries_to_geojson(input: Vec) -> geojson::GeoJson { - let mut features = Vec::new(); - for geom in input { - features.push(geojson::Feature::from(geom)); - } - geojson::GeoJson::from(geojson::FeatureCollection { - bbox: None, - features, - foreign_members: None, - }) -} - -#[cfg(test)] -mod tests { - use super::*; - use rand::{Rng, SeedableRng}; - - #[test] - fn f64_trimming() { - let mut rng = rand_xorshift::XorShiftRng::seed_from_u64(42); - for _ in 0..1_000 { - // Generate a random point - let input = Pt2D::new( - rng.gen_range(-214_000.00..214_000.0), - rng.gen_range(-214_000.00..214_000.0), - ); - // Round-trip to JSON and bincode - let json_roundtrip: Pt2D = - serde_json::from_slice(serde_json::to_string(&input).unwrap().as_bytes()).unwrap(); - let bincode_roundtrip: Pt2D = - bincode::deserialize(&bincode::serialize(&input).unwrap()).unwrap(); - - // Make sure everything is EXACTLY equal - if !exactly_eq(input, json_roundtrip) || !exactly_eq(input, bincode_roundtrip) { - panic!("Round-tripping mismatch!\ninput= {:?}\njson_roundtrip= {:?}\nbincode_roundtrip={:?}", input,json_roundtrip, bincode_roundtrip); - } - } - - // Hardcode a particular case, where we can hand-verify that it trims to 4 decimal places - let input = Pt2D::new(1.2345678, 9.876543); - let json_roundtrip: Pt2D = - serde_json::from_slice(serde_json::to_string(&input).unwrap().as_bytes()).unwrap(); - let bincode_roundtrip: Pt2D = - bincode::deserialize(&bincode::serialize(&input).unwrap()).unwrap(); - assert_eq!(input.x(), 1.2346); - assert_eq!(input.y(), 9.8765); - assert!(exactly_eq(input, json_roundtrip)); - assert!(exactly_eq(input, bincode_roundtrip)); - } - - // Don't use the PartialEq implementation, which does an epsilon check - fn exactly_eq(pt1: Pt2D, pt2: Pt2D) -> bool { - pt1.x() == pt2.x() && pt1.y() == pt2.y() - } -} diff --git a/geom/src/line.rs b/geom/src/line.rs deleted file mode 100644 index 804be032f2..0000000000 --- a/geom/src/line.rs +++ /dev/null @@ -1,279 +0,0 @@ -use std::fmt; - -use anyhow::Result; -use serde::{Deserialize, Serialize}; - -use crate::{Angle, Distance, PolyLine, Polygon, Pt2D, EPSILON_DIST}; - -/// A line segment. -// TODO Rename? -#[derive(Clone, Serialize, Deserialize, Debug, PartialEq)] -pub struct Line(Pt2D, Pt2D); - -impl Line { - /// Creates a line segment between two points, which must not be the same - pub fn new(pt1: Pt2D, pt2: Pt2D) -> Result { - if pt1.dist_to(pt2) <= EPSILON_DIST { - bail!("Line from {:?} to {:?} too small", pt1, pt2); - } - Ok(Line(pt1, pt2)) - } - - /// Equivalent to `Line::new(pt1, pt2).unwrap()`. Use this to effectively document an assertion - /// at the call-site. - pub fn must_new(pt1: Pt2D, pt2: Pt2D) -> Line { - Line::new(pt1, pt2).unwrap() - } - - /// Returns an infinite line passing through this line's two points. - pub fn infinite(&self) -> InfiniteLine { - InfiniteLine(self.0, self.1) - } - - /// Returns the first point in this line segment. - pub fn pt1(&self) -> Pt2D { - self.0 - } - - /// Returns the second point in this line segment. - pub fn pt2(&self) -> Pt2D { - self.1 - } - - /// Returns the two points in this line segment. - pub fn points(&self) -> Vec { - vec![self.0, self.1] - } - - /// Returns a polyline containing these two points. - pub fn to_polyline(&self) -> PolyLine { - PolyLine::must_new(self.points()) - } - - /// Returns a thick line segment. - pub fn make_polygons(&self, thickness: Distance) -> Polygon { - self.to_polyline().make_polygons(thickness) - } - - /// Length of the line segment - pub fn length(&self) -> Distance { - self.pt1().dist_to(self.pt2()) - } - - /// If two line segments intersect -- including endpoints -- return the point where they hit. - /// Undefined if the two lines have more than one intersection point! - // TODO Also return the distance along self - pub fn intersection(&self, other: &Line) -> Option { - // From http://bryceboe.com/2006/10/23/line-segment-intersection-algorithm/ - if is_counter_clockwise(self.pt1(), other.pt1(), other.pt2()) - == is_counter_clockwise(self.pt2(), other.pt1(), other.pt2()) - || is_counter_clockwise(self.pt1(), self.pt2(), other.pt1()) - == is_counter_clockwise(self.pt1(), self.pt2(), other.pt2()) - { - return None; - } - - let hit = self.infinite().intersection(&other.infinite())?; - if self.contains_pt(hit) { - // TODO and other contains pt, then we dont need ccw check thing - Some(hit) - } else { - // TODO Should be impossible, but I was hitting it somewhere - println!( - "{} and {} intersect, but first line doesn't contain_pt({})", - self, other, hit - ); - None - } - } - - /// Determine if two line segments intersect, but more so than just two endpoints touching. - pub fn crosses(&self, other: &Line) -> bool { - #[allow(clippy::suspicious_operation_groupings)] // false positive - if self.pt1() == other.pt1() - || self.pt1() == other.pt2() - || self.pt2() == other.pt1() - || self.pt2() == other.pt2() - { - return false; - } - self.intersection(other).is_some() - } - - /// If the line segment intersects with an infinite line -- including endpoints -- return the - /// point where they hit. Undefined if the segment and infinite line intersect at more than one - /// point! - // TODO Also return the distance along self - pub fn intersection_infinite(&self, other: &InfiniteLine) -> Option { - let hit = self.infinite().intersection(other)?; - if self.contains_pt(hit) { - Some(hit) - } else { - None - } - } - - /// Perpendicularly shifts the line over to the right. Width must be non-negative. - pub fn shift_right(&self, width: Distance) -> Line { - assert!(width >= Distance::ZERO); - let angle = self.angle().rotate_degs(90.0); - Line::must_new( - self.pt1().project_away(width, angle), - self.pt2().project_away(width, angle), - ) - } - - /// Perpendicularly shifts the line over to the left. Width must be non-negative. - pub fn shift_left(&self, width: Distance) -> Line { - assert!(width >= Distance::ZERO); - let angle = self.angle().rotate_degs(-90.0); - Line::must_new( - self.pt1().project_away(width, angle), - self.pt2().project_away(width, angle), - ) - } - - /// Perpendicularly shifts the line to the right if positive or left if negative. - pub fn shift_either_direction(&self, width: Distance) -> Line { - if width >= Distance::ZERO { - self.shift_right(width) - } else { - self.shift_left(-width) - } - } - - /// Returns a reversed line segment - pub fn reversed(&self) -> Line { - Line::must_new(self.pt2(), self.pt1()) - } - - /// The angle of the line segment, from the first to the second point - pub fn angle(&self) -> Angle { - self.pt1().angle_to(self.pt2()) - } - - /// Returns a point along the line segment, unless the distance exceeds the segment's length. - pub fn dist_along(&self, dist: Distance) -> Result { - let len = self.length(); - if dist < Distance::ZERO || dist > len { - bail!("dist_along({}) of a length {} line", dist, len); - } - self.percent_along(dist / len) - } - /// Equivalent to `self.dist_along(dist).unwrap()`. Use this to document an assertion at the - /// call-site. - pub fn must_dist_along(&self, dist: Distance) -> Pt2D { - self.dist_along(dist).unwrap() - } - - pub fn unbounded_dist_along(&self, dist: Distance) -> Pt2D { - self.unbounded_percent_along(dist / self.length()) - } - - pub fn unbounded_percent_along(&self, percent: f64) -> Pt2D { - Pt2D::new( - self.pt1().x() + percent * (self.pt2().x() - self.pt1().x()), - self.pt1().y() + percent * (self.pt2().y() - self.pt1().y()), - ) - } - pub fn percent_along(&self, percent: f64) -> Result { - if !(0.0..=1.0).contains(&percent) { - bail!("percent_along({}) of some line outside [0, 1]", percent); - } - Ok(self.unbounded_percent_along(percent)) - } - - pub fn slice(&self, from: Distance, to: Distance) -> Result { - if from < Distance::ZERO || to < Distance::ZERO || from >= to { - bail!("slice({}, {}) makes no sense", from, to); - } - Line::new(self.dist_along(from)?, self.dist_along(to)?) - } - - /// Returns a subset of this line, with two percentages along the line's length. - pub fn percent_slice(&self, from: f64, to: f64) -> Result { - self.slice(from * self.length(), to * self.length()) - } - - pub fn middle(&self) -> Result { - self.dist_along(self.length() / 2.0) - } - - pub fn contains_pt(&self, pt: Pt2D) -> bool { - self.dist_along_of_point(pt).is_some() - } - - pub fn dist_along_of_point(&self, pt: Pt2D) -> Option { - let dist1 = self.pt1().raw_dist_to(pt); - let dist2 = pt.raw_dist_to(self.pt2()); - let length = self.pt1().raw_dist_to(self.pt2()); - if (dist1 + dist2 - length).abs() < EPSILON_DIST.inner_meters() { - Some(Distance::meters(dist1)) - } else { - None - } - } - pub fn percent_along_of_point(&self, pt: Pt2D) -> Option { - let dist = self.dist_along_of_point(pt)?; - Some(dist / self.length()) - } -} - -impl fmt::Display for Line { - fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { - writeln!(f, "Line::new(")?; - writeln!(f, " Pt2D::new({}, {}),", self.0.x(), self.0.y())?; - writeln!(f, " Pt2D::new({}, {}),", self.1.x(), self.1.y())?; - write!(f, ")") - } -} - -fn is_counter_clockwise(pt1: Pt2D, pt2: Pt2D, pt3: Pt2D) -> bool { - (pt3.y() - pt1.y()) * (pt2.x() - pt1.x()) > (pt2.y() - pt1.y()) * (pt3.x() - pt1.x()) -} - -#[derive(Clone, Serialize, Deserialize, Debug)] -pub struct InfiniteLine(Pt2D, Pt2D); - -impl InfiniteLine { - /// Fails for parallel lines. - // https://stackoverflow.com/a/565282 by way of - // https://github.com/ucarion/line_intersection/blob/master/src/lib.rs - pub fn intersection(&self, other: &InfiniteLine) -> Option { - #![allow(clippy::many_single_char_names)] - fn cross(a: (f64, f64), b: (f64, f64)) -> f64 { - a.0 * b.1 - a.1 * b.0 - } - - let p = self.0; - let q = other.0; - let r = (self.1.x() - self.0.x(), self.1.y() - self.0.y()); - let s = (other.1.x() - other.0.x(), other.1.y() - other.0.y()); - - let r_cross_s = cross(r, s); - let q_minus_p = (q.x() - p.x(), q.y() - p.y()); - //let q_minus_p_cross_r = cross(q_minus_p, r); - - if r_cross_s == 0.0 { - // Parallel - None - } else { - let t = cross(q_minus_p, (s.0 / r_cross_s, s.1 / r_cross_s)); - //let u = cross(q_minus_p, Pt2D::new(r.x() / r_cross_s, r.y() / r_cross_s)); - Some(Pt2D::new(p.x() + t * r.0, p.y() + t * r.1)) - } - } - - pub fn from_pt_angle(pt: Pt2D, angle: Angle) -> InfiniteLine { - Line::must_new(pt, pt.project_away(Distance::meters(1.0), angle)).infinite() - } -} - -impl fmt::Display for InfiniteLine { - fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { - writeln!(f, "InfiniteLine::new(")?; - writeln!(f, " Pt2D::new({}, {}),", self.0.x(), self.0.y())?; - writeln!(f, " Pt2D::new({}, {}),", self.1.x(), self.1.y())?; - write!(f, ")") - } -} diff --git a/geom/src/percent.rs b/geom/src/percent.rs deleted file mode 100644 index 4d245eff5f..0000000000 --- a/geom/src/percent.rs +++ /dev/null @@ -1,28 +0,0 @@ -use std::fmt; - -/// Most of the time, [0, 1]. But some callers may go outside this range. -#[derive(Clone, Copy, PartialEq)] -pub struct Percent(f64); - -impl Percent { - pub fn inner(self) -> f64 { - self.0 - } - - pub fn int(x: usize) -> Percent { - if x > 100 { - panic!("Percent::int({}) too big", x); - } - Percent((x as f64) / 100.0) - } - - pub fn of(x: usize, total: usize) -> Percent { - Percent((x as f64) / (total as f64)) - } -} - -impl fmt::Display for Percent { - fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { - write!(f, "{:.2}%", self.0 * 100.0) - } -} diff --git a/geom/src/polygon.rs b/geom/src/polygon.rs deleted file mode 100644 index 44c82c6172..0000000000 --- a/geom/src/polygon.rs +++ /dev/null @@ -1,559 +0,0 @@ -use std::collections::BTreeMap; -use std::fmt; - -use anyhow::Result; -use geo::{ - Area, BooleanOps, Contains, ConvexHull, Intersects, MapCoordsInPlace, SimplifyVwPreserve, -}; -use serde::{Deserialize, Serialize}; - -use crate::{ - Angle, Bounds, CornerRadii, Distance, GPSBounds, LonLat, PolyLine, Pt2D, Ring, Tessellation, - Triangle, -}; - -#[derive(PartialEq, Serialize, Deserialize, Clone, Debug)] -pub struct Polygon { - // [0] is the outer/exterior ring, and the others are holes/interiors - pub(crate) rings: Vec, - // For performance reasons, some callers may want to generate the polygon's Rings and - // Tessellation at the same time, instead of using earcutr. - pub(crate) tessellation: Option, -} - -impl Polygon { - pub fn with_holes(outer: Ring, mut inner: Vec) -> Self { - inner.insert(0, outer); - Self { - rings: inner, - tessellation: None, - } - } - - pub fn from_rings(rings: Vec) -> Self { - assert!(!rings.is_empty()); - Self { - rings, - tessellation: None, - } - } - - pub(crate) fn pretessellated(rings: Vec, tessellation: Tessellation) -> Self { - Self { - rings, - tessellation: Some(tessellation), - } - } - - pub fn from_geojson(raw: &[Vec>]) -> Result { - let mut rings = Vec::new(); - for pts in raw { - let transformed: Vec = - pts.iter().map(|pair| Pt2D::new(pair[0], pair[1])).collect(); - rings.push(Ring::new(transformed)?); - } - Ok(Self::from_rings(rings)) - } - - pub fn from_triangle(tri: &Triangle) -> Self { - Ring::must_new(vec![tri.pt1, tri.pt2, tri.pt3, tri.pt1]).into_polygon() - } - - pub fn triangles(&self) -> Vec { - Tessellation::from(self.clone()).triangles() - } - - /// Does this polygon contain the point in its interior? - pub fn contains_pt(&self, pt: Pt2D) -> bool { - self.to_geo().contains(&geo::Point::from(pt)) - } - - pub fn get_bounds(&self) -> Bounds { - // Interiors should always be strictly contained within the polygon's exterior - Bounds::from(self.get_outer_ring().points()) - } - - /// Transformations must preserve Rings. - fn transform Pt2D>(&self, f: F) -> Result { - let mut rings = Vec::new(); - for ring in &self.rings { - rings.push(Ring::new(ring.points().iter().map(&f).collect())?); - } - Ok(Self { - rings, - tessellation: self.tessellation.clone().take().map(|mut t| { - t.transform(f); - t - }), - }) - } - - pub fn translate(&self, dx: f64, dy: f64) -> Self { - self.transform(|pt| pt.offset(dx, dy)) - .expect("translate shouldn't collapse Rings") - } - - /// When `factor` is small, this may collapse Rings and thus fail. - pub fn scale(&self, factor: f64) -> Result { - self.transform(|pt| Pt2D::new(pt.x() * factor, pt.y() * factor)) - } - - /// When `factor` is known to be over 1, then scaling can't fail. - pub fn must_scale(&self, factor: f64) -> Self { - if factor < 1.0 { - panic!("must_scale({factor}) might collapse Rings. Use scale()"); - } - self.transform(|pt| Pt2D::new(pt.x() * factor, pt.y() * factor)) - .expect("must_scale collapsed a Ring") - } - - pub fn rotate(&self, angle: Angle) -> Self { - self.rotate_around(angle, self.center()) - } - - pub fn rotate_around(&self, angle: Angle, pivot: Pt2D) -> Self { - self.transform(|pt| { - let origin_pt = Pt2D::new(pt.x() - pivot.x(), pt.y() - pivot.y()); - let (sin, cos) = angle.normalized_radians().sin_cos(); - Pt2D::new( - pivot.x() + origin_pt.x() * cos - origin_pt.y() * sin, - pivot.y() + origin_pt.y() * cos + origin_pt.x() * sin, - ) - }) - .expect("rotate_around shouldn't collapse Rings") - } - - pub fn centered_on(&self, center: Pt2D) -> Self { - let bounds = self.get_bounds(); - let dx = center.x() - bounds.width() / 2.0; - let dy = center.y() - bounds.height() / 2.0; - self.translate(dx, dy) - } - - pub fn get_outer_ring(&self) -> &Ring { - &self.rings[0] - } - - pub fn into_outer_ring(mut self) -> Ring { - self.rings.remove(0) - } - - /// Returns the arithmetic mean of the outer ring's points. The result could wind up inside a - /// hole in the polygon. Consider using `polylabel` too. - pub fn center(&self) -> Pt2D { - let mut pts = self.get_outer_ring().clone().into_points(); - pts.pop(); - Pt2D::center(&pts) - } - - /// Top-left at the origin. Doesn't take Distance, because this is usually pixels, actually. - pub fn maybe_rectangle(width: f64, height: f64) -> Result { - Ring::new(vec![ - Pt2D::new(0.0, 0.0), - Pt2D::new(width, 0.0), - Pt2D::new(width, height), - Pt2D::new(0.0, height), - Pt2D::new(0.0, 0.0), - ]) - .map(|ring| ring.into_polygon()) - } - - /// Top-left at the origin. Doesn't take Distance, because this is usually pixels, actually. - /// Note this will panic if `width` or `height` is 0. - pub fn rectangle(width: f64, height: f64) -> Self { - Self::maybe_rectangle(width, height).unwrap() - } - - pub fn rectangle_centered(center: Pt2D, width: Distance, height: Distance) -> Self { - Self::rectangle(width.inner_meters(), height.inner_meters()).translate( - center.x() - width.inner_meters() / 2.0, - center.y() - height.inner_meters() / 2.0, - ) - } - - pub fn rectangle_two_corners(pt1: Pt2D, pt2: Pt2D) -> Option { - if Pt2D::new(pt1.x(), 0.0) == Pt2D::new(pt2.x(), 0.0) - || Pt2D::new(0.0, pt1.y()) == Pt2D::new(0.0, pt2.y()) - { - return None; - } - - let (x1, width) = if pt1.x() < pt2.x() { - (pt1.x(), pt2.x() - pt1.x()) - } else { - (pt2.x(), pt1.x() - pt2.x()) - }; - let (y1, height) = if pt1.y() < pt2.y() { - (pt1.y(), pt2.y() - pt1.y()) - } else { - (pt2.y(), pt1.y() - pt2.y()) - }; - Some(Self::rectangle(width, height).translate(x1, y1)) - } - - /// Top-left at the origin. Doesn't take Distance, because this is usually pixels, actually. - pub fn maybe_rounded_rectangle>(w: f64, h: f64, r: R) -> Option { - let r = r.into(); - let max_r = r - .top_left - .max(r.top_right) - .max(r.bottom_right) - .max(r.bottom_left); - if 2.0 * max_r > w || 2.0 * max_r > h { - return None; - } - - let mut pts = vec![]; - - const RESOLUTION: usize = 5; - let mut arc = |r: f64, center: Pt2D, angle1_degs: f64, angle2_degs: f64| { - for i in 0..=RESOLUTION { - let angle = Angle::degrees( - angle1_degs + (angle2_degs - angle1_degs) * ((i as f64) / (RESOLUTION as f64)), - ); - pts.push(center.project_away(Distance::meters(r), angle.invert_y())); - } - }; - - arc(r.top_left, Pt2D::new(r.top_left, r.top_left), 180.0, 90.0); - arc( - r.top_right, - Pt2D::new(w - r.top_right, r.top_right), - 90.0, - 0.0, - ); - arc( - r.bottom_right, - Pt2D::new(w - r.bottom_right, h - r.bottom_right), - 360.0, - 270.0, - ); - arc( - r.bottom_left, - Pt2D::new(r.bottom_left, h - r.bottom_left), - 270.0, - 180.0, - ); - // Close it off - pts.push(Pt2D::new(0.0, r.top_left)); - - // If the radius was maximized, then some of the edges will be zero length. - pts.dedup(); - - Some(Ring::must_new(pts).into_polygon()) - } - - /// A rectangle, two sides of which are fully rounded half-circles. - pub fn pill(w: f64, h: f64) -> Self { - let r = w.min(h) / 2.0; - Self::maybe_rounded_rectangle(w, h, r).unwrap() - } - - /// Top-left at the origin. Doesn't take Distance, because this is usually pixels, actually. - /// If it's not possible to apply the specified radius, fallback to a regular rectangle. - pub fn rounded_rectangle>(w: f64, h: f64, r: R) -> Self { - Self::maybe_rounded_rectangle(w, h, r).unwrap_or_else(|| Self::rectangle(w, h)) - } - - /// Union all of the polygons into one geo::MultiPolygon - pub fn union_all_into_multipolygon(mut list: Vec) -> geo::MultiPolygon { - // TODO Not sure why this happened, or if this is really valid to construct... - if list.is_empty() { - return geo::MultiPolygon(Vec::new()); - } - - let mut result = geo::MultiPolygon(vec![list.pop().unwrap().into()]); - for p in list { - result = result.union(&p.into()); - } - result - } - - pub fn intersection(&self, other: &Self) -> Result> { - match std::panic::catch_unwind(std::panic::AssertUnwindSafe(|| { - from_multi(self.to_geo().intersection(&other.to_geo())) - })) { - Ok(result) => result, - Err(err) => { - println!("BooleanOps crashed: {err:?}"); - bail!("BooleanOps crashed: {err:?}"); - } - } - } - - pub fn difference(&self, other: &Self) -> Result> { - match std::panic::catch_unwind(std::panic::AssertUnwindSafe(|| { - from_multi(self.to_geo().difference(&other.to_geo())) - })) { - Ok(result) => result, - Err(err) => { - println!("BooleanOps crashed: {err:?}"); - bail!("BooleanOps crashed: {err:?}"); - } - } - } - - pub fn convex_hull(list: Vec) -> Result { - let mp: geo::MultiPolygon = list.into_iter().map(|p| p.to_geo()).collect(); - mp.convex_hull().try_into() - } - - pub fn concave_hull(points: Vec, concavity: u32) -> Result { - use geo::k_nearest_concave_hull::KNearestConcaveHull; - let points: Vec = points.iter().map(|p| geo::Point::from(*p)).collect(); - points.k_nearest_concave_hull(concavity).try_into() - } - - /// Find the "pole of inaccessibility" -- the most distant internal point from the polygon - /// outline - pub fn polylabel(&self) -> Pt2D { - let pt = polylabel::polylabel(&self.to_geo(), &1.0).unwrap(); - Pt2D::new(pt.x(), pt.y()) - } - - /// Do two polygons intersect at all? - pub fn intersects(&self, other: &Self) -> bool { - self.to_geo().intersects(&other.to_geo()) - } - - /// Does this polygon intersect a polyline? - pub fn intersects_polyline(&self, pl: &PolyLine) -> bool { - self.to_geo().intersects(&pl.to_geo()) - } - - pub(crate) fn get_rings(&self) -> &[Ring] { - &self.rings - } - - /// Creates the outline around the polygon (both the exterior and holes), with the thickness - /// half straddling the polygon and half of it just outside. - /// - /// Returns a `Tessellation` that may union together the outline from the exterior and multiple - /// holes. Callers that need a `Polygon` must call `to_outline` on the individual `Rings`. - pub fn to_outline(&self, thickness: Distance) -> Tessellation { - Tessellation::union_all( - self.rings - .iter() - .map(|r| Tessellation::from(r.to_outline(thickness))) - .collect(), - ) - } - - /// Usually m^2, unless the polygon is in screen-space - pub fn area(&self) -> f64 { - // Don't use signed_area, since we may work with polygons that have different orientations - self.to_geo().unsigned_area() - } - - /// Doesn't handle multiple crossings in and out. - pub fn clip_polyline(&self, input: &PolyLine) -> Option> { - let hits = self.get_outer_ring().all_intersections(input); - - if hits.is_empty() { - // All the points must be inside, or none - if self.contains_pt(input.first_pt()) { - Some(input.points().clone()) - } else { - None - } - } else if hits.len() == 1 { - // Which end? - if self.contains_pt(input.first_pt()) { - input - .get_slice_ending_at(hits[0]) - .map(|pl| pl.into_points()) - } else { - input - .get_slice_starting_at(hits[0]) - .map(|pl| pl.into_points()) - } - } else if hits.len() == 2 { - Some(input.trim_to_endpts(hits[0], hits[1]).into_points()) - } else { - // TODO Not handled - None - } - } - - // TODO Only handles a few cases - pub fn clip_ring(&self, input: &Ring) -> Option> { - let ring = self.get_outer_ring(); - let hits = ring.all_intersections(&PolyLine::unchecked_new(input.clone().into_points())); - - if hits.is_empty() { - // If the first point is inside, then all must be - if self.contains_pt(input.points()[0]) { - return Some(input.points().clone()); - } - } else if hits.len() == 2 { - let (pl1, pl2) = input.get_both_slices_btwn(hits[0], hits[1])?; - - // One of these should be partly outside the polygon. The endpoints won't be in the - // polygon itself, but they'll be on the ring. - if pl1 - .points() - .iter() - .all(|pt| self.contains_pt(*pt) || ring.contains_pt(*pt)) - { - return Some(pl1.into_points()); - } - if pl2 - .points() - .iter() - .all(|pt| self.contains_pt(*pt) || ring.contains_pt(*pt)) - { - return Some(pl2.into_points()); - } - // Huh? - } - - None - } - - /// Optionally map the world-space points back to GPS. - pub fn to_geojson(&self, gps: Option<&GPSBounds>) -> geojson::Geometry { - let mut geom: geo::Geometry = self.to_geo().into(); - if let Some(ref gps_bounds) = gps { - geom.map_coords_in_place(|c| { - let gps = Pt2D::new(c.x, c.y).to_gps(gps_bounds); - (gps.x(), gps.y()).into() - }); - } - - geojson::Geometry { - bbox: None, - value: geojson::Value::from(&geom), - foreign_members: None, - } - } - - /// Extracts all polygons from raw bytes representing a GeoJSON file, along with the string - /// key/value properties. Only the first polygon from multipolygons is returned. If - /// `require_in_bounds` is set, then the polygon must completely fit within the `gps_bounds`. - pub fn from_geojson_bytes( - raw_bytes: &[u8], - gps_bounds: &GPSBounds, - require_in_bounds: bool, - ) -> Result)>> { - let raw_string = std::str::from_utf8(raw_bytes)?; - let geojson = raw_string.parse::()?; - let features = match geojson { - geojson::GeoJson::Feature(feature) => vec![feature], - geojson::GeoJson::FeatureCollection(collection) => collection.features, - _ => bail!("Unexpected geojson: {:?}", geojson), - }; - - let mut results = Vec::new(); - for feature in features { - if let Some(geom) = &feature.geometry { - let raw_pts = match &geom.value { - geojson::Value::Polygon(pts) => pts, - // If there are multiple, just use the first - geojson::Value::MultiPolygon(polygons) => &polygons[0], - _ => { - continue; - } - }; - // TODO Handle holes - let gps_pts: Vec = raw_pts[0] - .iter() - .map(|pt| LonLat::new(pt[0], pt[1])) - .collect(); - let pts = if !require_in_bounds { - gps_bounds.convert(&gps_pts) - } else if let Some(pts) = gps_bounds.try_convert(&gps_pts) { - pts - } else { - continue; - }; - if let Ok(ring) = Ring::new(pts) { - let mut tags = BTreeMap::new(); - for (key, value) in feature.properties_iter() { - if let Some(value) = value.as_str() { - tags.insert(key.to_string(), value.to_string()); - } - } - results.push((ring.into_polygon(), tags)); - } - } - } - Ok(results) - } - - /// If simplification fails, just keep the original polygon - pub fn simplify(&self, epsilon: f64) -> Self { - self.to_geo() - .simplify_vw_preserve(&epsilon) - .try_into() - .unwrap_or_else(|_| self.clone()) - } - - /// An arbitrary placeholder value, when Option types aren't worthwhile - pub fn dummy() -> Self { - Self::rectangle(0.1, 0.1) - } - - // A less verbose way of invoking the From/Into impl. Note this hides a potentially expensive - // clone. - fn to_geo(&self) -> geo::Polygon { - self.clone().into() - } - - /// Convert to `geo` and also map from world-space to WGS84 - pub fn to_geo_wgs84(&self, gps_bounds: &GPSBounds) -> geo::Polygon { - let mut p = self.to_geo(); - p.map_coords_in_place(|c| { - let gps = Pt2D::new(c.x, c.y).to_gps(gps_bounds); - (gps.x(), gps.y()).into() - }); - p - } - - pub fn from_geo_wgs84(mut polygon: geo::Polygon, gps_bounds: &GPSBounds) -> Result { - polygon.map_coords_in_place(|c| { - let pt = LonLat::new(c.x, c.y).to_pt(gps_bounds); - (pt.x(), pt.y()).into() - }); - polygon.try_into() - } -} - -impl fmt::Display for Polygon { - fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { - writeln!(f, "Polygon with {} rings", self.rings.len())?; - for ring in &self.rings { - writeln!(f, " {}", ring)?; - } - Ok(()) - } -} - -impl TryFrom for Polygon { - type Error = anyhow::Error; - - fn try_from(poly: geo::Polygon) -> Result { - let (exterior, interiors) = poly.into_inner(); - let mut holes = Vec::new(); - for linestring in interiors { - holes.push(Ring::try_from(linestring)?); - } - Ok(Polygon::with_holes(Ring::try_from(exterior)?, holes)) - } -} - -impl From for geo::Polygon { - fn from(mut poly: Polygon) -> Self { - let exterior = poly.rings.remove(0); - let interiors: Vec = - poly.rings.into_iter().map(geo::LineString::from).collect(); - Self::new(exterior.into(), interiors) - } -} - -pub(crate) fn from_multi(multi: geo::MultiPolygon) -> Result> { - let mut result = Vec::new(); - for polygon in multi { - result.push(Polygon::try_from(polygon)?); - } - Ok(result) -} diff --git a/geom/src/polyline.rs b/geom/src/polyline.rs deleted file mode 100644 index e7f3179cb7..0000000000 --- a/geom/src/polyline.rs +++ /dev/null @@ -1,1225 +0,0 @@ -use std::collections::{BTreeMap, HashSet}; -use std::fmt; - -use anyhow::{Context, Result}; -use geo::{ClosestPoint, Winding}; -use serde::{Deserialize, Serialize}; - -use crate::{ - Angle, Bounds, Circle, Distance, GPSBounds, HashablePt2D, InfiniteLine, Line, LonLat, Polygon, - Pt2D, Ring, Tessellation, EPSILON_DIST, -}; - -// TODO How to tune this? -const MITER_THRESHOLD: f64 = 500.0; - -// TODO There used to be a second style that just has extra little hooks going out -pub enum ArrowCap { - Triangle, -} - -// TODO Document and enforce invariants: -// - at least two points -// - no duplicate points, whether adjacent or loops -// - no "useless" intermediate points with the same angle -#[derive(Clone, Debug, PartialEq, Serialize, Deserialize)] -pub struct PolyLine { - pts: Vec, - // TODO Note that caching length doesn't improve profiling results (by running - // small_spawn_completes test in release mode). May not be worth doing this. - length: Distance, -} - -impl PolyLine { - pub fn new(pts: Vec) -> Result { - if pts.len() < 2 { - bail!("Need at least two points for a PolyLine"); - } - let length = pts.windows(2).fold(Distance::ZERO, |so_far, pair| { - so_far + pair[0].dist_to(pair[1]) - }); - - if pts.windows(2).any(|pair| pair[0] == pair[1]) { - bail!( - "PL with total length {} and {} pts has ~dupe adjacent pts", - length, - pts.len(), - ); - } - - let result = PolyLine { pts, length }; - - // Can't have duplicates! If the polyline ever crosses back on itself, all sorts of things - // are broken. - let (_, dupes) = to_set(result.points()); - if !dupes.is_empty() { - bail!( - "PL with total length {} and {} pts has dupe non-adjacent pts", - result.length, - result.pts.len(), - ); - } - - Ok(result) - } - pub fn must_new(pts: Vec) -> PolyLine { - PolyLine::new(pts).unwrap() - } - - /// Doesn't check for duplicates. Use at your own risk. - pub fn unchecked_new(pts: Vec) -> PolyLine { - assert!(pts.len() >= 2); - let length = pts.windows(2).fold(Distance::ZERO, |so_far, pair| { - so_far + pair[0].dist_to(pair[1]) - }); - - PolyLine { pts, length } - } - - /// First dedupes adjacent points - pub fn deduping_new(mut pts: Vec) -> Result { - pts.dedup(); - PolyLine::new(pts) - } - - /// Like make_polygons, but make sure the points actually form a ring. - pub fn to_thick_ring(&self, width: Distance) -> Ring { - let mut side1 = self - .shift_with_sharp_angles(width / 2.0, MITER_THRESHOLD) - .unwrap(); - let mut side2 = self - .shift_with_sharp_angles(-width / 2.0, MITER_THRESHOLD) - .unwrap(); - side2.reverse(); - side1.extend(side2); - side1.push(side1[0]); - side1.dedup(); - Ring::must_new(side1) - } - - pub fn to_thick_boundary( - &self, - self_width: Distance, - boundary_width: Distance, - ) -> Option { - if self_width <= boundary_width || self.length() <= boundary_width + EPSILON_DIST { - return None; - } - // TODO exact_slice() used to work fine here, but the SUMO montlake map triggers a problem - // there - let slice = self - .maybe_exact_slice(boundary_width / 2.0, self.length() - boundary_width / 2.0) - .ok()?; - Some( - slice - .to_thick_ring(self_width - boundary_width) - .to_outline(boundary_width), - ) - } - - pub fn reversed(&self) -> PolyLine { - let mut pts = self.pts.clone(); - pts.reverse(); - PolyLine::must_new(pts) - } - - pub fn maybe_reverse(&self, reverse: bool) -> PolyLine { - if reverse { - self.reversed() - } else { - self.clone() - } - } - - /// Returns the quadrant where the overall angle of this polyline (pointing from the first to - /// last point) is in. Output between 0 and 3. - pub fn quadrant(&self) -> i64 { - let line_angle: f64 = self.overall_angle().normalized_radians(); - let line_angle = (line_angle / (std::f64::consts::PI / 2.0)) as i64; - line_angle.rem_euclid(4) + 1 - } - - /// Glue together two polylines in order. The last point of `self` must be the same as the - /// first point of `other`. This method handles removing unnecessary intermediate points if the - /// extension happens to be at the same angle as the last line segment of `self`. - pub fn extend(self, other: PolyLine) -> Result { - if *self.pts.last().unwrap() != other.pts[0] { - bail!("can't extend PL; last and first points don't match"); - } - - let mut self_pts = self.pts; - let mut other_pts = other.pts; - - loop { - let (pl1, _) = to_set(&self_pts); - let (pl2, _) = to_set(&other_pts[1..]); - - if pl1.intersection(&pl2).next().is_some() { - // Happens on some walking turns. Just clip out the loop. Start searching from the - // end of 'other'. - // TODO Measure the length of the thing being clipped out, to be sure this isn't - // running amok. - for (other_rev_idx, pt) in other_pts.iter().rev().enumerate() { - if pl1.contains(&pt.to_hashable()) { - while self_pts.last().unwrap() != pt { - self_pts.pop(); - } - other_pts = other_pts[other_pts.len() - 1 - other_rev_idx..].to_vec(); - break; - } - } - // Sanity check - assert_eq!(*self_pts.last().unwrap(), other_pts[0]); - } else { - break; - } - } - - // There's an exciting edge case: the next point to add is on self's last line. - if other_pts.len() >= 2 { - let same_line = self_pts[self_pts.len() - 2] - .angle_to(self_pts[self_pts.len() - 1]) - .approx_eq(other_pts[0].angle_to(other_pts[1]), 0.1); - if same_line { - self_pts.pop(); - } - } - self_pts.extend(other_pts.iter().skip(1)); - PolyLine::new(self_pts) - } - - /// Like `extend`, but panics on failure. - pub fn must_extend(self, other: PolyLine) -> PolyLine { - self.extend(other).unwrap() - } - - /// Extends `self` by a single point. If the new point is close enough to the last, dedupes. - /// Doesn't clean up any intermediate points. - pub fn optionally_push(self, pt: Pt2D) -> PolyLine { - let orig = self.clone(); - let mut pts = self.into_points(); - pts.push(pt); - match PolyLine::deduping_new(pts) { - Ok(pl) => pl, - // If the polyline loops back on itself and someone manages to exactly repeat an - // earlier point, just don't add this point - Err(_) => orig, - } - } - - /// Like `extend`, but handles the last and first point not matching by inserting that point. - /// Doesn't clean up any intermediate points. - pub fn force_extend(mut self, other: PolyLine) -> Result { - if *self.pts.last().unwrap() != other.pts[0] { - // TODO Blindly... what if we need to do the angle collapsing? - self.pts.push(other.pts[0]); - } - self.extend(other) - } - - /// One or both args might be empty. - pub fn append(first: Vec, second: Vec) -> Result> { - if second.is_empty() { - return Ok(first); - } - if first.is_empty() { - return Ok(second); - } - - Ok(PolyLine::new(first)? - .extend(PolyLine::new(second)?)? - .into_points()) - } - - pub fn points(&self) -> &Vec { - &self.pts - } - pub fn into_points(self) -> Vec { - self.pts - } - - pub fn lines(&self) -> impl Iterator + '_ { - self.pts - .windows(2) - .map(|pair| Line::must_new(pair[0], pair[1])) - } - - pub fn length(&self) -> Distance { - self.length - } - - /// Returns the excess distance left over from the end - pub fn slice(&self, start: Distance, end: Distance) -> Result<(PolyLine, Distance)> { - if start > end || start < Distance::ZERO || end < Distance::ZERO { - bail!("Can't get a polyline slice [{}, {}]", start, end); - } - if start > self.length() { - bail!( - "Can't get a polyline slice [{}, {}] on something of length {}", - start, - end, - self.length() - ); - } - if end - start < EPSILON_DIST { - bail!( - "Can't get a polyline slice [{}, {}] -- too small", - start, - end - ); - } - - let mut result: Vec = Vec::new(); - let mut dist_so_far = Distance::ZERO; - - for line in self.lines() { - let length = line.length(); - - // Does this line contain the first point of the slice? - if result.is_empty() && dist_so_far + length >= start { - result.push(line.must_dist_along(start - dist_so_far)); - } - - // Does this line contain the last point of the slice? - if dist_so_far + length >= end { - let last_pt = line.must_dist_along(end - dist_so_far); - if *result.last().unwrap() == last_pt { - result.pop(); - } - result.push(last_pt); - if result.len() == 1 { - // TODO Understand what happened here. - bail!("slice({}, {}) on {} did something weird", start, end, self); - } - return Ok((PolyLine::new(result)?, Distance::ZERO)); - } - - // If we're in the middle, just collect the endpoint. But not if it's too close to the - // previous point (namely, the start, which could be somewhere far along a line) - if !result.is_empty() && *result.last().unwrap() != line.pt2() { - result.push(line.pt2()); - } - - dist_so_far += length; - } - - if result.is_empty() { - bail!( - "Slice [{}, {}] has a start too big for polyline of length {}", - start, - end, - self.length() - ); - } - if result.len() == 1 { - bail!( - "Slice [{}, {}] on {} wound up a single point", - start, - end, - self - ); - } - - Ok((PolyLine::new(result)?, end - dist_so_far)) - } - - /// No excess leftover distance allowed. - // TODO Lot of callers of this. Make safer later. - pub fn exact_slice(&self, start: Distance, end: Distance) -> PolyLine { - self.maybe_exact_slice(start, end).unwrap() - } - pub fn maybe_exact_slice(&self, start: Distance, end: Distance) -> Result { - let (pl, leftover) = self - .slice(start, end) - .with_context(|| format!("exact_slice({}, {}) yielded empty slice", start, end))?; - if leftover > EPSILON_DIST { - bail!( - "exact_slice({}, {}) on a PL of length {} yielded leftover distance of {}", - start, - end, - self.length(), - leftover - ); - } - Ok(pl) - } - - pub fn first_half(&self) -> Result { - self.maybe_exact_slice(Distance::ZERO, self.length() / 2.0) - } - - pub fn second_half(&self) -> Result { - self.maybe_exact_slice(self.length() / 2.0, self.length()) - } - - pub fn dist_along(&self, dist_along: Distance) -> Result<(Pt2D, Angle)> { - if dist_along < Distance::ZERO { - bail!("dist_along {} is negative", dist_along); - } - if dist_along > self.length() { - bail!("dist_along {} is longer than {}", dist_along, self.length()); - } - if dist_along == self.length() { - return Ok((self.last_pt(), self.last_line().angle())); - } - - let mut dist_left = dist_along; - for (idx, l) in self.lines().enumerate() { - let length = l.length(); - let epsilon = if idx == self.pts.len() - 2 { - EPSILON_DIST - } else { - Distance::ZERO - }; - if dist_left <= length + epsilon { - // Floating point errors means sometimes we ask for something slightly longer than - // the line - let dist = l.dist_along(dist_left).unwrap_or_else(|_| l.pt2()); - return Ok((dist, l.angle())); - } - dist_left -= length; - } - // Leaving this panic, because I haven't seen this in ages, and something is seriously - // wrong if we get here - panic!( - "PolyLine dist_along of {} broke on length {}: {}", - dist_along, - self.length(), - self - ); - } - pub fn must_dist_along(&self, dist_along: Distance) -> (Pt2D, Angle) { - self.dist_along(dist_along).unwrap() - } - - pub fn middle(&self) -> Pt2D { - // If this fails, must be some super tiny line. Just return the first point in that case. - match self.dist_along(self.length() / 2.0) { - Ok((pt, _)) => pt, - Err(err) => { - println!( - "Guessing middle of PL with length {}: {}", - self.length(), - err - ); - self.first_pt() - } - } - } - - pub fn first_pt(&self) -> Pt2D { - self.pts[0] - } - pub fn last_pt(&self) -> Pt2D { - *self.pts.last().unwrap() - } - pub fn first_line(&self) -> Line { - Line::must_new(self.pts[0], self.pts[1]) - } - pub fn last_line(&self) -> Line { - Line::must_new(self.pts[self.pts.len() - 2], self.pts[self.pts.len() - 1]) - } - - pub fn shift_right(&self, width: Distance) -> Result { - self.shift_with_corrections(width) - } - pub fn must_shift_right(&self, width: Distance) -> PolyLine { - self.shift_right(width).unwrap() - } - - pub fn shift_left(&self, width: Distance) -> Result { - self.shift_with_corrections(-width) - } - pub fn must_shift_left(&self, width: Distance) -> PolyLine { - self.shift_left(width).unwrap() - } - - /// Perpendicularly shifts the polyline to the right if positive or left if negative. - pub fn shift_either_direction(&self, width: Distance) -> Result { - self.shift_with_corrections(width) - } - - /// `self` represents some center, with `total_width`. Logically this shifts left by - /// `total_width / 2`, then right by `width_from_left_side`, but without exasperating sharp - /// bends. - pub fn shift_from_center( - &self, - total_width: Distance, - width_from_left_side: Distance, - ) -> Result { - let half_width = total_width / 2.0; - if width_from_left_side < half_width { - self.shift_left(half_width - width_from_left_side) - } else { - self.shift_right(width_from_left_side - half_width) - } - } - - // Things to remember about shifting polylines: - // - the length before and after probably don't match up - // - the number of points may not match - fn shift_with_corrections(&self, width: Distance) -> Result { - let raw = self.shift_with_sharp_angles(width, MITER_THRESHOLD)?; - let result = PolyLine::deduping_new(raw)?; - if result.pts.len() == self.pts.len() { - fix_angles(self, result) - } else { - Ok(result) - } - } - - // If we start with a valid PolyLine, I'm not sure how we can ever possibly fail here, but it's - // happening. Avoid crashing. - fn shift_with_sharp_angles(&self, width: Distance, miter_threshold: f64) -> Result> { - if self.pts.len() == 2 { - let l = Line::new(self.pts[0], self.pts[1])?.shift_either_direction(width); - return Ok(vec![l.pt1(), l.pt2()]); - } - - let mut result: Vec = Vec::new(); - - let mut pt3_idx = 2; - let mut pt1_raw = self.pts[0]; - let mut pt2_raw = self.pts[1]; - - loop { - let pt3_raw = self.pts[pt3_idx]; - - let l1 = Line::new(pt1_raw, pt2_raw)?.shift_either_direction(width); - let l2 = Line::new(pt2_raw, pt3_raw)?.shift_either_direction(width); - - if pt3_idx == 2 { - result.push(l1.pt1()); - } - - if let Some(pt2_shift) = l1.infinite().intersection(&l2.infinite()) { - // Miter caps sometimes explode out to infinity. Hackily work around this. - let dist_away = l1.pt1().raw_dist_to(pt2_shift); - if dist_away < miter_threshold { - result.push(pt2_shift); - } else { - result.push(l1.pt2()); - } - } else { - // When the lines are perfectly parallel, it means pt2_shift_1st == pt2_shift_2nd - // and the original geometry is redundant. - result.push(l1.pt2()); - } - if pt3_idx == self.pts.len() - 1 { - result.push(l2.pt2()); - break; - } - - pt1_raw = pt2_raw; - pt2_raw = pt3_raw; - pt3_idx += 1; - } - - assert!(result.len() == self.pts.len()); - Ok(result) - } - - /// This produces a `Polygon` with a valid `Ring`. It may crash if this polyline was made with - /// `unchecked_new` - pub fn make_polygons(&self, width: Distance) -> Polygon { - let tessellation = self.thicken_tessellation(width); - let ring = Ring::deduping_new(tessellation.points.clone()) - .expect("PolyLine::make_polygons() failed"); - Polygon::pretessellated(vec![ring], tessellation) - } - - /// This produces a `Polygon` with a possibly invalid `Ring`. It shouldn't crash. Use sparingly. - pub fn unsafe_make_polygons(&self, width: Distance) -> Polygon { - let tessellation = self.thicken_tessellation(width); - let ring = Ring::unsafe_deduping_new(tessellation.points.clone()) - .expect("PolyLine::make_polygons() failed"); - Polygon::pretessellated(vec![ring], tessellation) - } - - /// Just produces a Tessellation - pub fn thicken_tessellation(&self, width: Distance) -> Tessellation { - // TODO Don't use the angle corrections yet -- they seem to do weird things. - let side1 = match self.shift_with_sharp_angles(width / 2.0, MITER_THRESHOLD) { - Ok(pl) => pl, - Err(err) => { - // TODO Circles will look extremely bizarre, but it emphasizes there's a bug - // without just crashing - println!("make_polygons({}) of {:?} failed: {}", width, self, err); - return Tessellation::from(Circle::new(self.first_pt(), width).to_polygon()); - } - }; - let mut side2 = match self.shift_with_sharp_angles(-width / 2.0, MITER_THRESHOLD) { - Ok(pl) => pl, - Err(err) => { - println!("make_polygons({}) of {:?} failed: {}", width, self, err); - return Tessellation::from(Circle::new(self.first_pt(), width).to_polygon()); - } - }; - assert_eq!(side1.len(), side2.len()); - - // Order the points so that they form a ring. No deduplication yet, though. - let len = 2 * side1.len(); - let mut points = side1; - side2.reverse(); - points.extend(side2); - points.push(points[0]); - let mut indices = Vec::new(); - - // Walk along the first side, making two triangles each step. This is easy to understand - // with a simple diagram, which I should eventually draw in ASCII art here. - for high_idx in 1..self.pts.len() { - indices.extend(vec![high_idx, high_idx - 1, len - high_idx]); - indices.extend(vec![len - high_idx, len - high_idx - 1, high_idx]); - } - - Tessellation::new(points.clone(), indices) - } - - /// This does the equivalent of make_polygons, returning the (start left, start right, end - /// left, end right). Fails rarely. - pub fn get_four_corners_of_thickened( - &self, - width: Distance, - ) -> Option<(Pt2D, Pt2D, Pt2D, Pt2D)> { - let mut side1 = self - .shift_with_sharp_angles(width / 2.0, MITER_THRESHOLD) - .ok()?; - let mut side2 = self - .shift_with_sharp_angles(-width / 2.0, MITER_THRESHOLD) - .ok()?; - Some(( - side1[0], - side2[0], - side1.pop().unwrap(), - side2.pop().unwrap(), - )) - } - - pub fn exact_dashed_polygons( - &self, - width: Distance, - dash_len: Distance, - dash_separation: Distance, - ) -> Vec { - let mut polygons: Vec = Vec::new(); - - let total_length = self.length(); - - let mut start = Distance::ZERO; - loop { - if start + dash_len >= total_length { - break; - } - - polygons.push( - self.exact_slice(start, start + dash_len) - .make_polygons(width), - ); - start += dash_len + dash_separation; - } - - polygons - } - - /// Don't draw the dashes too close to the ends. - pub fn dashed_lines( - &self, - width: Distance, - dash_len: Distance, - dash_separation: Distance, - ) -> Vec { - if self.length() <= dash_separation * 2.0 + EPSILON_DIST { - return vec![self.make_polygons(width)]; - } - self.exact_slice(dash_separation, self.length() - dash_separation) - .exact_dashed_polygons(width, dash_len, dash_separation) - } - - /// Fail if the length is too short. - pub fn maybe_make_arrow(&self, thickness: Distance, cap: ArrowCap) -> Option { - let head_size = thickness * 2.0; - let triangle_height = head_size / 2.0_f64.sqrt(); - - let slice = self - .maybe_exact_slice(Distance::ZERO, self.length() - triangle_height) - .ok()?; - - let angle = slice.last_pt().angle_to(self.last_pt()); - let corner1 = self - .last_pt() - .project_away(head_size, angle.rotate_degs(-135.0)); - let corner2 = self - .last_pt() - .project_away(head_size, angle.rotate_degs(135.0)); - - let mut pts = slice - .shift_with_sharp_angles(thickness / 2.0, MITER_THRESHOLD) - .ok()?; - match cap { - ArrowCap::Triangle => { - pts.push(corner2); - pts.push(self.last_pt()); - pts.push(corner1); - } - } - let mut side2 = slice - .shift_with_sharp_angles(-thickness / 2.0, MITER_THRESHOLD) - .ok()?; - side2.reverse(); - pts.extend(side2); - pts.push(pts[0]); - pts.dedup(); - Some(Ring::must_new(pts).into_polygon()) - } - - /// If the length is too short, just give up and make the thick line - pub fn make_arrow(&self, thickness: Distance, cap: ArrowCap) -> Polygon { - if let Some(p) = self.maybe_make_arrow(thickness, cap) { - p - } else { - // Just give up and make the thick line. - self.make_polygons(thickness) - } - } - - pub fn make_double_arrow(&self, thickness: Distance, cap: ArrowCap) -> Polygon { - let head_size = thickness * 2.0; - let triangle_height = head_size / 2.0_f64.sqrt(); - - if self.length() < triangle_height * 2.0 + EPSILON_DIST { - // Just give up and make the thick line. - return self.make_polygons(thickness); - } - let slice = self.exact_slice(triangle_height, self.length() - triangle_height); - - let angle = slice.last_pt().angle_to(self.last_pt()); - let corner1 = self - .last_pt() - .project_away(head_size, angle.rotate_degs(-135.0)); - let corner2 = self - .last_pt() - .project_away(head_size, angle.rotate_degs(135.0)); - - let mut pts = match slice.shift_with_sharp_angles(thickness / 2.0, MITER_THRESHOLD) { - Ok(pl) => pl, - Err(_) => { - return self.make_polygons(thickness); - } - }; - match cap { - ArrowCap::Triangle => { - pts.push(corner2); - pts.push(self.last_pt()); - pts.push(corner1); - } - } - let mut side2 = match slice.shift_with_sharp_angles(-thickness / 2.0, MITER_THRESHOLD) { - Ok(pl) => pl, - Err(_) => { - return self.make_polygons(thickness); - } - }; - side2.reverse(); - pts.extend(side2); - - let angle = self.first_pt().angle_to(slice.first_pt()); - let corner3 = self - .first_pt() - .project_away(head_size, angle.rotate_degs(-45.0)); - let corner4 = self - .first_pt() - .project_away(head_size, angle.rotate_degs(45.0)); - match cap { - ArrowCap::Triangle => { - pts.push(corner3); - pts.push(self.first_pt()); - pts.push(corner4); - } - } - - pts.push(pts[0]); - pts.dedup(); - Ring::must_new(pts).into_polygon() - } - - pub fn dashed_arrow( - &self, - width: Distance, - dash_len: Distance, - dash_separation: Distance, - cap: ArrowCap, - ) -> Vec { - let mut polygons = self.exact_dashed_polygons(width, dash_len, dash_separation); - // And a cap on the arrow. In case the last line is long, trim it to be the dash - // length. - let last_line = self.last_line(); - let last_len = last_line.length(); - let arrow_line = if last_len <= dash_len { - last_line - } else { - Line::must_new( - last_line.must_dist_along(last_len - dash_len), - last_line.pt2(), - ) - }; - polygons.push(arrow_line.to_polyline().make_arrow(width, cap)); - polygons - } - - /// Also return the angle of the line where the hit was found - // TODO Also return distance along self of the hit - pub fn intersection(&self, other: &PolyLine) -> Option<(Pt2D, Angle)> { - if self == other { - // This is likely a much biggger problem for the caller. Previously we would crash, but - // that's often inappropriate. - println!("intersection() between two PolyLines that're exactly the same"); - return None; - } - - // There could be several collisions. Pick the "first" from self's perspective. - let mut closest_intersection: Option<(Pt2D, Angle)> = None; - let mut closest_intersection_distance: Option = None; - - for l1 in self.lines() { - for l2 in other.lines() { - if let Some(pt) = l1.intersection(&l2) { - if let Some(new_distance) = self.get_slice_ending_at(pt).map(|pl| pl.length()) { - match closest_intersection_distance { - None => { - closest_intersection = Some((pt, l1.angle())); - closest_intersection_distance = Some(new_distance); - } - Some(existing_distance) if existing_distance > new_distance => { - closest_intersection = Some((pt, l1.angle())); - closest_intersection_distance = Some(new_distance); - } - _ => {} - } - } - } - } - } - - // TODO Why is any of this necessary? Found a test case at the intersection geometry for - // https://www.openstreetmap.org/node/274088813 where this made a huge difference! - if closest_intersection.is_none() && self.last_pt() == other.last_pt() { - return Some((self.last_pt(), self.last_line().angle())); - } - - closest_intersection - } - - // TODO Also distance along - pub fn intersection_infinite(&self, other: &InfiniteLine) -> Option { - for l in self.lines() { - if let Some(hit) = l.intersection_infinite(other) { - return Some(hit); - } - } - None - } - - /// Panics if the pt is not on the polyline. Returns None if the point is the first point - /// (meaning the slice is empty). - pub fn get_slice_ending_at(&self, pt: Pt2D) -> Option { - if self.first_pt() == pt { - return None; - } - - if let Some(idx) = self.lines().position(|l| l.contains_pt(pt)) { - let mut pts = self.pts.clone(); - pts.truncate(idx + 1); - // Make sure the last line isn't too tiny - if *pts.last().unwrap() == pt { - pts.pop(); - } - pts.push(pt); - if pts.len() == 1 { - return None; - } - Some(PolyLine::must_new(pts)) - } else { - panic!("Can't get_slice_ending_at: {} doesn't contain {}", self, pt); - } - } - - /// Returns None if the point is the last point. - pub fn get_slice_starting_at(&self, pt: Pt2D) -> Option { - if self.last_pt() == pt { - return None; - } - - if let Some(idx) = self.lines().position(|l| l.contains_pt(pt)) { - let mut pts = self.pts.clone(); - pts = pts.split_off(idx + 1); - if pt != pts[0] { - pts.insert(0, pt); - } - Some(PolyLine::must_new(pts)) - } else { - panic!( - "Can't get_slice_starting_at: {} doesn't contain {}", - self, pt - ); - } - } - - /// Same as get_slice_ending_at, but returns None if the point isn't on the polyline. - // TODO Switch everything to this, after better understanding why this is happening at all. - pub fn safe_get_slice_ending_at(&self, pt: Pt2D) -> Option { - if self.first_pt() == pt { - return None; - } - - if let Some(idx) = self.lines().position(|l| l.contains_pt(pt)) { - let mut pts = self.pts.clone(); - pts.truncate(idx + 1); - // Make sure the last line isn't too tiny - if *pts.last().unwrap() == pt { - pts.pop(); - } - pts.push(pt); - if pts.len() == 1 { - return None; - } - Some(PolyLine::must_new(pts)) - } else { - None - } - } - - /// Same as get_slice_starting_at, but returns None if the point isn't on the polyline. - pub fn safe_get_slice_starting_at(&self, pt: Pt2D) -> Option { - if self.last_pt() == pt { - return None; - } - - if let Some(idx) = self.lines().position(|l| l.contains_pt(pt)) { - let mut pts = self.pts.clone(); - pts = pts.split_off(idx + 1); - if pt != pts[0] { - pts.insert(0, pt); - } - Some(PolyLine::must_new(pts)) - } else { - None - } - } - - pub fn dist_along_of_point(&self, pt: Pt2D) -> Option<(Distance, Angle)> { - let mut dist_along = Distance::ZERO; - for l in self.lines() { - if let Some(dist) = l.dist_along_of_point(pt) { - return Some((dist_along + dist, l.angle())); - } else { - dist_along += l.length(); - } - } - None - } - - pub fn trim_to_endpts(&self, pt1: Pt2D, pt2: Pt2D) -> PolyLine { - assert!(pt1 != pt2); - let mut dist1 = self.dist_along_of_point(pt1).unwrap().0; - let mut dist2 = self.dist_along_of_point(pt2).unwrap().0; - if dist1 > dist2 { - std::mem::swap(&mut dist1, &mut dist2); - } - self.exact_slice(dist1, dist2) - } - - pub fn get_bounds(&self) -> Bounds { - Bounds::from(&self.pts) - } - - /// If the current line is at least this long, return it. Otherwise, extend the end of it, - /// following the angle of the last line. - pub fn extend_to_length(&self, min_len: Distance) -> PolyLine { - let need_len = min_len - self.length(); - if need_len <= Distance::ZERO { - return self.clone(); - } - let line = self.last_line(); - // We might be extending a very tiny amount - if let Ok(extension) = PolyLine::new(vec![ - line.pt2(), - line.pt2().project_away(need_len, line.angle()), - ]) { - self.clone().must_extend(extension) - } else { - let mut pts = self.clone().into_points(); - pts.pop(); - pts.push(line.pt2().project_away(need_len, line.angle())); - PolyLine::must_new(pts) - } - } - - /// Produces a GeoJSON linestring, optionally mapping the world-space points back to GPS. - pub fn to_geojson(&self, gps: Option<&GPSBounds>) -> geojson::Geometry { - let mut pts = Vec::new(); - if let Some(gps) = gps { - for pt in gps.convert_back(&self.pts) { - pts.push(vec![pt.x(), pt.y()]); - } - } else { - for pt in &self.pts { - pts.push(vec![pt.x(), pt.y()]); - } - } - geojson::Geometry::new(geojson::Value::LineString(pts)) - } - - pub fn from_geojson(feature: &geojson::Feature, gps: Option<&GPSBounds>) -> Result { - if let Some(geojson::Geometry { - value: geojson::Value::LineString(ref pts), - .. - }) = feature.geometry - { - let mut points = Vec::new(); - for pt in pts { - let x = pt[0]; - let y = pt[1]; - if let Some(ref gps) = gps { - points.push(LonLat::new(x, y).to_pt(gps)); - } else { - points.push(Pt2D::new(x, y)); - } - } - PolyLine::new(points) - } else { - bail!("Input isn't a LineString") - } - } - - /// Returns the point on the polyline closest to the query. - pub fn project_pt(&self, query: Pt2D) -> Pt2D { - match self - .to_geo() - .closest_point(&geo::Point::new(query.x(), query.y())) - { - geo::Closest::Intersection(hit) | geo::Closest::SinglePoint(hit) => { - Pt2D::new(hit.x(), hit.y()) - } - geo::Closest::Indeterminate => unreachable!(), - } - } - - /// Returns the angle from the start to end of this polyline. - pub fn overall_angle(&self) -> Angle { - self.first_pt().angle_to(self.last_pt()) - } - - pub(crate) fn to_geo(&self) -> geo::LineString { - let pts: Vec = self - .pts - .iter() - .map(|pt| geo::Point::new(pt.x(), pt.y())) - .collect(); - pts.into() - } - - /// Walk along the PolyLine, starting `buffer_ends` from the start and ending `buffer_ends` - /// before the end. Advance in increments of `step_size`. Returns the point and angle at each - /// step. - pub fn step_along(&self, step_size: Distance, buffer_ends: Distance) -> Vec<(Pt2D, Angle)> { - self.step_along_start_end(step_size, buffer_ends, buffer_ends) - } - - /// Walk along the PolyLine, from `start_buffer` to `length - end_buffer`. Advance in - /// increments of `step_size`. Returns the point and angle at each step. - pub fn step_along_start_end( - &self, - step_size: Distance, - start_buffer: Distance, - end_buffer: Distance, - ) -> Vec<(Pt2D, Angle)> { - let mut result = Vec::new(); - let mut dist_along = start_buffer; - let length = self.length(); - while dist_along < length - end_buffer { - result.push(self.must_dist_along(dist_along)); - dist_along += step_size; - } - result - } - - /// - /// ``` - /// use geom::{PolyLine, Pt2D, Distance}; - /// - /// let polyline = PolyLine::must_new(vec![ - /// Pt2D::new(0.0, 0.0), - /// Pt2D::new(0.0, 10.0), - /// Pt2D::new(10.0, 20.0), - /// ]); - /// - /// assert_eq!( - /// polyline.interpolate_points(Distance::meters(20.0)).points(), - /// &vec![ - /// Pt2D::new(0.0, 0.0), - /// Pt2D::new(0.0, 10.0), - /// Pt2D::new(10.0, 20.0), - /// ] - /// ); - /// - /// assert_eq!( - /// polyline.interpolate_points(Distance::meters(10.0)).points(), - /// &vec![ - /// Pt2D::new(0.0, 0.0), - /// Pt2D::new(0.0, 10.0), - /// Pt2D::new(5.0, 15.0), - /// Pt2D::new(10.0, 20.0), - /// ] - /// ); - /// - /// ``` - pub fn interpolate_points(&self, max_step: Distance) -> PolyLine { - if self.pts.len() < 2 { - return self.clone(); - } - - let mut output = vec![]; - for line in self.lines() { - let points = (line.length() / max_step).ceil(); - let step_size = line.length() / points; - for i in 0..(points as usize) { - output.push(line.must_dist_along(step_size * i as f64)); - } - } - - output.push(*self.pts.last().unwrap()); - - PolyLine::new(output).unwrap() - } - - /// An arbitrary placeholder value, when Option types aren't worthwhile - pub fn dummy() -> PolyLine { - PolyLine::must_new(vec![Pt2D::new(0.0, 0.0), Pt2D::new(0.1, 0.1)]) - } - - /// Extracts all linestrings from raw bytes representing a GeoJSON file, along with the string - /// key/value properties. Uses `PolyLine::unchecked_new`, which doesn't dedupe or handle any - /// invalid cases. - pub fn from_geojson_bytes( - raw_bytes: &[u8], - gps_bounds: &GPSBounds, - require_in_bounds: bool, - ) -> Result)>> { - let raw_string = std::str::from_utf8(raw_bytes)?; - let geojson = raw_string.parse::()?; - let features = match geojson { - geojson::GeoJson::Feature(feature) => vec![feature], - geojson::GeoJson::FeatureCollection(collection) => collection.features, - _ => bail!("Unexpected geojson: {:?}", geojson), - }; - - let mut results = Vec::new(); - for feature in features { - if let Some(geom) = &feature.geometry { - let raw_pts = match &geom.value { - geojson::Value::LineString(pts) => pts, - // If there are multiple, just use the first - geojson::Value::MultiLineString(strings) => &strings[0], - _ => { - continue; - } - }; - let gps_pts: Vec = - raw_pts.iter().map(|pt| LonLat::new(pt[0], pt[1])).collect(); - let pts = if !require_in_bounds { - gps_bounds.convert(&gps_pts) - } else if let Some(pts) = gps_bounds.try_convert(&gps_pts) { - pts - } else { - continue; - }; - let mut tags = BTreeMap::new(); - for (key, value) in feature.properties_iter() { - if let Some(value) = value.as_str() { - tags.insert(key.to_string(), value.to_string()); - } - } - results.push((Self::unchecked_new(pts), tags)); - } - } - Ok(results) - } - - // Note this being false does not necessarily imply counter-clockwise; it might be neither - pub fn is_clockwise(&self) -> bool { - self.to_geo().is_cw() - } -} - -impl fmt::Display for PolyLine { - fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { - writeln!(f, "PolyLine::new(vec![ // length {}", self.length)?; - for (idx, pt) in self.pts.iter().enumerate() { - write!(f, " Pt2D::new({}, {}),", pt.x(), pt.y())?; - if idx > 0 { - let line = Line::must_new(self.pts[idx - 1], *pt); - write!( - f, - " // {}, {} (+ {} @ {})", - pt.x() - self.pts[idx - 1].x(), - pt.y() - self.pts[idx - 1].y(), - line.length(), - line.angle(), - )?; - } - writeln!(f)?; - } - write!(f, "])") - } -} - -fn fix_angles(orig: &PolyLine, result: PolyLine) -> Result { - let mut pts = result.pts.clone(); - - // Check that the angles roughly match up between the original and shifted line - for (idx, (orig_l, shifted_l)) in orig.lines().zip(result.lines()).enumerate() { - let orig_angle = orig_l.angle(); - let shifted_angle = shifted_l.angle(); - - if !orig_angle.approx_eq(shifted_angle, 1.0) { - // When this happens, the rotation is usually right around 180 -- so try swapping - // the points! - /*println!( - "Points changed angles from {} to {} (rot {})", - orig_angle, shifted_angle, rot - );*/ - pts.swap(idx, idx + 1); - // TODO Start the fixing over. but make sure we won't infinite loop... - //return fix_angles(orig, result); - } - } - - // When we swap points, length of the entire PolyLine may change! Recalculating is vital. - PolyLine::new(pts) -} - -// Also returns the duplicates. -fn to_set(pts: &[Pt2D]) -> (HashSet, HashSet) { - let mut deduped = HashSet::new(); - let mut dupes = HashSet::new(); - for pt in pts { - let pt = pt.to_hashable(); - if deduped.contains(&pt) { - dupes.insert(pt); - } else { - deduped.insert(pt); - } - } - (deduped, dupes) -} - -impl From<&PolyLine> for geo::LineString { - fn from(poly_line: &PolyLine) -> Self { - let coords: Vec = poly_line - .pts - .iter() - .map(|pt| geo::coord! { x: pt.x(), y: pt.y() }) - .collect(); - geo::LineString::new(coords) - } -} diff --git a/geom/src/pt.rs b/geom/src/pt.rs deleted file mode 100644 index 01ad36e8b8..0000000000 --- a/geom/src/pt.rs +++ /dev/null @@ -1,192 +0,0 @@ -use std::fmt; - -use geo::Simplify; -use ordered_float::NotNan; -use serde::{Deserialize, Serialize}; - -use crate::conversions::pts_to_line_string; -use crate::{ - deserialize_f64, serialize_f64, trim_f64, Angle, Distance, GPSBounds, LonLat, EPSILON_DIST, -}; - -/// This represents world-space in meters. -#[derive(Clone, Copy, Debug, Serialize, Deserialize)] -pub struct Pt2D { - #[serde(serialize_with = "serialize_f64", deserialize_with = "deserialize_f64")] - x: f64, - #[serde(serialize_with = "serialize_f64", deserialize_with = "deserialize_f64")] - y: f64, -} - -impl std::cmp::PartialEq for Pt2D { - fn eq(&self, other: &Pt2D) -> bool { - self.approx_eq(*other, EPSILON_DIST) - } -} - -impl Pt2D { - pub fn new(x: f64, y: f64) -> Pt2D { - if !x.is_finite() || !y.is_finite() { - panic!("Bad Pt2D {}, {}", x, y); - } - - // TODO enforce >=0 - - Pt2D { - x: trim_f64(x), - y: trim_f64(y), - } - } - - pub fn zero() -> Self { - Self::new(0.0, 0.0) - } - - // TODO This is a small first step... - pub fn approx_eq(self, other: Pt2D, threshold: Distance) -> bool { - self.dist_to(other) <= threshold - } - - /// Can go out of bounds. - pub fn to_gps(self, b: &GPSBounds) -> LonLat { - b.convert_back_xy(self.x(), self.y()) - } - - pub fn x(self) -> f64 { - self.x - } - - pub fn y(self) -> f64 { - self.y - } - - /// If distance is negative, this projects a point in theta.opposite() - pub fn project_away(self, dist: Distance, theta: Angle) -> Pt2D { - let (sin, cos) = theta.normalized_radians().sin_cos(); - Pt2D::new( - self.x() + dist.inner_meters() * cos, - self.y() + dist.inner_meters() * sin, - ) - } - - pub(crate) fn raw_dist_to(self, to: Pt2D) -> f64 { - ((self.x() - to.x()).powi(2) + (self.y() - to.y()).powi(2)).sqrt() - } - - pub fn dist_to(self, to: Pt2D) -> Distance { - Distance::meters(self.raw_dist_to(to)) - } - - /// Pretty meaningless units, for comparing distances very roughly - pub fn fast_dist(self, other: Pt2D) -> NotNan { - NotNan::new((self.x() - other.x()).powi(2) + (self.y() - other.y()).powi(2)).unwrap() - } - - pub fn angle_to(self, to: Pt2D) -> Angle { - // DON'T invert y here - Angle::new_rads((to.y() - self.y()).atan2(to.x() - self.x())) - } - - pub fn offset(self, dx: f64, dy: f64) -> Pt2D { - Pt2D::new(self.x() + dx, self.y() + dy) - } - - pub fn center(pts: &[Pt2D]) -> Pt2D { - if pts.is_empty() { - panic!("Can't find center of 0 points"); - } - let mut x = 0.0; - let mut y = 0.0; - for pt in pts { - x += pt.x(); - y += pt.y(); - } - let len = pts.len() as f64; - Pt2D::new(x / len, y / len) - } - - // Temporary until Pt2D has proper resolution. - pub fn approx_dedupe(pts: Vec, threshold: Distance) -> Vec { - // Just use dedup() on the Vec. - assert_ne!(threshold, EPSILON_DIST); - let mut result: Vec = Vec::new(); - for pt in pts { - if result.is_empty() || !result.last().unwrap().approx_eq(pt, threshold) { - result.push(pt); - } - } - result - } - - pub fn to_hashable(self) -> HashablePt2D { - HashablePt2D { - x_nan: NotNan::new(self.x()).unwrap(), - y_nan: NotNan::new(self.y()).unwrap(), - } - } - - /// Simplifies a list of points using Ramer-Douglas-Peuckr - pub fn simplify_rdp(pts: Vec, epsilon: f64) -> Vec { - let mut pts = pts_to_line_string(&pts) - .simplify(&epsilon) - .into_points() - .into_iter() - .map(|pt| pt.into()) - .collect::>(); - // TODO Not sure why, but from geo 0.23 to 0.24, this became necessary? - pts.dedup(); - pts - } - - pub fn to_geojson(self, gps: Option<&GPSBounds>) -> geojson::Geometry { - if let Some(gps) = gps { - self.to_gps(gps).to_geojson() - } else { - geojson::Geometry::new(geojson::Value::Point(vec![self.x(), self.y()])) - } - } -} - -impl fmt::Display for Pt2D { - fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { - write!(f, "Pt2D({0}, {1})", self.x(), self.y()) - } -} - -/// This represents world space, NOT LonLat. -// TODO So rename it HashablePair or something -#[derive(Clone, Copy, Debug, Eq, Hash, PartialEq, PartialOrd, Ord)] -pub struct HashablePt2D { - x_nan: NotNan, - y_nan: NotNan, -} - -impl HashablePt2D { - pub fn to_pt2d(self) -> Pt2D { - Pt2D::new(self.x_nan.into_inner(), self.y_nan.into_inner()) - } -} - -impl From for geo::Coord { - fn from(pt: Pt2D) -> Self { - geo::Coord { x: pt.x, y: pt.y } - } -} - -impl From for geo::Point { - fn from(pt: Pt2D) -> Self { - geo::Point::new(pt.x, pt.y) - } -} - -impl From for Pt2D { - fn from(coord: geo::Coord) -> Self { - Pt2D::new(coord.x, coord.y) - } -} - -impl From for Pt2D { - fn from(point: geo::Point) -> Self { - Pt2D::new(point.x(), point.y()) - } -} diff --git a/geom/src/ring.rs b/geom/src/ring.rs deleted file mode 100644 index 6ecc0be3a1..0000000000 --- a/geom/src/ring.rs +++ /dev/null @@ -1,331 +0,0 @@ -use std::collections::HashSet; -use std::fmt; -use std::fmt::Write; - -use anyhow::Result; -use serde::{Deserialize, Serialize}; - -use crate::{Distance, GPSBounds, Line, PolyLine, Polygon, Pt2D, Tessellation}; - -/// Maybe a misnomer, but like a PolyLine, but closed. -#[derive(Clone, Debug, PartialEq, Serialize, Deserialize)] -pub struct Ring { - // first equals last - pts: Vec, -} - -impl Ring { - pub fn new(pts: Vec) -> Result { - if pts.len() < 3 { - bail!("Can't make a ring with < 3 points"); - } - if pts[0] != *pts.last().unwrap() { - bail!("Can't make a ring with mismatching first/last points"); - } - - if let Some(pair) = pts.windows(2).find(|pair| pair[0] == pair[1]) { - bail!("Ring has duplicate adjacent points near {}", pair[0]); - } - - let result = Ring { pts }; - - let mut seen_pts = HashSet::new(); - for pt in result.pts.iter().skip(1) { - if seen_pts.contains(&pt.to_hashable()) { - bail!("Ring has repeat non-adjacent points near {}", pt); - } - seen_pts.insert(pt.to_hashable()); - } - - Ok(result) - } - - /// Use with caution. Ignores duplicate points - pub fn unsafe_deduping_new(mut pts: Vec) -> Result { - pts.dedup(); - if pts.len() < 3 { - bail!("Can't make a ring with < 3 points"); - } - if pts[0] != *pts.last().unwrap() { - bail!("Can't make a ring with mismatching first/last points"); - } - - if let Some(pair) = pts.windows(2).find(|pair| pair[0] == pair[1]) { - bail!("Ring has duplicate adjacent points near {}", pair[0]); - } - - let result = Ring { pts }; - - let mut seen_pts = HashSet::new(); - for pt in result.pts.iter().skip(1) { - if seen_pts.contains(&pt.to_hashable()) { - // Just spam logs instead of crashing - println!("Ring has repeat non-adjacent points near {}", pt); - } - seen_pts.insert(pt.to_hashable()); - } - - Ok(result) - } - - pub fn must_new(pts: Vec) -> Ring { - Ring::new(pts).unwrap() - } - - /// First dedupes adjacent points - pub fn deduping_new(mut pts: Vec) -> Result { - pts.dedup(); - Self::new(pts) - } - - pub fn as_polyline(&self) -> PolyLine { - PolyLine::unchecked_new(self.pts.clone()) - } - - /// Draws the ring with some thickness, with half of it straddling the interor of the ring, and - /// half on the outside. - pub fn to_outline(&self, thickness: Distance) -> Tessellation { - // TODO Has a weird corner. Use the polygon offset thing instead? - self.as_polyline().thicken_tessellation(thickness) - } - - pub fn into_polygon(self) -> Polygon { - Polygon::with_holes(self, Vec::new()) - } - - pub fn points(&self) -> &Vec { - &self.pts - } - pub fn into_points(self) -> Vec { - self.pts - } - - /// Be careful with the order of results. Hits on an earlier line segment of other show up - /// first, but if the ring hits a line segment at multiple points, who knows. Dedupes. - pub fn all_intersections(&self, other: &PolyLine) -> Vec { - let mut hits = Vec::new(); - let mut seen = HashSet::new(); - for l1 in other.lines() { - for l2 in self - .pts - .windows(2) - .map(|pair| Line::must_new(pair[0], pair[1])) - { - if let Some(pt) = l1.intersection(&l2) { - if !seen.contains(&pt.to_hashable()) { - hits.push(pt); - seen.insert(pt.to_hashable()); - } - } - } - } - hits - } - - pub(crate) fn get_both_slices_btwn( - &self, - pt1: Pt2D, - pt2: Pt2D, - ) -> Option<(PolyLine, PolyLine)> { - assert!(pt1 != pt2); - let pl = self.as_polyline(); - - let mut dist1 = pl.dist_along_of_point(pt1)?.0; - let mut dist2 = pl.dist_along_of_point(pt2)?.0; - if dist1 > dist2 { - std::mem::swap(&mut dist1, &mut dist2); - } - if dist1 == dist2 { - return None; - } - - // TODO If we reversed the points, we need to reverse these results! Argh - let candidate1 = pl.maybe_exact_slice(dist1, dist2).ok()?; - let candidate2 = pl - .maybe_exact_slice(dist2, pl.length()) - .ok()? - .must_extend(pl.maybe_exact_slice(Distance::ZERO, dist1).ok()?); - Some((candidate1, candidate2)) - } - - /// Assuming both points are somewhere along the ring, return the points in between the two, by - /// tracing along the ring in the longer or shorter direction (depending on `longer`). If both - /// points are the same, returns `None`. The result is oriented from `pt1` to `pt2`. - pub fn get_slice_between(&self, pt1: Pt2D, pt2: Pt2D, longer: bool) -> Option { - if pt1 == pt2 { - return None; - } - let (candidate1, candidate2) = self.get_both_slices_btwn(pt1, pt2)?; - let slice = if longer == (candidate1.length() > candidate2.length()) { - candidate1 - } else { - candidate2 - }; - if slice.first_pt() == pt1 { - Some(slice) - } else { - // TODO Do we want to be paranoid here? Or just do the fix in get_both_slices_btwn - // directly? - Some(slice.reversed()) - } - } - - /// Assuming both points are somewhere along the ring, return the points in between the two, by - /// tracing along the ring in the shorter direction. If both points are the same, returns - /// `None`. The result is oriented from `pt1` to `pt2`. - pub fn get_shorter_slice_between(&self, pt1: Pt2D, pt2: Pt2D) -> Option { - self.get_slice_between(pt1, pt2, false) - } - - // TODO Rmove this one, fix all callers - pub fn get_shorter_slice_btwn(&self, pt1: Pt2D, pt2: Pt2D) -> Option { - let (candidate1, candidate2) = self.get_both_slices_btwn(pt1, pt2)?; - if candidate1.length() < candidate2.length() { - Some(candidate1) - } else { - Some(candidate2) - } - } - - /// Extract all PolyLines and Rings. Doesn't handle crazy double loops and stuff. - pub fn split_points(pts: &[Pt2D]) -> Result<(Vec, Vec)> { - let mut seen = HashSet::new(); - let mut intersections = HashSet::new(); - for pt in pts { - let pt = pt.to_hashable(); - if seen.contains(&pt) { - intersections.insert(pt); - } else { - seen.insert(pt); - } - } - intersections.insert(pts[0].to_hashable()); - intersections.insert(pts.last().unwrap().to_hashable()); - - let mut polylines = Vec::new(); - let mut rings = Vec::new(); - let mut current = Vec::new(); - for pt in pts.iter().cloned() { - current.push(pt); - if intersections.contains(&pt.to_hashable()) && current.len() > 1 { - if current[0] == pt && current.len() >= 3 { - rings.push(Ring::new(current.drain(..).collect())?); - } else { - polylines.push(PolyLine::new(current.drain(..).collect())?); - } - current.push(pt); - } - } - Ok((polylines, rings)) - } - - pub fn contains_pt(&self, pt: Pt2D) -> bool { - self.as_polyline().dist_along_of_point(pt).is_some() - } - - /// Produces a GeoJSON polygon, optionally mapping the world-space points back to GPS. - pub fn to_geojson(&self, gps: Option<&GPSBounds>) -> geojson::Geometry { - let mut pts = Vec::new(); - if let Some(gps) = gps { - for pt in gps.convert_back(&self.pts) { - pts.push(vec![pt.x(), pt.y()]); - } - } else { - for pt in &self.pts { - pts.push(vec![pt.x(), pt.y()]); - } - } - geojson::Geometry::new(geojson::Value::Polygon(vec![pts])) - } - - /// Translates the ring by a fixed offset. - pub fn translate(mut self, dx: f64, dy: f64) -> Ring { - for pt in &mut self.pts { - *pt = pt.offset(dx, dy); - } - self - } - - /// Find the "pole of inaccessibility" -- the most distant internal point from the polygon - /// outline - pub fn polylabel(&self) -> Pt2D { - // TODO Refactor to_geo? - let polygon = geo::Polygon::new( - geo::LineString::from( - self.pts - .iter() - .map(|pt| geo::Point::new(pt.x(), pt.y())) - .collect::>(), - ), - Vec::new(), - ); - let pt = polylabel::polylabel(&polygon, &1.0).unwrap(); - Pt2D::new(pt.x(), pt.y()) - } - - /// Look for "bad" rings that double back on themselves. These're likely to cause many - /// downstream problems. "Bad" means the order of points doesn't match the order when sorting - /// by angle from the center. - /// - /// TODO I spot many false positives. Look for better definitions -- maybe self-intersecting - /// polygons? - pub fn doubles_back(&self) -> bool { - // Skip the first=last point - let mut orig = self.pts.clone(); - orig.pop(); - // Polylabel is better than center; sometimes the center is very close to an edge - let center = self.polylabel(); - - let mut sorted = orig.clone(); - sorted.sort_by_key(|pt| pt.angle_to(center).normalized_degrees() as i64); - - // Align things again - while sorted[0] != orig[0] { - sorted.rotate_right(1); - } - - // Do they match up? - orig != sorted - } - - /// Print the coordinates of this ring as a `geo::LineString` for easy bug reports - pub fn as_geo_linestring(&self) -> String { - let mut output = String::new(); - writeln!(output, "let line_string = geo_types::line_string![").unwrap(); - for pt in &self.pts { - writeln!(output, " (x: {}, y: {}),", pt.x(), pt.y()).unwrap(); - } - writeln!(output, "];").unwrap(); - output - } -} - -impl fmt::Display for Ring { - fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result { - writeln!(f, "Ring::new(vec![")?; - for pt in &self.pts { - writeln!(f, " Pt2D::new({}, {}),", pt.x(), pt.y())?; - } - write!(f, "])") - } -} - -impl From for geo::LineString { - fn from(ring: Ring) -> Self { - let coords = ring - .pts - .into_iter() - .map(geo::Coord::from) - .collect::>(); - Self(coords) - } -} - -impl TryFrom for Ring { - type Error = anyhow::Error; - - fn try_from(line_string: geo::LineString) -> Result { - let pts: Vec = line_string.0.into_iter().map(Pt2D::from).collect(); - Self::deduping_new(pts) - } -} diff --git a/geom/src/speed.rs b/geom/src/speed.rs deleted file mode 100644 index 80af829808..0000000000 --- a/geom/src/speed.rs +++ /dev/null @@ -1,123 +0,0 @@ -use std::{cmp, ops}; - -use serde::{Deserialize, Serialize}; - -use crate::{deserialize_f64, serialize_f64, trim_f64, Distance, Duration, UnitFmt}; - -/// In meters per second. Can be negative. -#[derive(Clone, Copy, Debug, PartialEq, PartialOrd, Serialize, Deserialize)] -pub struct Speed( - #[serde(serialize_with = "serialize_f64", deserialize_with = "deserialize_f64")] f64, -); - -// By construction, Speed is a finite f64 with trimmed precision. -impl Eq for Speed {} - -#[allow(clippy::derive_ord_xor_partial_ord)] // false positive -impl Ord for Speed { - fn cmp(&self, other: &Speed) -> cmp::Ordering { - self.partial_cmp(other).unwrap() - } -} - -impl Speed { - pub const ZERO: Speed = Speed::const_meters_per_second(0.0); - - pub fn meters_per_second(value: f64) -> Speed { - if !value.is_finite() { - panic!("Bad Speed {}", value); - } - - Speed(trim_f64(value)) - } - - pub const fn const_meters_per_second(value: f64) -> Speed { - Speed(value) - } - - pub fn miles_per_hour(value: f64) -> Speed { - Speed::meters_per_second(0.44704 * value) - } - - pub fn km_per_hour(value: f64) -> Speed { - Speed::meters_per_second(0.277778 * value) - } - - pub fn from_dist_time(d: Distance, t: Duration) -> Speed { - Speed::meters_per_second(d.inner_meters() / t.inner_seconds()) - } - - // TODO Remove if possible. - pub fn inner_meters_per_second(self) -> f64 { - self.0 - } - - pub fn to_miles_per_hour(self) -> f64 { - self.0 * 2.23694 - } - - /// Describes the speed according to formatting rules. - pub fn to_string(self, fmt: &UnitFmt) -> String { - if fmt.metric { - format!("{} km/h", (self.0 * 3.6).round()) - } else { - format!("{} mph", self.to_miles_per_hour().round()) - } - } -} - -impl ops::Add for Speed { - type Output = Speed; - - fn add(self, other: Speed) -> Speed { - Speed::meters_per_second(self.0 + other.0) - } -} - -impl ops::Sub for Speed { - type Output = Speed; - - fn sub(self, other: Speed) -> Speed { - Speed::meters_per_second(self.0 - other.0) - } -} - -impl ops::Div for Speed { - type Output = f64; - - fn div(self, other: Speed) -> f64 { - self.0 / other.0 - } -} - -impl ops::Neg for Speed { - type Output = Speed; - - fn neg(self) -> Speed { - Speed::meters_per_second(-self.0) - } -} - -impl ops::Mul for Speed { - type Output = Speed; - - fn mul(self, scalar: f64) -> Speed { - Speed::meters_per_second(self.0 * scalar) - } -} - -impl ops::Mul for f64 { - type Output = Speed; - - fn mul(self, other: Speed) -> Speed { - Speed::meters_per_second(self * other.0) - } -} - -impl ops::Mul for Speed { - type Output = Distance; - - fn mul(self, other: Duration) -> Distance { - Distance::meters(self.0 * other.inner_seconds()) - } -} diff --git a/geom/src/stats.rs b/geom/src/stats.rs deleted file mode 100644 index 394625b0b7..0000000000 --- a/geom/src/stats.rs +++ /dev/null @@ -1,198 +0,0 @@ -use serde::{Deserialize, Serialize}; - -use crate::{Distance, Duration}; - -#[derive(Clone, Copy, PartialEq, Eq, PartialOrd, Ord, Serialize, Deserialize)] -pub enum Statistic { - Min, - Mean, - P50, - P90, - P99, - Max, -} - -impl Statistic { - pub fn all() -> Vec { - vec![ - Statistic::Min, - Statistic::Mean, - Statistic::P50, - Statistic::P90, - Statistic::P99, - Statistic::Max, - ] - } -} - -impl std::fmt::Display for Statistic { - fn fmt(&self, f: &mut std::fmt::Formatter) -> std::fmt::Result { - match self { - Statistic::Min => write!(f, "minimum"), - Statistic::Mean => write!(f, "mean"), - Statistic::P50 => write!(f, "50%ile"), - Statistic::P90 => write!(f, "90%ile"), - Statistic::P99 => write!(f, "99%ile"), - Statistic::Max => write!(f, "maximum"), - } - } -} - -pub trait HgramValue: Copy + std::cmp::Ord + std::fmt::Display { - // TODO Weird name because I can't figure out associated type mess in FanChart - fn hgram_zero() -> T; - fn to_u64(self) -> u64; - fn from_u64(x: u64) -> T; -} - -impl HgramValue for Duration { - fn hgram_zero() -> Duration { - Duration::ZERO - } - fn to_u64(self) -> u64 { - self.to_u64() - } - fn from_u64(x: u64) -> Duration { - Duration::from_u64(x) - } -} - -impl HgramValue for Distance { - fn hgram_zero() -> Distance { - Distance::ZERO - } - fn to_u64(self) -> u64 { - self.to_u64() - } - fn from_u64(x: u64) -> Distance { - Distance::from_u64(x) - } -} - -impl HgramValue for u16 { - fn hgram_zero() -> u16 { - 0 - } - fn to_u64(self) -> u64 { - self as u64 - } - fn from_u64(x: u64) -> u16 { - u16::try_from(x).unwrap() - } -} - -impl HgramValue for usize { - fn hgram_zero() -> usize { - 0 - } - fn to_u64(self) -> u64 { - self as u64 - } - fn from_u64(x: u64) -> usize { - x as usize - } -} - -// TODO Maybe consider having a type-safe NonEmptyHistogram. -#[derive(Clone)] -pub struct Histogram> { - count: usize, - histogram: histogram::Histogram, - min: T, - max: T, -} - -impl> Default for Histogram { - fn default() -> Histogram { - Histogram { - count: 0, - histogram: Default::default(), - min: T::hgram_zero(), - max: T::hgram_zero(), - } - } -} - -impl> Histogram { - pub fn new() -> Histogram { - Default::default() - } - - pub fn add(&mut self, x: T) { - if self.count == 0 { - self.min = x; - self.max = x; - } else { - self.min = self.min.min(x); - self.max = self.max.max(x); - } - self.count += 1; - self.histogram - .increment(x.to_u64()) - .map_err(|err| format!("Can't add {}: {}", x, err)) - .unwrap(); - } - - pub fn remove(&mut self, x: T) { - // TODO This doesn't update min/max! Why are we tracking that ourselves? Do we not trust - // the lossiness of the underlying histogram? - self.count -= 1; - self.histogram - .decrement(x.to_u64()) - .map_err(|err| format!("Can't remove {}: {}", x, err)) - .unwrap(); - } - - pub fn describe(&self) -> String { - if self.count == 0 { - return "no data yet".to_string(); - } - - format!( - "{} count, 50%ile {}, 90%ile {}, 99%ile {}, min {}, mean {}, max {}", - crate::utils::prettyprint_usize(self.count), - self.select(Statistic::P50).unwrap(), - self.select(Statistic::P90).unwrap(), - self.select(Statistic::P99).unwrap(), - self.select(Statistic::Min).unwrap(), - self.select(Statistic::Mean).unwrap(), - self.select(Statistic::Max).unwrap(), - ) - } - - /// None if empty - pub fn percentile(&self, p: f64) -> Option { - if self.count == 0 { - return None; - } - Some(T::from_u64(self.histogram.percentile(p).unwrap())) - } - - pub fn select(&self, stat: Statistic) -> Option { - if self.count == 0 { - return None; - } - let raw = match stat { - Statistic::P50 => self.histogram.percentile(50.0).unwrap(), - Statistic::P90 => self.histogram.percentile(90.0).unwrap(), - Statistic::P99 => self.histogram.percentile(99.0).unwrap(), - Statistic::Min => { - return Some(self.min); - } - Statistic::Mean => self.histogram.mean().unwrap(), - Statistic::Max => { - return Some(self.max); - } - }; - Some(T::from_u64(raw)) - } - - pub fn count(&self) -> usize { - self.count - } - - /// Could implement PartialEq, but be a bit more clear how approximate this is - pub fn seems_eq(&self, other: &Histogram) -> bool { - self.describe() == other.describe() - } -} diff --git a/geom/src/tessellation.rs b/geom/src/tessellation.rs deleted file mode 100644 index fbd4ed13b2..0000000000 --- a/geom/src/tessellation.rs +++ /dev/null @@ -1,264 +0,0 @@ -use anyhow::Result; -use serde::{Deserialize, Serialize}; - -use crate::{Angle, Bounds, GPSBounds, Polygon, Pt2D}; - -// Only serializable for Polygons that precompute a tessellation -/// A tessellated polygon, ready for rendering. -#[derive(PartialEq, Serialize, Deserialize, Clone, Debug)] -pub struct Tessellation { - /// These points aren't in any meaningful order. It's not generally possible to reconstruct a - /// `Polygon` from this. - pub(crate) points: Vec, - /// Groups of three indices make up the triangles - indices: Vec, -} - -#[derive(Clone, Debug)] -pub struct Triangle { - pub pt1: Pt2D, - pub pt2: Pt2D, - pub pt3: Pt2D, -} - -impl From for Tessellation { - fn from(mut polygon: Polygon) -> Self { - if let Some(tessellation) = polygon.tessellation.take() { - return tessellation; - } - - let geojson_style: Vec>> = polygon - .rings - .iter() - .map(|ring| { - ring.points() - .iter() - .map(|pt| vec![pt.x(), pt.y()]) - .collect() - }) - .collect(); - let (vertices, holes, dims) = earcutr::flatten(&geojson_style); - let indices = downsize(earcutr::earcut(&vertices, &holes, dims).unwrap()); - - Self { - points: vertices - .chunks(2) - .map(|pair| Pt2D::new(pair[0], pair[1])) - .collect(), - indices, - } - } -} - -impl From for Tessellation { - fn from(poly: geo::Polygon) -> Self { - // geo::Polygon -> geom::Polygon may fail, so we can't just do two hops. We can tessellate - // something even if it has invalid Rings. - let (exterior, mut interiors) = poly.into_inner(); - interiors.insert(0, exterior); - - let geojson_style: Vec>> = interiors - .into_iter() - .map(|ring| { - ring.into_inner() - .into_iter() - .map(|pt| vec![pt.x, pt.y]) - .collect() - }) - .collect(); - let (vertices, holes, dims) = earcutr::flatten(&geojson_style); - let indices = earcutr::earcut(&vertices, &holes, dims).unwrap(); - - let points = vertices - .chunks(2) - .map(|pair| Pt2D::new(pair[0], pair[1])) - .collect(); - - Self::new(points, indices) - } -} - -impl Tessellation { - pub fn new(points: Vec, indices: Vec) -> Self { - Tessellation { - points, - indices: downsize(indices), - } - } - - /// The `points` are not necessarily a `Ring`, which has strict requirements about no duplicate - /// points. We can render various types of invalid polygon. - pub fn from_ring(points: Vec) -> Self { - assert!(points.len() >= 3); - - let mut vertices = Vec::new(); - for pt in &points { - vertices.push(pt.x()); - vertices.push(pt.y()); - } - let indices = downsize(earcutr::earcut(&vertices, &Vec::new(), 2).unwrap()); - - Self { points, indices } - } - - /// Returns (points, indices) for rendering - pub fn consume(self) -> (Vec, Vec) { - (self.points, self.indices) - } - - pub fn triangles(&self) -> Vec { - let mut triangles: Vec = Vec::new(); - for slice in self.indices.chunks_exact(3) { - triangles.push(Triangle { - pt1: self.points[slice[0] as usize], - pt2: self.points[slice[1] as usize], - pt3: self.points[slice[2] as usize], - }); - } - triangles - } - - pub fn get_bounds(&self) -> Bounds { - Bounds::from(&self.points) - } - - pub fn center(&self) -> Pt2D { - self.get_bounds().center() - } - - pub(crate) fn transform Pt2D>(&mut self, f: F) { - for pt in &mut self.points { - *pt = f(pt); - } - } - - pub fn translate(&mut self, dx: f64, dy: f64) { - self.transform(|pt| pt.offset(dx, dy)); - } - - pub fn scale(&mut self, factor: f64) { - self.transform(|pt| Pt2D::new(pt.x() * factor, pt.y() * factor)); - } - - pub fn scale_xy(&mut self, x_factor: f64, y_factor: f64) { - self.transform(|pt| Pt2D::new(pt.x() * x_factor, pt.y() * y_factor)) - } - - pub fn rotate(&mut self, angle: Angle) { - self.rotate_around(angle, self.center()) - } - - pub fn rotate_around(&mut self, angle: Angle, pivot: Pt2D) { - self.transform(|pt| { - let origin_pt = Pt2D::new(pt.x() - pivot.x(), pt.y() - pivot.y()); - let (sin, cos) = angle.normalized_radians().sin_cos(); - Pt2D::new( - pivot.x() + origin_pt.x() * cos - origin_pt.y() * sin, - pivot.y() + origin_pt.y() * cos + origin_pt.x() * sin, - ) - }) - } - - /// Equivalent to `self.scale(scale).translate(translate_x, translate_y).rotate_around(rotate, - /// pivot)`, but modifies the polygon in-place and is faster. - pub fn inplace_multi_transform( - &mut self, - scale: f64, - translate_x: f64, - translate_y: f64, - rotate: Angle, - pivot: Pt2D, - ) { - let (sin, cos) = rotate.normalized_radians().sin_cos(); - - for pt in &mut self.points { - // Scale - let x = scale * pt.x(); - let y = scale * pt.y(); - // Translate - let x = x + translate_x; - let y = y + translate_y; - // Rotate - let origin_pt = Pt2D::new(x - pivot.x(), y - pivot.y()); - *pt = Pt2D::new( - pivot.x() + origin_pt.x() * cos - origin_pt.y() * sin, - pivot.y() + origin_pt.y() * cos + origin_pt.x() * sin, - ); - } - } - - pub fn union(self, other: Self) -> Self { - let mut points = self.points; - let mut indices = self.indices; - let offset = points.len() as u16; - points.extend(other.points); - for idx in other.indices { - indices.push(offset + idx); - } - Self { points, indices } - } - - pub fn union_all(mut list: Vec) -> Self { - let mut result = list.pop().unwrap(); - for p in list { - result = result.union(p); - } - result - } - - /// Produces a GeoJSON multipolygon consisting of individual triangles. Optionally map the - /// world-space points back to GPS. - pub fn to_geojson(&self, gps: Option<&GPSBounds>) -> geojson::Geometry { - let mut polygons = Vec::new(); - for triangle in self.triangles() { - let raw_pts = vec![triangle.pt1, triangle.pt2, triangle.pt3, triangle.pt1]; - let mut pts = Vec::new(); - if let Some(gps) = gps { - for pt in gps.convert_back(&raw_pts) { - pts.push(vec![pt.x(), pt.y()]); - } - } else { - for pt in raw_pts { - pts.push(vec![pt.x(), pt.y()]); - } - } - polygons.push(vec![pts]); - } - - geojson::Geometry::new(geojson::Value::MultiPolygon(polygons)) - } - - // TODO This only makes sense for something vaguely Ring-like - fn to_geo(&self) -> geo::Polygon { - let exterior = crate::conversions::pts_to_line_string(&self.points); - geo::Polygon::new(exterior, Vec::new()) - } - - // TODO After making to_outline return a real Polygon, get rid of this - pub fn difference(&self, other: &Tessellation) -> Result> { - use geo::BooleanOps; - - // TODO Remove after https://github.com/georust/geo/issues/913 - match std::panic::catch_unwind(std::panic::AssertUnwindSafe(|| { - crate::polygon::from_multi(self.to_geo().difference(&other.to_geo())) - })) { - Ok(result) => result, - Err(err) => { - println!("BooleanOps crashed: {err:?}"); - bail!("BooleanOps crashed: {err:?}"); - } - } - } -} - -fn downsize(input: Vec) -> Vec { - let mut output = Vec::new(); - for x in input { - if let Ok(x) = u16::try_from(x) { - output.push(x); - } else { - panic!("{} can't fit in u16, some polygon is too huge", x); - } - } - output -} diff --git a/geom/src/time.rs b/geom/src/time.rs deleted file mode 100644 index 54d947aeb2..0000000000 --- a/geom/src/time.rs +++ /dev/null @@ -1,226 +0,0 @@ -use std::{cmp, ops}; - -use anyhow::Result; -use ordered_float::NotNan; -use serde::{Deserialize, Serialize}; - -use crate::{deserialize_f64, serialize_f64, trim_f64, Duration}; - -/// In seconds since midnight. Can't be negative. -#[derive(Clone, Copy, Debug, Default, PartialEq, PartialOrd, Serialize, Deserialize)] -pub struct Time( - #[serde(serialize_with = "serialize_f64", deserialize_with = "deserialize_f64")] f64, -); - -// By construction, Time is a finite f64 with trimmed precision. -impl Eq for Time {} - -#[allow(clippy::derive_ord_xor_partial_ord)] // false positive -impl Ord for Time { - fn cmp(&self, other: &Time) -> cmp::Ordering { - self.partial_cmp(other).unwrap() - } -} - -#[allow(clippy::derive_hash_xor_eq)] // false positive -impl std::hash::Hash for Time { - fn hash(&self, state: &mut H) { - NotNan::new(self.0).unwrap().hash(state); - } -} - -impl Time { - pub const START_OF_DAY: Time = Time(0.0); - - // No direct public constructors. Explicitly do Time::START_OF_DAY + duration. - fn seconds_since_midnight(value: f64) -> Time { - if !value.is_finite() || value < 0.0 { - panic!("Bad Time {}", value); - } - - Time(trim_f64(value)) - } - - /// (hours, minutes, seconds, centiseconds) - fn get_parts(self) -> (usize, usize, usize, usize) { - let mut remainder = self.0; - let hours = (remainder / 3600.0).floor(); - remainder -= hours * 3600.0; - let minutes = (remainder / 60.0).floor(); - remainder -= minutes * 60.0; - let seconds = remainder.floor(); - remainder -= seconds; - let centis = (remainder / 0.1).floor(); - - ( - hours as usize, - minutes as usize, - seconds as usize, - centis as usize, - ) - } - /// Rounded down. 6:59:00 is hour 6. - pub fn get_hours(self) -> usize { - self.get_parts().0 - } - - pub fn ampm_tostring(self) -> String { - let (mut hours, minutes, seconds, _) = self.get_parts(); - let next_day = if hours >= 24 { - let days = hours / 24; - hours %= 24; - format!(" (+{} days)", days) - } else { - "".to_string() - }; - let suffix = if hours < 12 { "AM" } else { "PM" }; - if hours == 0 { - hours = 12; - } else if hours >= 24 { - hours -= 24; - } else if hours > 12 { - hours -= 12; - } - - format!( - "{:02}:{:02}:{:02} {}{}", - hours, minutes, seconds, suffix, next_day - ) - } - - pub fn as_filename(self) -> String { - let (hours, minutes, seconds, remainder) = self.get_parts(); - format!( - "{0:02}h{1:02}m{2:02}.{3:01}s", - hours, minutes, seconds, remainder - ) - } - - pub fn parse(string: &str) -> Result