Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Unify line measures #1216

Merged
merged 11 commits into from
Oct 8, 2024
Merged
Show file tree
Hide file tree
Changes from 6 commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 3 additions & 3 deletions .github/workflows/test.yml
Original file line number Diff line number Diff line change
Expand Up @@ -77,7 +77,7 @@ jobs:
# giving us about 6 months of coverage.
#
# Minimum supported rust version (MSRV)
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.74"
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.75"
# Two most recent releases - we omit older ones for expedient CI
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.78"
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.79"
Expand Down Expand Up @@ -105,7 +105,7 @@ jobs:
# giving us about 6 months of coverage.
#
# Minimum supported rust version (MSRV)
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.74"
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.75"
# Two most recent releases - we omit older ones for expedient CI
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.78"
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.79"
Expand All @@ -132,7 +132,7 @@ jobs:
# giving us about 6 months of coverage.
#
# Minimum supported rust version (MSRV)
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.74"
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.75"
# Two most recent releases - we omit older ones for expedient CI
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.78"
- "ghcr.io/georust/geo-ci:proj-9.4.0-rust-1.79"
Expand Down
2 changes: 1 addition & 1 deletion geo-postgis/Cargo.toml
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@ documentation = "https://docs.rs/geo-postgis/"
readme = "../README.md"
keywords = ["gis", "geo", "geography", "geospatial", "postgis"]
description = "Conversion between `geo-types` and `postgis` types."
rust-version = "1.74"
rust-version = "1.75"
edition = "2021"

[dependencies]
Expand Down
2 changes: 1 addition & 1 deletion geo-types/Cargo.toml
Original file line number Diff line number Diff line change
Expand Up @@ -7,7 +7,7 @@ documentation = "https://docs.rs/geo-types/"
readme = "../README.md"
keywords = ["gis", "geo", "geography", "geospatial"]
description = "Geospatial primitive data types"
rust-version = "1.74"
rust-version = "1.75"
edition = "2021"

[features]
Expand Down
20 changes: 20 additions & 0 deletions geo/CHANGES.md
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,26 @@
* <https://github.com/georust/geo/pull/1201>
* Add `StitchTriangles` trait which implements a new kind of combining algorithm for `Triangle`s
* <https://github.com/georust/geo/pull/1087>
* BREAKING: Remove deprecated `Bearing` trait
* Unify various line measurements under new `line_measures::{Bearing, Distance, Destination, InterpolatePoint}` traits
Before:
```
use geo::{GeodesicBearing, HaversineBearing, GeodesicDistance, HaversineDistance};
GeodesicBearing::geodesic_bearing(p1, p2)
HaversineBearing::haversine_bearing(p1, p2)
GeodesicDistance::geodesic_distance(p1, p2)
HaversineDistance::haversine_distance(p1, p2)
```

After:
```
use geo::{Geodesic, Haversine, Bearing, Distance};
Geodesic::bearing(p1, p2)
Haversine::bearing(p1, p2)
Geodesic::distance(p1, p2)
Haversine::distance(p1, p2)
```
* <https://github.com/georust/geo/pull/1216>

## 0.28.0

Expand Down
2 changes: 1 addition & 1 deletion geo/Cargo.toml
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@ readme = "../README.md"
keywords = ["gis", "geo", "geography", "geospatial"]
autobenches = true
edition = "2021"
rust-version = "1.74"
rust-version = "1.75"
categories = ["science::geo"]

[features]
Expand Down
33 changes: 0 additions & 33 deletions geo/src/algorithm/bearing.rs

This file was deleted.

2 changes: 0 additions & 2 deletions geo/src/algorithm/geodesic_intermediate.rs
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,6 @@ use crate::{CoordFloat, Point};
use geographiclib_rs::{DirectGeodesic, Geodesic, InverseGeodesic};

/// Returns a new Point along a route between two existing points on an ellipsoidal model of the earth

