Dateline "no-wrap" odd in some scenarios

[corporealfunk]

I'm experiencing an issue with wrapping that you can recreate on the blog demo:

https://blog.cyclemap.link/Leaflet.Geodesic/nowrap-interactive.html

Simply drag the marker that is in southern CA of the United States and drag it to, say, Portugal.

I get lines like this by creating a line with the following two points and setting wrap to false:

const lax = [33.9416, -118.4085];
const nrt = [35.7720, 140.3929];

JSBin: https://jsbin.com/lamihiloda/edit?html,js,output

[henrythasler]

Yeah, I know. I was hoping nobody would notice... :-) I will think of something to address this issue.

[henrythasler]

There could be 2 solutions:

Which one would make more sense to you?

[corporealfunk]

Shifting either way is OK for my needs, though I was hoping there might be a solution in the library. The shifting solution requires that we keep track of when points cross the antimeridian and then shift all subsequent points in that line:

(the following are some of the GPS waypoints from an international flight plan from LAX to NRT):

const lax = [33.9416, -118.4085];
const rizon = [53, -170];
const eyeti = [53, 180 - 360]; // crossed antimeridian, shift

// all following points must shift:
const opake = [53.001222, 169.89275 - 360];
const olcot = [51.430194, 165.55575 - 360]
const nrt = [35.7720, 140.3929 - 360];

It's not a difficult solution to do this in client code, though it would be nice if the library were able to sort this out for us instead of us having to keep track as we add points to a line.

JS bin example of the above: https://jsbin.com/duyivevike/1/edit?html,js,output

[henrythasler]

I haven't had a good idea yet how to resolve this issue w/o introducing lots of if-statements which I don't want to do. If you can live with your workaround then I'd like to leave it for now. Maybe someone has a good idea in the future.

[stale[bot]]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. Thank you for your contributions.

[cherbert]

Is it rude of me to reopen this old wound? It's a real shame there isn't a fix for this. This problem is plaguing my project. I'm trying to migrate from Google Maps which doesn't have this glitch.

It's creating a bit of a mess on this aircraft flight plan map. https://i.imgur.com/uonfJuI.png

[henrythasler]

Absolutely not. Let's see if I can help you with that one.

[bugfindr]

@henrythasler I am also facing this issue. Any update on this?

[henrythasler]

I have just started to work on that issue. Let's see how is goes.

[zakjan]

If you find a solution, this is exactly one of the cases I was looking for sharing at https://github.com/henrythasler/Leaflet.Geodesic/issues/64 for my https://github.com/zakjan/mapbox-gl-draw-geodesic (WIP)

[zakjan]

Update: I have used https://github.com/springmeyer/arc.js for the calculation, then I shift the longitudes past the antimeridian crossing in https://github.com/zakjan/mapbox-gl-draw-geodesic/blob/705032a18f6af3a3b44c9ae898c16eb95336a5d7/src/utils/create_geodesic_line.js as @corporealfunk suggested.

Demo is available at https://zakjan.github.io/mapbox-gl-draw-geodesic/

@cherbert Do you mind to share your dataset, so that I could test it with?

[bugfindr]

@zakjan Can you please provided simple Geodesic Polyline example code as per your solution JavaScript so I can try that with my data?

[zakjan]

@bugfindr Source code for the geodesic line computation is in the file I linked above ^^ https://github.com/zakjan/mapbox-gl-draw-geodesic/blob/811c55ec32d6fc9c0fcb06a43c9b8d7e35c782c3/src/utils/create_geodesic_line.js

I'm sorry for hijacking this thread. If anyone has questions, please fill them in mapbox-gl-draw-geodesic repo.

[henrythasler]

The interesting thing is, that the features are repeated over the visible map if more than 360° are shown. Is this on purpose?

[zakjan]

@henrythasler This is the default Mapbox behavior for all features. Even with repeatWorldCopies: false (Mapbox equivalent of Leaflet noWrap: true), features crossing antimeridian are rendered on both sides of the single world.

I wonder why Leaflet core doesn't do the same. I only found Leaflet.RepeatedMarkers and L.VectorGrid.Slicer plugins for it.

[henrythasler]

I found a limited solution that covers corporealfunk's lax-nrt-testcase (red line). It still has some limitations if source and target a too far apart (blue).

If this is acceptable for now I'd like to close the issue. Please let me know what you think.

[henrythasler]

Need feedback and a release before closing.

[bugfindr]

