Interior straight skeleton leaks on the exterior with inexact constructions

[strk]

As reported in #42 (comment) the interior straight skeleton of polygon can end up being outside of the input polygon due to robustness issues when created with exact predicates and inexact constructions kernel.

I've reduced the testcase to this (in WKT form):

POLYGON((1687891.03971877 5401915.74731365,1684030.810026 5402682.520688,
1684000 5402672.31399816,1681298.83486802 5401963.62058474,
1681281.134295 5403807.881624,1683915.838219 5402676.741864,
1684000 5402711.41802501,1684020.618162 5402719.913077,
1687892.11090561 5404508.53561586,1687891.03971877 5401915.74731365))

The artifact happens near segments of the boundary which are relatively short compared to the overall bounding box.

Running the skeleton builder with exact constructions kernel gives the correct answer.
Is it expected ?

[strk]

I've noticed that translating the input polygon to have it's lower-left point at 0,0 before calling the straight skeleton computer fixes the issue, which clearly identifies the problem as a robustness one.

[strk]

I realized that constructing the input polygon from WKT does not reproduce the problem due to roundup. This WKB input (in HEX form) has all the precision it takes to trigger it:

010300002091080000010000000A00000063022B0A53C139419EFCD3EF4E9B54412DDD5DCF3EB23941C4F352A10E9C54410400000020B23941C08B18140C9C544116E9B9D592A739410DA9B7E75A9B54413B28612281A7394113876CF8279D54413C8595D6CBB1394125B37A2F0D9C54410400000020B23941F9EBC0DA159C544168DD3F9E34B2394184DA6FFA179C5441524F641C54C13941C3874722D79D544163022B0A53C139419EFCD3EF4E9B5441

See https://en.wikipedia.org/wiki/Well-known_text#Well-known_binary for interpretation of the WKB

[lrineau]

Fernando? @fcacciola

[fcacciola]

Hello. Indeed, the results you're seeing are to be expected when using inexact constructions.
Have in mind that a Straight Skeleton has both topological and geometrical components. The topological component is the graph of nodes that link all skeleton vertices. Using exact predicates guarantees that graph is correct. The geometrical component is the actual coordinates of skeleton vertices in the output. You can only guarantee correctly rounded vertex coordinates by using exact constructions, as you discovered.
But as you also discovered, unfortunately, exact constructions take up way much more time, so they are better avoided whenever possible.

Robustness issues are always directly dependent on the concrete numerical values involved in the computations, which is why you discovered that a "perturbation" of the input polygon (the translation) changes things. Unfortunately, that change doesn't mean the problem goes away, just that it is hidden on this test case, so it won't work a general solution.
As you apparently figured out, the translation to the origin reduces the effective magnitud of the numerical values involved, reducing the chances of an incorrect output, but it won't work in most cases.

The input modification which does do the trick in almost all cases is, as vmora found on Oslandia/SFCGAL#89, to merge almost collinear edges with a proper tolerance. So is merging very very small edges as you found here.

So, the general solution is not to translate (nor rotate nor scale) the input, but to "clean" it, removing almost collinear and extremely tiny edges.

In fact, in ALL practical cases, I always run a polyline simplification (for both collinearity and zero-length edges) with an application-specific threshold (which is usually MUCH larger than epsion) on every input to the straight skeleton. This not only significantly reduces the chances of incorrect vertex coordinates in the output (when using inexact constructions), but also reduces the running time by several factors in my real cases since it usually happens that lots of almost collinear vertices are found in the input due to curve sampling of the "real" input.

[lrineau]

@strk I want to point out that your bug report was well documented. Thanks.

After the explanation from Fernando, do you accept that the bug is closed?

[strk]

I understand it is by policy that invalid results are accepted with "inexact constructions" kernel. If you confirm, ok with me to close the bug. For comparison, in GEOS (http://github.com/libgeos/libgeos) robustness issues are expected to throw an exception (at the expense of validation steps).

[strk]

Oh, and thanks a lot Fernando for your exhaustive explanation !