pub trait GeodesicIntermediate<T: CoordFloat> {
/// Returns a new Point along a route between two existing points on an ellipsoidal model of the earth
///
Expand All @@ -25,7 +24,6 @@ pub trait GeodesicIntermediate<T: CoordFloat> {
/// assert_relative_eq!(i50, i50_should, epsilon = 1.0e-6);
/// assert_relative_eq!(i80, i80_should, epsilon = 1.0e-6);
/// ```

fn geodesic_intermediate(&self, other: &Point<T>, f: T) -> Point<T>;
fn geodesic_intermediate_fill(
&self,
Expand Down
13 changes: 13 additions & 0 deletions geo/src/algorithm/line_measures/bearing.rs
Original file line number Diff line number Diff line change
@@ -0,0 +1,13 @@
use geo_types::{CoordFloat, Point};

/// Calculate the bearing between two points
pub trait Bearing<F: CoordFloat> {
/// Calculate the bearing from `origin` to `destination` in degrees.
///
/// See [specific implementations](#implementors) for details.
///
/// # Units
/// - `origin`, `destination`: Point where the units of x/y depend on the [trait implementation](#implementors).
/// - returns: degrees, where: North: 0°, East: 90°, South: 180°, West: 270°
fn bearing(origin: Point<F>, destination: Point<F>) -> F;
}
19 changes: 19 additions & 0 deletions geo/src/algorithm/line_measures/destination.rs
Original file line number Diff line number Diff line change
@@ -0,0 +1,19 @@
use geo_types::{CoordFloat, Point};

/// Calculate the destination point from an origin point, a bearing and a distance.
pub trait Destination<F: CoordFloat> {
/// Returns a new point having travelled the `distance` along a line
/// from the `origin` point with the given `bearing`.
///
/// See [specific implementations](#implementors) for details.
///
/// # Units
///
/// - `origin`: Point where the units of x/y depend on the [trait implementation](#implementors).
/// - `bearing`: degrees, where: North: 0°, East: 90°, South: 180°, West: 270°
/// - `distance`: depends on the [trait implementation](#implementors).
/// - returns: Point where the units of x/y depend on the [trait implementation](#implementors).
///
/// [`metric_spaces`]: super::metric_spaces
fn destination(origin: Point<F>, bearing: F, distance: F) -> Point<F>;
}
12 changes: 12 additions & 0 deletions geo/src/algorithm/line_measures/distance.rs
Original file line number Diff line number Diff line change
@@ -0,0 +1,12 @@
/// Calculate the distance between the `Origin` and `Destination` geometry.
pub trait Distance<F, Origin, Destination> {
/// Note that not all implementations support all geometry combinations, but at least `Point` to `Point`
/// is supported.
/// See [specific implementations](#implementors) for details.
///
/// # Units
///
/// - `origin`, `destination`: geometry where the units of x/y depend on the trait implementation.
/// - returns: depends on the trait implementation.
fn distance(origin: Origin, destination: Destination) -> F;
}
34 changes: 34 additions & 0 deletions geo/src/algorithm/line_measures/interpolate_point.rs
Original file line number Diff line number Diff line change
@@ -0,0 +1,34 @@
use crate::{CoordFloat, Point};

// REVIEW: Naming alternatives:
// - LinearReferencing
// - PointAlongLine
// - LineInterpolatePoint (postgis)
// - Interpolate (shapely)
// - Position (geographiclib)
// - Intermediate (georust::geo)
pub trait InterpolatePoint<F: CoordFloat> {
/// Returns a new Point along a line between two existing points
///
/// See [specific implementations](#implementors) for details.
fn point_at_ratio_between(start: Point<F>, end: Point<F>, ratio_from_start: F) -> Point<F>;

// TODO:
// fn point_at_distance_between(start: Point<F>, end: Point<F>, distance_from_start: F) -> Point<F>;

/// Interpolates `Point`s along a line between `start` and `end`.
///
/// See [specific implementations](#implementors) for details.
///
/// As many points as necessary will be added such that the distance between points
/// never exceeds `max_distance`. If the distance between start and end is less than
/// `max_distance`, no additional points will be included in the output.
///
/// `include_ends`: Should the start and end points be included in the output?
fn points_along_line(
start: Point<F>,
end: Point<F>,
max_distance: F,
include_ends: bool,
) -> impl Iterator<Item = Point<F>>;
}
77 changes: 77 additions & 0 deletions geo/src/algorithm/line_measures/metric_spaces/euclidean.rs
Original file line number Diff line number Diff line change
@@ -0,0 +1,77 @@
use super::super::Distance;
use crate::{GeoFloat, Point};

/// Operations on the [Euclidean plane] measure distance with the pythagorean formula -
/// what you'd measure with a ruler.
///
/// If you have lon/lat points, use the [`Haversine`], [`Geodesic`], or other [metric spaces] -
/// Euclidean methods will give nonsense results.
///
/// Alternatively, you *can* use lon/lat points with Euclidean methods if you first [`Transform`]
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I understand the point you're trying to make, but if you transform lon/lat coordinates they're not lon/lat coordinates anymore. Could we say something like "If you wish to use Euclidean operations with lon/lat coordinates these must first be be transformed using the transform / transform_crs_to_crs methods or their immutable variants. Use of these requires the proj feature"?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm OK with your phrasing — I've updated it.

For my own understanding, are you trying to avoid the usage of the term "projection"?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Not consciously, but yes: If you're working with non-Euclidean coordinates, you're probably not that interested in projection. If you were, you wouldn't be using measures and algorithms that are slower and more complex than their 2D counterparts.

This is all part of a larger endeavour, in the sense that geo is trying to be a practical tool with which to do real work. A great deal of that work can be done by transforming geographic coordinates into geodetic coordinates, giving you access to the simpler, more robust abstractions in the Euclidean plane.

But we also know (because the abstractions and algorithms now exist in geo because people have asked for / contributed them) that some work is more suited to spherical geometry.

As it stands, I feel like our docs and the way that some of our algorithms are organised don't demarcate that relationship very well; In my mind, we're "Euclidean and Cartesian by default, but if you need spherical geometry and algorithms we've got some of those".

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In my mind, we're "Euclidean and Cartesian by default, but if you need spherical geometry and algorithms we've got some of those".

Yeah, I think that accurately reflects georust/geo. I think it also pretty accurately reflects most historical GIS software.

This is all part of a larger endeavour, in the sense that geo is trying to be a practical tool with which to do real work.

💯 ❤️

If you were [interested in projection], you wouldn't be using measures and algorithms that are slower and more complex than their 2D counterparts.

With "The Web" as a publishing platform (e.g. webmaps, geojson) and more global datasets like OSM, using lon/lat coordinates in algorithms is increasingly not a conscious decision of the analyst/cartographer/person as "the best", but rather just happens to be the format that their source data is in, or that their presentation layer demands.

Having talked through this a bit, I think that's the group I'm trying to connect with, and often am implicitly biased towards in my documentation and API design.

/// your points to an appropriate projection.
///
/// [Euclidean plane]: https://en.wikipedia.org/wiki/Euclidean_plane
/// [`Transform`]: crate::Transform
/// [`Haversine`]: super::Haversine
/// [`Geodesic`]: super::Geodesic
/// [metric spaces]: super
pub struct Euclidean;

/// Calculate the Euclidean distance (a.k.a. pythagorean distance) between two Points
impl<F: GeoFloat> Distance<F, Point<F>, Point<F>> for Euclidean {
/// Calculate the Euclidean distance (a.k.a. pythagorean distance) between two Points
///
/// # Units
/// - `origin`, `destination`: Point where the units of x/y represent non-angular units
/// — e.g. meters or miles, not lon/lat. For lon/lat points, use the
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In a similar vein, I don't think we should be referring to "meters or miles" here: non-Euclidean distances are often measured in meters / miles, so this is confusing. "Not lon/lat" is OK, but I'd prefer that this is either implicit (you can't measure lon / lat distances in Euclidean space, and we are in Euclidean space, therefore…) or that we decline to explain with physical-world examples completely.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

non-Euclidean distances are often measured in meters / miles, so this is confusing.

I believe you, but I honestly didn't know that was true. Can you give me an example?

you can't measure lon / lat distances in Euclidean space

Oh but you can! And it's a common error that gives you a meaningless result. Look no further than our own example code on the legacy euclidean distance:

let p1 = point!(x: -72.1235, y: 42.3521);
let p2 = point!(x: -72.1260, y: 42.45);
let distance = p1.euclidean_distance(&p2);

assert_relative_eq!(distance, 0.09793191512474639);

I think it's important to keep simple relatable language and examples to help mitigate that.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I believe you, but I honestly didn't know that was true. Can you give me an example?

Any physical measurement of distance on the earth's surface (or a simplified abstraction of it) is definitionally non-Euclidean because it's a curved surface.

Copy link
Member Author

@michaelkirk michaelkirk Sep 25, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ah I see, that makes sense. Thank you for bearing with me. I don't have a formal GIS background, so a lot of things that may seem obvious, I truly don't know.

I've included the link to [Euclidean Plane], which I hope addresses your ambiguity concern.

I've also linked to the Transform trait in reference to what we mean by "projection". While also not 100% precise, I think that it will be helpful to those who know a little, without confusing those who know a lot.

I still think it's important that our documentation aim to help people avoid the most likely errors by saying something like:

this method takes lng/lat

But the negation:

this method does not take lng/lat

Feels awkward and unnecessarily mysterious without an example. I now understand how Euclidean space is not the only space that could use meters as a unit, but as an example of what we mean by "not lng/lat", I think it is clarifying, and unlikely to confuse.

So if I've said something that is incorrect, let me know and I'll try again to fix it. But if it's more so that you think it reads like it was written for dummies, well then I'd please beg of you to have mercy on us dummies.

/// [`Haversine`] or [`Geodesic`] [metric spaces].
/// - returns: distance in the same units as the `origin` and `destination` points
///
/// # Example
/// ```
/// use geo::{Euclidean, Distance};
/// use geo::Point;
/// // web mercator
/// let new_york_city = Point::new(-8238310.24, 4942194.78);
/// // web mercator
/// let london = Point::new(-14226.63, 6678077.70);
/// let distance: f64 = Euclidean::distance(new_york_city, london);
///
/// assert_eq!(
/// 8_405_286., // meters in web mercator
/// distance.round()
/// );
/// ```
///
/// [`Haversine`]: super::Haversine
/// [`Geodesic`]: super::Geodesic
/// [metric spaces]: super
fn distance(origin: Point<F>, destination: Point<F>) -> F {
crate::EuclideanDistance::euclidean_distance(&origin, &destination)
}
}

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

type MetricSpace = Euclidean;

mod distance {
use super::*;

#[test]
fn new_york_to_london() {
// web mercator
let new_york_city = Point::new(-8238310.24, 4942194.78);
// web mercator
let london = Point::new(-14226.63, 6678077.70);
let distance: f64 = MetricSpace::distance(new_york_city, london);

assert_relative_eq!(
8_405_286., // meters in web mercator
distance.round()
);
}
}
}
Loading
Loading