Refactor connect edges

[rowanwins]

Hi @w8r

Well I've made a fairly good start on refactoring the connect edges routine.

I need to do some more testing and check why the various tests aren't passing, it looks like part of it is just to do with whether the contours are closed or not so fingers crossed that should be a simple fix.

And I also need to check against some of the issues in the repo.

Anyway feels like reasonable progress.

[rowanwins]

Resolves

• Self intersecting edges in union #103
• Misattributed hole with intersection #96
• intersection returning invalid polygon #93
• Union show hole as multipolygon #71
• Incorrect inner ring added to diff result. #69
• outer ring added to wrong polygon in diff() result. #68
• Invalid polygon in expected result of test 'disjoint union testing' #61

[rowanwins]

I think this is now good @w8r

I made some minor modifications to the test suite where there were dodgy outputs specified.

Things are perhaps a touch slower but you might be able to work your magic to speed things up again now that they are slightly more stable again.

With this PR

Hole_Hole
Martinez x 36,537 ops/sec ±2.05% (87 runs sampled)
JSTS x 2,746 ops/sec ±2.31% (85 runs sampled)
- Fastest is Martinez

Asia union
Martinez x 9.11 ops/sec ±6.72% (25 runs sampled)
JSTS x 9.66 ops/sec ±6.27% (28 runs sampled)
- Fastest is JSTS

States clip
Martinez x 224 ops/sec ±1.16% (87 runs sampled)
JSTS x 134 ops/sec ±2.00% (78 runs sampled)
- Fastest is Martinez

In Master

Hole_Hole
Martinez x 36,993 ops/sec ±1.06% (90 runs sampled)
JSTS x 2,675 ops/sec ±1.52% (89 runs sampled)
- Fastest is Martinez

Asia union
Martinez x 11.57 ops/sec ±5.51% (33 runs sampled)
JSTS x 9.79 ops/sec ±4.80% (29 runs sampled)
- Fastest is Martinez

States clip
Martinez x 282 ops/sec ±3.07% (84 runs sampled)
JSTS x 125 ops/sec ±2.07% (80 runs sampled)
- Fastest is Martinez

[rowanwins]

One possible improvement is that we could ensure contours are going the correct way - clockwise or counterwise, I think that's what the precomputedCC field is for on the contour class, but I haven't touched that yet.

[bluenote10]

Nice, I have checked this PR against an issue I recently found in the Rust version, and with this change the result is now correct. On master the output is still two separate polygons (i.e., it missed to report a ring):

With this PR the ring is properly assigned to the polygon:

If desired I can do a more thorough review in the next days, testing further edge cases etc.

[rowanwins]

Awesome - thanks for checking it out @bluenote10
I think it would be worth having a few tests in your suite that looked specifically at nesting of contours as it's quite complicated and we don't want to introduce regressions down the track, so perhaps two or three of the issues listed above could be made into tests once your PR get's merged.

[bluenote10]

I have already prepared something ;)

[bluenote10]

I've looked a bit into the PR, and of first glance the refactoring looks pretty good to me. There are just a few things that still need to be investigated.

Diffing with hole

When creating the diff/xor of a polygon with hole to a larger outer polygon, the following happens:

The result correctly contains the outer ring, but there should also be an interior polygon in the result.

Test Data

``` { "features": [ { "geometry": { "coordinates": [ [ [0, 0], [1, 0], [1, 3], [0, 3], [0, 0] ], [ [0.1, 0.1], [0.1, 0.9], [0.9, 0.9], [0.9, 0.1], [0.1, 0.1] ], [ [0.1, 2.1], [0.1, 2.9], [0.9, 2.9], [0.9, 2.1], [0.1, 2.1] ] ], "type": "Polygon" }, "properties": {}, "type": "Feature" }, { "geometry": { "coordinates": [ [ [0.1, 1.1], [0.1, 1.9], [0.9, 1.9], [0.9, 1.1], [0.1, 1.1] ], [ [0.2, 1.2], [0.2, 1.8], [0.8, 1.8], [0.8, 1.2], [0.2, 1.2] ] ], "type": "Polygon" }, "properties": {}, "type": "Feature" } ], "type": "FeatureCollection" } ```

Union of nested polygons

This test simply takes the union of two multi polygons which each consists of an outer and inner "wall" -- one bigger than the other (blue shading is multi polygon A, red shading is multi polygon B):

The result should simply contain all 4 walls, but stops after the outer one.

Test Data