@henrythasler I have tested it and still seeing issue. Please try below example latlng. From LatLng: (34.5975,-117.383) To LatLng: (28.5661,77.1061)

Btw I am using version 2.5.4 as its seems to me latest. Thanks!

[henrythasler]

I haven't created a release yet, so you won't be able to see any changes. I will create a pre-release for testing this week.

[bugfindr]

Ok. Let me know when we ready to test. Thanks!

[henrythasler]

You can now try with 2.5.5-0 pre-release

[bugfindr]

We have issue with this pre-release as well If plotting more than two points together. I will provide example soon here.

[bugfindr]

@henrythasler sorry for delay response. Here's issue I am facing with pre-release in case of multiline. I am attaching demo html file to for your test as well. Thanks!

Multilineline-test.zip

Output image.

[nolatron]

I tried the prerelease 2.5.5-0 and still getting the cut for flight path line over the pacific.

Here's what I'm currently using to create the polyline:

var start3 = new L.LatLng(33.9425, -119.593); var end3 = new L.LatLng(37.4625, 126.439); const polyline3a = L.geodesic([start3, end3], { color: '#22A7F0', weight: 3, dashArray: "15 10", opacity: 1 });

In a previous 2017 version that I had still, the line would sorta work if I added 360 to the -119 coordinate, but that doesn't seem to work in the 2.5.4 or 2.5.5-0.

[stale[bot]]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. Thank you for your contributions.

[stale[bot]]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. Thank you for your contributions.

[stale[bot]]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. Thank you for your contributions.

[bugfindr]

@henrythasler Any update on this?

[henrythasler]

No, not yet... Sorry.

[stale[bot]]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. Thank you for your contributions.

[stale[bot]]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. Thank you for your contributions.

[matafokka]

Hello, seems like I've managed to find a fix for this issue.

As suggested above, I looked at arc.js and source with formulas behind it. Algorithm that finds intermediate points depends on distance between two points along great circle. However, algorithm that calculates that distance will wrap when difference between lngs is greater than 180 degrees, and give distance along small part of a circle. Thus, first algorithm wraps too. The fix is to use a distance along big part of a circle (360 - d).

I've written this code, and it seems to work fine with different kinds of points (i.e. looks identical to what this library produces when it doesn't wrap):

const turfHelpers = require("@turf/helpers"); // I need it for other stuff, so I use it here to convert degs to rads and back

// Original points
let lng1Orig = -120, lat1Orig = -50, lng2Orig = 80, lat2Orig = 40,
    lng1 = turfHelpers.degreesToRadians(lng1Orig),
    lat1 = turfHelpers.degreesToRadians(lat1Orig),
    lng2 = turfHelpers.degreesToRadians(lng2Orig),
    lat2 = turfHelpers.degreesToRadians(lat2Orig),

    // Calculate great circle distance between these points
    w = lng2 - lng1,
    h = lat2 - lat1,
    z = Math.sin(h / 2) ** 2 + Math.cos(lat1) * Math.cos(lat2) * Math.sin(w / 2) ** 2,
    d = 2 * Math.asin(Math.sqrt(z)),
    coords = [];

// These formulas will "wrap" (calculate distance along small part of a circle between two points) when lng difference
// is greater than 180 degrees. Algorithm below is using that distance and wraps too.
// The fix is to use a distance along big part of a circle (360 - d).

if (Math.abs(w) > Math.PI)
    d = Math.PI * 2 - d;

for (let f = 0; f <= 1; f += 0.1) {
    let A = Math.sin((1 - f) * d) / Math.sin(d),
        B = Math.sin(f * d) / Math.sin(d),
        x = A * Math.cos(lat1) * Math.cos(lng1) + B * Math.cos(lat2) * Math.cos(lng2),
        y = A * Math.cos(lat1) * Math.sin(lng1) + B * Math.cos(lat2) * Math.sin(lng2),
        z = A * Math.sin(lat1) + B * Math.sin(lat2),
        lat = turfHelpers.radiansToDegrees(Math.atan2(z, Math.sqrt(x ** 2 + y ** 2))),
        lng = turfHelpers.radiansToDegrees(Math.atan2(y, x));

    let coord = [lat, lng];
    L.marker(coord).bindPopup(`${Math.round(f * 10)}: ${lng}, ${lat}`).addTo(map);
    coords.push(coord);
}

L.polyline(coords).addTo(map);

Since my knowledge in geodesy is little to none, I have almost no idea how these formulas work, and I found this solution by banging my head against my keyboard, please, verify it on some kind of dataset before using it.

