Avoid generating self-intersecting isochrone polygons that would otherwise result from ::Generalize() #3026

ktatterso · 2021-04-21T23:09:55Z

Issue

When generating isochrones, customers can specify a numeric value, generalize or gen_factor, which simplifies the returned polygon. This is particularly useful for reducing the size of the polygon, and hence, the returned payload.

The polyline/polygon simplification method employed by Valhalla is the Douglas-Peucker algorithm. On its own, it has no facilities that would prevent it from generating a self-intersecting polygon during its simplification step. Based on my studies, using a gen_factor=55, auto costing, and a 30 min contour - the probability of us returning a self-intersecting polygon is quite high, I'd estimate 80-90% (urban areas).

This is problematic for customers/users who take these isochrone polygons and input them to other GIS systems for analysis - only to find their analysis tools don't work with self-intersecting polygons.

My solution to avoid generating self-intersections all boils down to this question: will simplifying this polyline result in a self-intersection? To determine the answer to that question poses a spatial problem: are there any polygon points in the "vicinity" of the polyline that would cause a self-intersection with the polyline I might simplify?

To answer the question "are there any points near this polyline I'm about to simplify?" I have created the PointTileIndex class. This class can discretize the lat/lon space into whatever width is desired - in this case, gen_factor determines the width. The PointTileIndex constructor will "bin" all the given points into the index. Then, within the Douglas-Peucker algorithm, we can use the PointTileIndex to answer that very question: "what points are near this polyline?"

Armed with a set of "points near my polyline", we can very accurately ask, "will any of these nearby points cause a self-intersection if I simplify?" Geometrically, you can think of line simplification as the geometric collapse of many triangles. The triangles in question are formed by the start position, the end position, and each point along the polyline. Hence, if a stray polygon point is found in any of the triangles we're about to collapse, we would cause a self-intersection if we simplified.

Example 1

Location: 47.412426,9.334023
Costing: Auto
Contour minutes: 30

Generalize: 0 (no self-intersection)

Generalize: 55 (without fix, self-intersects)

Generalize: 55 (after fix, no self-intersection)

Example 2

Location: 34.061657,-118.268494
Costing: Auto
Contour minutes: 30

Generalize: 0 (no self-intersection)

Generalize: 55 (without fix, self-intersects)

Generalize: 55 (after fix, no self-intersection)

Performance

The ::tile() method is O(n). Asking for nearby points is k*O(1) where k is the number of nearby tiles. Fortunately, most segments considered for simplification aren't very long, so k is small.

For example, the 30-min auto isochrone (gen_factor=55) for Orlando previously took 560 ms (not including serialization). The ::Generalize() step took 6 ms of that 560 ms. After my changes, ::Generalize() takes 22 ms, increasing the total time to 578 ms.

Thoughts & other ideas

I tried using Boost's boost::geometry::buffer() method with a distance_strategy=0.0, thinking if I gave it a self-intersecting polygon it would fix it. Unfortunately, it did not. Given a self-intersecting polygon it would return an empty result. In contrast, it worked fine for non-self-intersecting polygons.

I also wondered about other routines that might heal/repair/fix self-intersecting polygons. In the end, I decided that this approach would keep with the spirit of the original algorithm, only avoiding self-intersections as needed.

Tasklist

Add tests
Add #fixes with the issue number that this PR addresses
Update the docs with any new request parameters or changes to behavior described
Update the changelog

…ersections.

…tion test cases that were self intersecting. Lots of cleanup.

…tions to have new results without self-intersections. Much tidying.

…2 tests.

gknisely · 2021-04-22T12:11:31Z

@ktatterso Can we get some before and after screen shots of the results.

merkispavel

The first round of review

test/isochrone.cc

valhalla/midgard/polyline2.h

valhalla/midgard/point_tile_index.h

src/midgard/point_tile_index.cc

src/midgard/polyline2.cc

valhalla/midgard/point_tile_index.h

src/midgard/point_tile_index.cc

src/midgard/polyline2.cc

valhalla/midgard/point_tile_index.h

ktatterso · 2021-04-22T15:21:28Z

@ktatterso Can we get some before and after screen shots of the results.

I added some examples. thanks.

as not to overflow the underlying Tile id (int32_t).

ktatterso · 2021-04-27T21:55:01Z

valhalla/midgard/point_tile_index.h