``` { "features": [ { "geometry": { "coordinates": [ [ [[-10, -10], [-10, 10], [10, 10], [10, -10], [-10, -10]], [[-9.5, -9.5], [-9.5, 9.5], [9.5, 9.5], [9.5, -9.5], [-9.5, -9.5]] ], [ [[-9, -9], [-9, 9], [9, 9], [9, -9], [-9, -9]], [[-8.5, -8.5], [-8.5, 8.5], [8.5, 8.5], [8.5, -8.5], [-8.5, -8.5]] ] ], "type": "MultiPolygon" }, "properties": {}, "type": "Feature" }, { "geometry": { "coordinates": [ [ [[-7, -7], [-7, 7], [7, 7], [7, -7], [-7, -7]], [[-6.5, -6.5], [-6.5, 6.5], [6.5, 6.5], [6.5, -6.5], [-6.5, -6.5]] ], [ [[-6, -6], [-6, 6], [6, 6], [6, -6], [-6, -6]], [[-5.5, -5.5], [-5.5, 5.5], [5.5, 5.5], [5.5, -5.5], [-5.5, -5.5]] ] ], "type": "MultiPolygon" }, "properties": {}, "type": "Feature" } ], "type": "FeatureCollection" } ```

[bluenote10]

Just a note about terminology: In the a paper a "contour" really is just any "linear ring / closed line segment", i.e., it can be either the exterior of polygon or polygon holes:

In that sense the class Contour is a bit confusing, because it doesn't seem to represent a single ring, but rather an entire polygon. Renaming the class to Polygon would be more consistent. But I don't fully understand the semantics of contourId yet, is it actually a "polygon id"?

[bluenote10]

I think I get it now -- the contours aren't storing the points of holes like a polygon would, they are merely references to other contour ids.

[bluenote10]

Why is this check needed (and the one in line 96)? Should the rings produced by the algorithm not be closed in any case? Or is this merely needed for the currently broken test cases that sometimes don't close the input polygons properly? If that's the case we should rather take care of that upfront.

[bluenote10]

What happens of contour.external is false? I think this leads to problems of relevant things not being added to the result.

[bluenote10]

This line doesn't need to be in the other scope, could be moved into line 93 to make it const outCoords = [contour.points].

[w8r]

true

[bluenote10]

I have looked into these two cases and think they have the same cause, and it is not related to the change in the PR:

For these nested linear rings, the expression resultEvents[i].prevInResult.resultInOut should basically alternate between outside-inside and inside-outside transition for each contour. However in the debugger I see that all nested contours read a value resultEvents[i].prevInResult.resultInOut == false, i.e., they all think they are inside (which is why they can't create new nested polygons). It looks like this is actually a bug in the first stage of the algorithm.

EDIT: Wait, resultInOut isn't computed in the first stage? I thought it is... If it is only computed at the bottom of connectEdges it might be related after all.

I think the problem is line 139 if (event.left) { event.resultInOut = false; ... }. It is too simple to make every left-event an outside-inside transition. In the Rust implementation I'm currently experimenting with a version that moves the entire resultInOut computation into compute_fields -- I think this is the best place to take care of it because that is where we also compute in_result. Results look promising, but it isn't fully finished yet.

[bluenote10]

Okay I think I've figured out now how the whole connect edges logic is supposed to work. For the Rust version of the algorithm I have all test cases working now, see 21re/rust-geo-booleanop#10 .

[w8r]

Is there a logical mistake somewhere in my version or have I just left something out?

[bluenote10]

Basically the logic around Figure 4 from the paper was missing, see the PR description for a high level explanation. Fixing it requires two minor changes:

• During the computeFields function I not only compute the boolean information whether an edge is inResult but rather the ternary information "not in result", "in result and in-out transition", "in result and out-in transition". The Rust version uses an enum for that, on JS we could simply encode it as a number. The relevant logic is:

// within compute_fields
let in_result = in_result(event, operation);
let result_transition = if !in_result {
    ResultTransition::None
} else {
    determine_result_transition(&event, operation)
};
event.set_result_transition(result_transition);

// using this helper function
fn determine_result_transition<F>(event: &SweepEvent<F>, operation: Operation) -> ResultTransition
where
    F: Float,
{
    let this_in = !event.is_in_out();
    let that_in = !event.is_other_in_out();
    let is_in = match operation {
        Operation::Intersection => this_in && that_in,
        Operation::Union        => this_in || that_in,
        Operation::Xor          => this_in ^  that_in,
        Operation::Difference   =>
            // Difference is assymmetric, so subject vs clipping matters.
            if event.is_subject {
                this_in && !that_in
            } else {
                that_in && !this_in
            }
    };
    if is_in { ResultTransition::OutIn } else { ResultTransition::InOut }
}

• When constructing a contour we have to consider whether the contour below is in-out vs out-in and work out the 4 cases of Figure 4. In Rust for instance this is used to construct a contour:

if let Some(prev_in_result) = event.get_prev_in_result() {
    // Note that it is valid to query the "previous in result" for its output contour id,
    // because we must have already processed it (i.e., assigned an output contour id)
    // in an earlier iteration, otherwise it wouldn't be possible that it is "previous in
    // result".
    let lower_contour_id = prev_in_result.get_output_contour_id();
    if prev_in_result.get_result_transition() == ResultTransition::OutIn {
        // We are inside. Now we have to check if the thing below us is another hole or
        // an exterior contour.
        let lower_contour = &contours[lower_contour_id as usize];
        if let Some(parent_contour_id) = lower_contour.hole_of {
            // The lower contour is a hole => Connect the new contour as a hole to its parent,
            // and use same depth.
            contours[parent_contour_id as usize].hole_ids.push(contour_id);
            let hole_of = Some(parent_contour_id);
            let depth = contours[lower_contour_id as usize].depth;
            Contour::new(hole_of, depth)
        } else {
            // The lower contour is an exterior contour => Connect the new contour as a hole,
            // and increment depth.
            contours[lower_contour_id as usize].hole_ids.push(contour_id);
            let hole_of = Some(lower_contour_id);
            let depth = contours[lower_contour_id as usize].depth + 1;
            Contour::new(hole_of, depth)
        }
    } else {
        // We are outside => this contour is an exterior contour of same depth.
        let depth = contours[lower_contour_id as usize].depth;
        Contour::new(None, depth)
    }
} else {
    // There is no lower/previous contour => this contour is an exterior contour of depth 0.
    Contour::new(None, 0)
}

@rowanwins version stores the holeOf and depth information in an array/hashmap, but I decided to directly attach it to the contour for better performance.

[rowanwins]

Can we do anything here to help move this one forward @w8r .

With the things you mention above @bluenote10 are any changes required to my PR or are you just saying that your approach was slightly different (perhaps slightly more performant).

[bluenote10]

Yes, I'm pretty certain that the changes mentioned are needed for the algorithm to work correctly. This is the core problem:

[...] the problem is line 139 if (event.left) { event.resultInOut = false; ... }. It is too simple to make every left-event an outside-inside transition.

The test cases I've posted above fail because of that and work fine with the changes applied. Basically the changes are just an exact modelling of the logic proposed in the paper (Fig. 4, see the Rust PR for an explanation).

[w8r]

Shall I try and add this logic from your Rust implementation to see how it plays?

[bluenote10]

Shall I try and add this logic from your Rust implementation to see how it plays?

I also thought about contributing it, and was mainly waiting for #117 to get merged, which would allow to simply copy over all the new test cases from the Rust version. Having it implemented by a more experienced JS developer won't hurt though ;)

[rowanwins]

Hi @bluenote10

I took a look but my rust skills aren't that flash so I had some trouble interpreting, I think I understand the main point of your proposal, I'm just not very sure how to code it up.

Perhaps @w8r we can make @bluenote10 a contributor to the repo so he can work on this branch also, rather than another branch via his fork?

Cheers

[w8r]

@rowanwins, I added @bluenote10 as a collaborator

[w8r]

@bluenote10 I saw you have successfully merged master branch in your PR, did you also have failing test cases? I have 8 of them

[bluenote10]

I didn't check that, I only ran the tests after I completed the refactoring. But on #121 all tests are passing.

[rowanwins]

I think all the test cases are passing @w8r

[w8r]

Guys, I cannot find time to take a closer look at this PR, I am sorry
Tests are passing, benchmarks look fine to me, I say I can release and not "let the better become the enemy of the good"
Thank you very much for your help! I will release now and see if we have any performance problems or anything

[bluenote10]

👍 In this context 21re/rust-geo-booleanop#15 might be interesting as well. It is basically an addon to this PR which produces nicer contours and can improve performance in cases where the next_pos searches unnecessarily in the wrong direction repeatedly.

[underbluewaters]

Thank you for this folks. I ran into problems with @turf/difference misplacing interior rings and using v0.7.0 solved it. Upgrading also seemed to give a big performance boost when clipping against big multi-megabyte features.