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

Rhumb destinations do not wrap angles after the first pole intersection #1152

Open
joe-saronic opened this issue Feb 21, 2024 · 8 comments · May be fixed by #1153
Open

Rhumb destinations do not wrap angles after the first pole intersection #1152

joe-saronic opened this issue Feb 21, 2024 · 8 comments · May be fixed by #1153

Comments

@joe-saronic
Copy link
Contributor

Given that a rhumb line has finite length to a given pole, and that RhumbDestination appears to wrap correctly in one direction, one would expect that latitude would be bounded to +/-90 and longitude to +/-180 as the path length increased. Instead, we see the following:

image

Here is a rust program that dumps the data:

use geo::{Point, RhumbDestination};

use std::fs::OpenOptions;
use std::io::{BufWriter, Error, Write};

fn main() -> Result<(), Error> {
    let write = OpenOptions::new().write(true).create(true).open("./rhumb.csv");
    let mut writer = BufWriter::new(write.unwrap());
    let origin = Point::new(0.0, 0.0);
    for i in 0..1000000 {
        let dist = 100.0 * (i as f64);
        let dest = origin.rhumb_destination(45.0, dist);
        let lon = dest.x();
        let lat = dest.y();
        writer.write(format!("{dist:.16},{lon:.16},{lat:.16}\n").as_bytes())?;
    }
    Ok(())
}

and some python code to load and plot:

import pandas as pd
from matplotlib import pyplot as plt

df = pd.read_csv('rhumb.csv', header=None)
fig, ax = plt.subplots(2, 1, constrained_layout=True)
ax[0].plot(df[0], df[1])
ax[1].plot(df[0], df[2])
ax[0].set_title('Lon vs Dist')
ax[1].set_title('Lat vs Dist')

The expected fix would be that latitude would be clamped, either by taking the modulo of the path length from pole to pole, or empirically.

@joe-saronic joe-saronic changed the title Rhumb destinations are unreliable for large distances Rhumb destinations do not wrap angles after the first pole intersection Feb 21, 2024
@urschrei
Copy link
Member

@apendleton this is your feature. I don't think this is a bug, but I agree that clamping lon and lat would be preferable. WDYT, and do you have capacity to add it?

@joe-saronic
Copy link
Contributor Author

Returning a latitude <-90 is pretty close to being a bug :) I will take a look if I have time myself.

@apendleton
Copy link
Contributor

For reference, here's a quick port of this test program to Javascript so I could test it with Turf (which is what the Rust implementation was ported from):

const fs = require('fs');
const turf = require('@turf/turf');

function main() {
    const out = [];
    const origin = [0.0, 0.0];
    for (let i = 0; i < 1000000; i++) {
        let dist = 100.0 * i;
        let dest = turf.rhumbDestination(origin, dist, 45.0, {units: 'meters'}).geometry.coordinates;
        out.push(`${dist},${dest[0]},${dest[1]}`);
    }
    fs.writeFileSync("./turf_rhumb.csv", out.join('\n'));
}

main();

And it produces broadly similar output with respect to latitude, but probably-saner output with respect to longitude:

image

With respect to longitude, we do already attempt to ensure that we're wrapping correctly, but it looks like there's maybe some corner-case weirdness that results in us producing some occasional out-of-bounds values as evidenced by this graph. I think wrapping (not clamping) is the preferred behavior for longitude, so probably we need to figure out what's going on that results in us occasionally producing longitudes <-180 in the Rust version.

For latitude, I'm actually not sure what we want the behavior to be for very long lines. I think clamping at +/- 90 is defensible, but it's not clear if we just allow the longitude number to continue to spin meaninglessly in that case or what. As I understand it, conceptually, for non-horizontal lines, we should be spiraling towards the pole for awhile and then eventually reach it at some fixed/calculable distance, after which further movement along the line is undefined.

I can pitch in however is most useful, but it might be good to see if any other implementations exist out in the wild to see if there's a standard/expected behavior in this circumstance? I don't have a clear frame of reference here except for Turf (which seems to also be arguably wrong).

@joe-saronic
Copy link
Contributor Author

Returning an Option<...> which is None when you've gone too far seems like a reasonable approach. The maximum distance is proportional to the difference in latitudes and the secant of the bearing, so is not hard to calculate up-front from a given latitude. For bearing = +/- 90, the option would always be Some, with longitude determined by the local radius.

@apendleton
Copy link
Contributor

An Option seems reasonable to me too, though that'd be a breaking API change, and would make rhumb destination inconsistent with haversine destination in terms of signature, so I could go either way.

As to how to calculate: yeah, makes sense. My use of geo tends to be in a pretty performance-sensitive context and trig is slow, so given that in the typical case I think we'd probably expect not to hit the pole, I might sort of be inclined to avoid doing that secant every time if it usually won't matter, vs. checking after the fact. I suppose that's empirically testable though.

@urschrei
Copy link
Member

I would say API inconsistency isn't a big deal; it's accounting for a different detail after all. The breaking change is a concern. On the other hand, the API isn't in wide use, and the change is certainly justifiable.
I can't speak to the perf concerns – we should look at the numbers first.

@apendleton
Copy link
Contributor

Makes sense. I'll probably have time to dig in in more detail either Friday or this weekend, but happy to cede this to someone else if there's interest in it being done sooner.

@joe-saronic joe-saronic linked a pull request Feb 22, 2024 that will close this issue
2 tasks
@joe-saronic
Copy link
Contributor Author

Turns out to be a lot simpler than I expected. Early return is actually the faster/simpler approach here since the conditional is already being checked anyway.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging a pull request may close this issue.

3 participants