I did more testing with segments length (length of each segment should be equal to line length / f), and this thing seems to be working just fine. This feature actually is pretty useful for aerial photography.

[henrythasler]

Thanks for this! Let me check it in the next days/weeks when I have some free time.

[matafokka]

Thank you, I'll be waiting for your feedback :)

Meanwhile, I'll do more testing and post an update when I'm done. I might even try to make a PR, but no promises here.

By the way, are you ok with ditching Vincenty's formula for this?

[matafokka]

Hello again, have I misunderstood the question or I just need a break right now? I thought that, despite glitches, a line should follow big part of a circle, i.e. stretch across the map, isn't it the case? Because I started testing and got this:

LInes through the same points ([-120, 4], [120, 4]), red line is produced by your library, blue - by my code above. Both lines indeed follow a great circle, just different parts of it. So, do I understand the problem correctly?

Anyway, even if your library should always follow the small part:

1. Code above still should work, just remove if (Math.abs(w) > Math.PI), and you'll basically get what arc.js does.
2. I have some speculations to why glitches occur. I'll clone the repo tomorrow and see what I can do.

Since I mentioned it, can I request an option that will make lines follow a big part of a circle (i.e. do what my code does)? I'd love to implement it, but I don't know how Vincenty's formula work and can't modify it.

[matafokka]

Yuppie, I found a solution!

The idea behind your wrapping algorithm is right, but it should be implemented a bit different. Check this code, I wrote an explanation on how it works in the comments:

/**
 * Finds in which "map" lies the first point, i.e. what number revolutions should we do to get to given lng.
 * @param lng {number} Lng to work with
 * @return Signed number revolutions
 */
private static numberOfRevolutions(lng: number) {
    return Math.sign(lng) * Math.floor(Math.abs(lng) / 180);
}

wrapMultiLineString(multilinestring: L.LatLng[][]): L.LatLng[][] {
    const result: L.LatLng[][] = [];

    multilinestring.forEach((linestring) => {
        let resultLine: L.LatLng[] = [], firstLng = linestring[0].lng, prevLng = firstLng,
                mapNumber = GeodesicGeometry.numberOfRevolutions(firstLng);

        for (let point of linestring) {
            // Find map number difference and move current point to the map containing first point
            let lng = point.lng, mapNumberDiff = GeodesicGeometry.numberOfRevolutions(lng) - mapNumber;

            if (mapNumberDiff !== 0)
                lng -= mapNumberDiff * 180;

                // If difference between lngs is greater than 90, line will be split. We should shift points with
                // such difference depending on their position on the map after split (left or right).

                // Thus, we wrapped lines to one map, but if original coordinates do more than one revolution,
                // we have to adjust lngs accordingly using the same principle.
                // With more revolutions, problem occur more times, so we'll adjust until we get acceptable difference.
                while (Math.abs(lng - prevLng) > 90)
                    lng -= 180 * Math.sign(lng - prevLng);

                prevLng = lng;
                resultLine.push(new L.LatLng(point.lat, lng));
        }
        result.push(resultLine);
    });
    return result;
}

This is before and after changes:

And polyline (red is wrapped, blue is split, both line are using same coordinates, don't mind markers, the behavior is what you've intended):

I also have a couple suggestions not related to this problem:

1. We can add method I described above as an alternative. It has a number of cool features I mentioned, and it's faster than yours, though, it might be less precise.
2. We can replace forEach() with for loops which will result in increased performance on underpowered devices and will be good for rerendering (i.e. when moving markers like in your demos). It happens because function calls are ~2x (if I remember correctly) times slower than iteration of the loop. You can google a comparison.

So anyway, my algorithm seems fine, but I'll further test it and create a PR with all my suggestions in a couple of days. I hope I'll have time, haha.

Edits: don't mind those, I'm just done for today. Meaningful edit 1: updated code to work with more than one revolution.

[matafokka]

Hey, another update! We can't wrap every segment of a polyline to only one map, if that polyline creates more than one revolution:

The lowest segment of a red line obviously should be shifted to the left. To move along shown path, we'll have to do two revolutions around the Earth. The more revolutions line do, the more shifts are required. I still don't fully understand why is that so, but I fixed it.

A working example with three revolutions, up to 4 shifts has been required to correctly wrap the line:

I've updated the code with these changes.

[matafokka]

Another update: wrapMultiLineString() now depends on midpoints, thus, it fails tests where lines are given without midpoints.