+ */
+class PointTileIndex {
+  // This guy can tell us which tile a points belongs to.
+  std::unique_ptr<Tiles<PointLL>> tiles;


Switched to using the Tiles class for creating a gridded space over the given points-of-interest. Works well and performs equal to my previous solution.

merkispavel

LGTM

merkispavel · 2021-04-28T14:01:15Z

valhalla/midgard/polyline2.h

@@ -6,6 +6,7 @@
 #include <valhalla/midgard/point2.h>
 #include <valhalla/midgard/pointll.h>

+#include <list>


nit: unused

merkispavel · 2021-04-28T15:08:39Z

src/midgard/point_tile_index.cc

+  for (auto iter = polyline.begin(); iter != polyline.end(); iter++) {
+    const PointLL& p = *iter;


nit:

Suggested change

for (auto iter = polyline.begin(); iter != polyline.end(); iter++) {

const PointLL& p = *iter;

for (const auto& p : polyline) {

merkispavel · 2021-04-28T15:10:16Z

src/midgard/point_tile_index.cc

+    if (p.lat() < min_lat)
+      min_lat = p.lat();
+    if (p.lat() > max_lat)
+      max_lat = p.lat();
+    if (p.lng() < min_lng)
+      min_lng = p.lng();
+    if (p.lng() > max_lng)
+      max_lng = p.lng();


nit

Suggested change

if (p.lat() < min_lat)

min_lat = p.lat();

if (p.lat() > max_lat)

max_lat = p.lat();

if (p.lng() < min_lng)

min_lng = p.lng();

if (p.lng() > max_lng)

max_lng = p.lng();

min_lat = std::min(min_lat, p.lat());

max_lat = std::min(max_lat, p.lat());

min_lng = std::min(min_lng, p.lng());

max_lng = std::max(max_lng, p.lng());

merkispavel · 2021-04-28T21:04:38Z

src/midgard/CMakeLists.txt

@@ -1,9 +1,11 @@
 file(GLOB headers ${VALHALLA_SOURCE_DIR}/valhalla/midgard/*.h)

-set_source_files_properties(
+if ((UNIX OR APPLE) AND ENABLE_SINGLE_FILES_WERROR)


I believe APPLE is unix https://cmake.org/cmake/help/v3.6/variable/UNIX.html

merkispavel

Re-approve

ktatterso added 14 commits April 12, 2021 13:12

Whole bunch of early ideas on our Normalize() logic avoiding self-int…

8c70eca

…ersections.

Cleaning up, performance tweaks, generalizing.

e20f28f

More comments.

6a03e03

Modularizing logic. Templating.

17be104

Handling tile-id under/overflow. Fixing existing isochrone simplifica…

be8eea3

…tion test cases that were self intersecting. Lots of cleanup.

Updating existing isochrone tests which were generating self-intersec…

ae626ca

…tions to have new results without self-intersections. Much tidying.

Added triangle containment test. Moved Generalize tests into polyline…

5b06d22

…2 tests.

formatting

4589995

Merge remote-tracking branch 'origin/master' into intersecting-polys

94c2fa8

Tidying.

b49a882

Adding basic unit tests for PointTileIndex.

6f9328f

Merge remote-tracking branch 'origin/master' into intersecting-polys

df91c33

formatting

ed74fde

changelog update

2f3e2d7

ktatterso requested review from mandeepsandhu, danpat, merkispavel, purew, gknisely and mookerji April 21, 2021 23:13

merkispavel requested changes Apr 22, 2021

View reviewed changes

gknisely reviewed Apr 22, 2021

View reviewed changes

src/midgard/polyline2.cc Outdated Show resolved Hide resolved

gknisely reviewed Apr 22, 2021

View reviewed changes

valhalla/midgard/point_tile_index.h Outdated Show resolved Hide resolved

ktatterso added 4 commits April 27, 2021 11:38

Switching from my tile implementation to Tiles.

7fe2905

More cleanup.

dc3f27f

Placing an upper limit on the number of divisions of the tile space so

2cbff32

as not to overflow the underlying Tile id (int32_t).

Merge remote-tracking branch 'origin/master' into intersecting-polys

92de0e8

ktatterso commented Apr 27, 2021

View reviewed changes

Adding include for <memory>.

81ad33a

merkispavel previously approved these changes Apr 28, 2021

View reviewed changes

Sigh, moved the #include <memory> to point_tile_index.h.

bb7895d

ktatterso dismissed merkispavel’s stale review via bb7895d April 28, 2021 15:40

ktatterso added 2 commits April 28, 2021 10:10

Initializing a variable due to stricter compiler check.

5f9d037

Conditinalizing -Werror to Unix and Apple targets only.

0ae129d

ktatterso requested a review from merkispavel April 28, 2021 19:36

merkispavel reviewed Apr 28, 2021

View reviewed changes

merkispavel approved these changes Apr 28, 2021

View reviewed changes

ktatterso merged commit 3af3781 into master Apr 28, 2021

ktatterso deleted the intersecting-polys branch April 28, 2021 22:28

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Avoid generating self-intersecting isochrone polygons that would otherwise result from ::Generalize() #3026

Avoid generating self-intersecting isochrone polygons that would otherwise result from ::Generalize() #3026

ktatterso commented Apr 21, 2021 •

edited

Loading

gknisely commented Apr 22, 2021

merkispavel left a comment •

edited

Loading

ktatterso commented Apr 22, 2021

ktatterso Apr 27, 2021

merkispavel left a comment

merkispavel Apr 28, 2021

merkispavel Apr 28, 2021

merkispavel Apr 28, 2021

merkispavel Apr 28, 2021

merkispavel left a comment

		for (auto iter = polyline.begin(); iter != polyline.end(); iter++) {
		const PointLL& p = *iter;

	for (auto iter = polyline.begin(); iter != polyline.end(); iter++) {
	const PointLL& p = *iter;
	for (const auto& p : polyline) {

Avoid generating self-intersecting isochrone polygons that would otherwise result from ::Generalize() #3026

Avoid generating self-intersecting isochrone polygons that would otherwise result from ::Generalize() #3026

Conversation

ktatterso commented Apr 21, 2021 • edited Loading

Issue

Example 1

Example 2

Performance

Thoughts & other ideas

Tasklist

gknisely commented Apr 22, 2021

merkispavel left a comment • edited Loading

Choose a reason for hiding this comment

ktatterso commented Apr 22, 2021

ktatterso Apr 27, 2021

Choose a reason for hiding this comment

merkispavel left a comment

Choose a reason for hiding this comment

merkispavel Apr 28, 2021

Choose a reason for hiding this comment

merkispavel Apr 28, 2021

Choose a reason for hiding this comment

merkispavel Apr 28, 2021

Choose a reason for hiding this comment

merkispavel Apr 28, 2021

Choose a reason for hiding this comment

merkispavel left a comment

Choose a reason for hiding this comment

ktatterso commented Apr 21, 2021 •

edited

Loading

merkispavel left a comment •

edited

Loading