georust / geo

Geospatial primitives and algorithms for Rust
https://crates.io/crates/geo
Other
1.49k stars 194 forks source link

HaversineDestination results don't wrap at the antimeridian as expected #1074

Closed apendleton closed 9 months ago

apendleton commented 10 months ago

During a recent refactor, we switched some logic that had previously used HaversineIntermediate to use HaversineDestination instead, in a way that should have been equivalent (calculating great-circle intermediate values between two points). The results do seem in general to be equivalent, but HaversineDestination sometimes produces result values that are outside [-180,180] longitude even when the inputs are in-bounds, which I found somewhat surprising.

Looking at the docs, it doesn't seem like geo actually commits to any particular longitude representation ([-180,180] vs [0,360] in particular), so maybe this is left undefined on purpose and applications should just normalize? But notably, the equivalent HaversineIntermediate operation does seem to wrap longitudes correctly as I'd expect. Seems like there's some value in consistency, or at least in documenting expected behavior.

Minimal rust-script reproduction that compares in-theory-equivalent HaversineDestination and HaversineIntermediate operations calculating points along an AM-crossing line:

#!/usr/bin/env rust-script
//! ```cargo
//! [dependencies]
//! geo = "0.26"
//! ```

use geo::{ point, HaversineDistance, HaversineDestination, HaversineBearing, HaversineIntermediate };

fn main() {
    let pt1 = point! { x: 170.0, y: -30.0 };
    let pt2 = point! { x: -170.0, y: -30.0 };

    for (start, end) in [(pt1, pt2), (pt2, pt1)] {
        println!("Going from {:?} to {:?}...", start.x_y(), end.x_y());

        let bearing = start.haversine_bearing(end);
        let total_distance = start.haversine_distance(&end);

        let mut dist = 0.;
        while dist < total_distance {
            let dest_point = start.haversine_destination(bearing, dist);
            let inter_point = start.haversine_intermediate(&end, dist / total_distance);

            println!("destination: {:?}, intermediate: {:?}", dest_point.x_y(), inter_point.x_y());

            dist += 250_000.;
        }
    }
}

prints:

Going from (170.0, -30.0) to (-170.0, -30.0)...
destination: (170.0, -29.999999999999996), intermediate: (170.0, -29.999999999999996)
destination: (172.5908050903293, -30.172086178256627), intermediate: (172.5908050903293, -30.172086178256635)
destination: (175.18932441403777, -30.293113311281495), intermediate: (175.18932441403777, -30.29311331128149)
destination: (177.79289734601838, -30.362706690084703), intermediate: (177.79289734601838, -30.362706690084703)
destination: (180.39880075944293, -30.380649596934468), intermediate: (-179.6011992405571, -30.380649596934468)
destination: (183.00428252156786, -30.346886007528365), intermediate: (-176.99571747843217, -30.346886007528365)
destination: (185.6065957562343, -30.261521293335424), intermediate: (-174.39340424376573, -30.261521293335424)
destination: (188.2030329670417, -30.124820901911836), intermediate: (-171.79696703295832, -30.12482090191183)
Going from (-170.0, -30.0) to (170.0, -30.0)...
destination: (-170.0, -29.999999999999996), intermediate: (-170.0, -29.999999999999996)
destination: (-172.5908050903293, -30.172086178256627), intermediate: (-172.5908050903293, -30.172086178256635)
destination: (-175.18932441403777, -30.293113311281495), intermediate: (-175.18932441403777, -30.29311331128149)
destination: (-177.79289734601838, -30.362706690084703), intermediate: (-177.79289734601838, -30.362706690084703)
destination: (-180.39880075944293, -30.380649596934468), intermediate: (179.6011992405571, -30.380649596934468)
destination: (-183.00428252156786, -30.346886007528365), intermediate: (176.99571747843217, -30.346886007528365)
destination: (-185.6065957562343, -30.261521293335424), intermediate: (174.39340424376573, -30.261521293335424)
destination: (-188.2030329670417, -30.124820901911836), intermediate: (171.79696703295832, -30.12482090191183)

Happy to help out with a PR or whatever if there's some consensus about what, if anything, ought to be done here, but just wanted to flag since it caught me off-guard.

Also FWIW, I compared the behavior to turf.destination in turf.js, and again somewhat to my surprise, it exhibits the same past-180 wrapping behavior, so maybe this is normal?

> turf.destination([-170.0, -30.0], 1750000.0, -95.03836877329749, {units: 'meters'})
{
  type: 'Feature',
  properties: {},
  geometry: {
    type: 'Point',
    coordinates: [ -188.2030329670417, -30.124820901911832 ]
  }
}
urschrei commented 10 months ago

The matching turf behaviour makes sense as HaversineDestination is inspired by the turf function: https://github.com/georust/geo/pull/124.

I'm in favour of modifying the method to match HaversineIntermediate. Anyone else have thoughts?