I think, you should update tests to use entire Geodesic class and somehow hide this function, so it won't be misused. It might hurt compatibility, I can't think of a clear solution for this.

Anyway, I won't touch both tests and API in my PR, it's up to you to resolve these issues.

The other method I wanted to add can work fine with only two points. I suppose, we can make wrapMultiLineString() call it as a compromise. Results might not be 100% same, but it'll do the thing, although, with less precision.

[henrythasler]

By the way, are you ok with ditching Vincenty's formula for this?

No I'm not. Vincenty's formula is the reason why I created this lib. :-D

[henrythasler]

The feature/better-wrap-next branch now works with simple lines and the test-cases. I hope I can make also make it work for multi-linestrings.

[matafokka]

@henrythasler, you have to wrap by segments between two source points and then shift them. The case with antimeridians is really hard to deal with when there's multiple points. I've been trying to do that in one loop for hours and failed. I'm currently doing segments.

Your solution also does better with this issue, but still doesn't work sometimes. Try [[25.799891182088334, -1546.875], [37.71859032558816, -286.87500000000006]] points.

Combination of previous points and half-working segments seems to work fine in this case:

Please, see this code:

wrapMultiLineString(multilinestring: L.LatLng[][]): L.LatLng[][] {
    const result: L.LatLng[][] = [], pointsInSegment = 2 ** (this.steps + 1);

    for (let linestring of multilinestring) {
        let segmentsCount = (linestring.length - 1) / pointsInSegment,
                resultLine: L.LatLng[] = [], segments: L.LatLng[][] = [];

        for (let currentSegment = 0; currentSegment < segmentsCount; currentSegment++) {
            let segment = [],startIndex = currentSegment * pointsInSegment,
                endIndex = (currentSegment + 1) * pointsInSegment,
                firstLng = linestring[startIndex].lng, prevLng = firstLng,
                shouldCopy = Math.abs(firstLng - linestring[endIndex].lng) <= 180;

            console.log(currentSegment, startIndex, endIndex);

            for (let i = startIndex; i <= endIndex; i++) {
                let point = linestring[i], lng = point.lng;

                // If difference between start point and end point is less than 180 degrees, line won't be split,
                // so just copy the whole string. It solves edge case when last point is close to the antimeridian
                // from the "left" side.
                if (shouldCopy) {
                    segment.push(point);
                    continue;
                }

                let untilDiff = i === endIndex ? 90 : 180;
                while (Math.abs(lng - prevLng) >= untilDiff) {
                    lng -= 180 * Math.sign(lng - prevLng);
                }
                prevLng = lng;
                segment.push(new L.LatLng(point.lat, lng));
            }

            segments.push(segment);
        }

        let moveLngTo;
        for (let segment of segments) {

            if (moveLngTo === undefined) {
                resultLine.push(...segment);
                moveLngTo = segment[segment.length - 1].lng;
                continue;
            }

            let lngDiff = moveLngTo - segment[0].lng;

            for (let i = 1; i < segment.length; i++) {
                let point = segment[i];
                point.lng += lngDiff;
                resultLine.push(point);
            }

            moveLngTo = segment[segment.length - 1].lng;
        }

        console.log("_____")
        result.push(resultLine);
    }
    return result;
}

It still sometimes fails on antimeridians (I suppose, we need to shift end point one more time, if difference is still <= 90) and doesn't count segments correctly.

---

We also can try to use your solution, shift the last point to solve this issue and combine it with the segments.

Please, try to do so. I'm gonna try to fix my code. Let me know, if you'll succeed. I'm too gonna post more comments, if I'll get something working.

[matafokka]

Updates!

Both this stuff and solution in my PR fails miserably on some near-antimeridian cases when small number of points is given.

I guess, it's better to use your formula and just shift the last point, if difference is insane.

Split too fails on antimeridians:

Points on screenshot:

[
    {
        "lat": 43,
        "lng": -468
    },
    {
        "lat": 52.908902047770255,
        "lng": -290.39062500000006
    },
    {
        "lat": -80.05804956215623,
        "lng": -470.39062500000006
    }
]

[henrythasler]

I think I've got it now. See feature/better-wrap-next. Unit-tests are still a bit weird but the interactive demo testing-multiline-nosplit.html works.

[matafokka]

Hey, works like a charm! It also passes all edge cases I've encountered! Good job!

[henrythasler]

A new release with this fix will be available within a few days.