[Boost.Geometry] - Union generates polygon with punctures

[Olivier Heriveaux]

Hello !

I'm trying to use Boost.Geometry library to implement Minkowsky
algorithm (i.e. polygon convolution). I try to implement it in the
same way Boost.Polygon does, but during the algorithm execution, a
boost::geometry::overlay_invalid_input_exception is raised.

It seems that one of the input polygons to the union_ function is
invalid. To understand my mistake, I outputed the operand polygons to
SVG files and I discovered punctures in one of them. I believe
Boost.Geometry thinks my polygon is self-intersecting (may be due to
precision error with floats ?) - but this polygon itself is the result
of a previous union operation.

For info, I attached the SVG file (use inkscape to see it completely
since segments are out of the box due to negative coordinates).

I don't know if solutions to this problem exists or if the library is
faulty (the bad polygon is the result of an union operation). One
workaround would be to detect punctures and remove them, but it's kind
of ugly imho...

Some more details:
- point coordinates type is float,
- I use custom class instead of point_xy model,
- polygon winding is counter-clockwise.

Thank you in advance for any help,
Olivier

[Barend Gehrels]

Hi,

Post by Olivier Heriveaux
Hello !
I'm trying to use Boost.Geometry library to implement Minkowsky
algorithm (i.e. polygon convolution). I try to implement it in the
same way Boost.Polygon does, but during the algorithm execution, a
boost::geometry::overlay_invalid_input_exception is raised.
It seems that one of the input polygons to the union_ function is
invalid. To understand my mistake, I outputed the operand polygons to
SVG files and I discovered punctures in one of them. I believe
Boost.Geometry thinks my polygon is self-intersecting (may be due to
precision error with floats ?) - but this polygon itself is the result
of a previous union operation.
For info, I attached the SVG file (use inkscape to see it completely
since segments are out of the box due to negative coordinates).
I don't know if solutions to this problem exists or if the library is
faulty (the bad polygon is the result of an union operation). One
workaround would be to detect punctures and remove them, but it's kind
of ugly imho...
- point coordinates type is float,
- I use custom class instead of point_xy model,
- polygon winding is counter-clockwise.
Thank you in advance for any help,
Olivier

Thanks for the report - can you give more details, especially WKT of the
two input polygons? Or the sequence of algorithms you are using? These
spikes normally never occur in unions, however in difference they might
occur in some circumstances. Is this also happening with double typed
data (in case you know)?

Regards, Barend

[Olivier Heriveaux]

Hi,

Thanks for your support !

I exported the two inputs which produce the error during the algorithm
execution and I wrote a simple program appart with thoose polygons. I also
simplified a little the inputs and i was able to reproduce the problem (see
attached source file main.cpp).

The toy program has two input polygons and generates the union. Here are
the two input polygons:

boost::geometry::read_wkt
(
"POLYGON(("
"-0.973360121250153 0.194701105356216,"
"-0.923463344573975 0.184775948524475,"
"-0.858578622341156 0.141421362757683,"
"0.18477588891983 2.07653665542603,"
"0.141421347856522 2.1414213180542,"
"-1.18477594852448 -0.076536700129509,"
"-1.1414213180542 -0.141421377658844,"
"-0.973360121250153 0.194701105356216"
"))",
pa
);

boost::geometry::read_wkt
(
"POLYGON(("
"-0.14142133295536 1.8585786819458,"
"-1.1414213180542 -0.141421377658844,"
"-1.07653665542603 -0.184775933623314,"
"-0.076536625623703 1.81522405147552,"
"-0.14142133295536 1.8585786819458"
"))",
pb
);

And here is the output polygon, printed to stdout:

-0.19351507723331451416 1.5812671184539794922
-0.07653662562370300293 1.8152240514755249023
-0.1414213329553604126 1.8585786819458007812
-0.72200763225555419922 0.69740593433380126953
-1.1847759485244750977 -0.076536700129508972168
-1.1414213180541992188 -0.14142137765884399414 <- A
-0.9733600616455078125 0.19470109045505523682 <- here is the puncture
-1.1414213180541992188 -0.14142137765884399414 <- point equal to A
-1.0765366554260253906 -0.18477593362331390381
-0.89969980716705322266 0.1688976585865020752
-0.85857862234115600586 0.14142136275768280029
0.18477588891983032227 2.0765366554260253906
0.14142134785652160645 2.1414213180541992188
-0.19351507723331451416 1.5812671184539794922

I'll attach the generated SVG for inputs and output.
As you suggested, I also tried to use double instead of float. The result
is different (note: not all digits printed for doubles) :

-0.19351506716472754999 1.5812671771090451855
-0.07653662562370300293 1.8152240514755200174
-0.14142133295535999626 1.8585786819457998931
-0.72200777157880025037 0.69740575279044780821
-1.1847759485244799826 -0.076536700129508999924
-1.1414213180542001069 -0.14142137765884399414 <- A
-0.97336012125015303198 0.19470110535621598657 <- one new extra-point
-0.97336008742385127235 0.19470109862769818809 <-
-1.1414213180542001069 -0.14142137765884399414 <- point equal to A (I
checked till last digit, not shown here)
-1.0765366554260300536 -0.18477593362331398708
-0.89969984067827302177 0.16889768269682844948
-0.85857862234115600586 0.14142136275768299458
0.1847758889198299892 2.0765366554260298315
0.14142134785652199502 2.1414213180542001069
-0.19351506716472754999 1.5812671771090451855

With double precision, a kind of hole is created instead of a puncture, but
I don't know if this should be accepted as a good polygon...

One more thing I forgot to mention: I'm using boost 1.53 - I see another
version has been released very recently and I think i'll give a try with it
asap.

I'm not very familiar with Boost.Geometry, so maybe I'm doing something
wrong or maybe I'm expecting too much.

Regards

Post by Barend Gehrels
Hi,

Post by Olivier Heriveaux
Hello !
I'm trying to use Boost.Geometry library to implement Minkowsky
algorithm (i.e. polygon convolution). I try to implement it in the
same way Boost.Polygon does, but during the algorithm execution, a
boost::geometry::overlay_**invalid_input_exception is raised.
It seems that one of the input polygons to the union_ function is
invalid. To understand my mistake, I outputed the operand polygons to
SVG files and I discovered punctures in one of them. I believe
Boost.Geometry thinks my polygon is self-intersecting (may be due to
precision error with floats ?) - but this polygon itself is the result
of a previous union operation.
For info, I attached the SVG file (use inkscape to see it completely
since segments are out of the box due to negative coordinates).
I don't know if solutions to this problem exists or if the library is
faulty (the bad polygon is the result of an union operation). One
workaround would be to detect punctures and remove them, but it's kind
of ugly imho...
- point coordinates type is float,
- I use custom class instead of point_xy model,
- polygon winding is counter-clockwise.
Thank you in advance for any help,
Olivier

Thanks for the report - can you give more details, especially WKT of the
two input polygons? Or the sequence of algorithms you are using? These
spikes normally never occur in unions, however in difference they might
occur in some circumstances. Is this also happening with double typed data
(in case you know)?
Regards, Barend
______________________________**_________________
Boost-users mailing list
http://lists.boost.org/**mailman/listinfo.cgi/boost-**users<http://lists.boost.org/mailman/listinfo.cgi/boost-users>

[Barend Gehrels]

Hi Olivier,


The policy of this list is: don't top-post. I reordered your mail to
read from top to bottom, my answers in between.

Post by Olivier Heriveaux
It seems that one of the input polygons to the union_ function is
invalid. To understand my mistake, I outputed the operand polygons to
SVG files and I discovered punctures in one of them. I believe
Boost.Geometry thinks my polygon is self-intersecting (may be due to
precision error with floats ?) - but this polygon itself is the result
of a previous union operation.
I don't know if solutions to this problem exists or if the library is
faulty (the bad polygon is the result of an union operation). One
workaround would be to detect punctures and remove them, but it's kind
of ugly imho...
- point coordinates type is float,
- I use custom class instead of point_xy model,
- polygon winding is counter-clockwise.
Thank you in advance for any help,
Olivier
I exported the two inputs which produce the error during the algorithm
execution and I wrote a simple program appart with thoose polygons. I
also simplified a little the inputs and i was able to reproduce the
problem (see attached source file main.cpp).

Thanks for your sample program and your tests.

Post by Olivier Heriveaux
I'll attach the generated SVG for inputs and output.
As you suggested, I also tried to use double instead of float. The
-1.1414213180542001069 -0.14142137765884399414 <- A
-0.97336012125015303198 0.19470110535621598657 <- one new extra-point
-0.97336008742385127235 0.19470109862769818809 <-
-1.1414213180542001069 -0.14142137765884399414 <- point equal to A
(I checked till last digit, not shown here)

My version with high precision (ttmath):
-1.1414213180542 -0.141421377658844 (point of both pa and pb)
-0.973360121250153 0.194701105356216 (point of pa)
-0.9733600874238512128115932867815276942
0.19470109862769819504261957771496882692 (generated point)
-1.1414213180542 -0.141421377658844 (point of both pa and pb)

Post by Olivier Heriveaux
One more thing I forgot to mention: I'm using boost 1.53 - I see
another version has been released very recently and I think i'll give
a try with it asap.

This behaviour is not changed.

Post by Olivier Heriveaux
I'm not very familiar with Boost.Geometry, so maybe I'm doing
something wrong or maybe I'm expecting too much.

I tried it with float, double and ttmath and can reproduce your case
with all of them. ttmath is expected to be high precision so it is
(probably) not really a floating point error. I also tried to make the
union with both SQL Server and PostGIS and both give a similar output,
including the pseudo-spike. SQL Server makes it one polygon, PostGIS
makes it a polygon with a hole. So this case is indeed quite
interesting. They both report the output polygon to be valid. A polygon
may contain a self-tangency.

The whole-like triangle has an area of
-0.000000006250291651472502798247346743117038847 so it is in fact not a
spike, at least not in ttmath (and according to your report, not in double).

So only in float we have the real spike. The output is indeed not
checked on these cases. They should not be generated but in float they
apparently are, because of the lower resolution.

Did you try the complete program (so including the steps after this
union) with double? Is it then also reporting the self-intersection? You
probably cannot use float for this case... We can filter out the spike
in a post-process, as you suggest, but it is (at least not now) the
intention to do this by default.

Regards, Barend

[Olivier Heriveaux]

The policy of this list is: don't top-post. I reordered your mail to read
from top to bottom, my answers in between.

Sorry for the inconvenient, I knew it but gmail hid your answer in my reply
form (have to click a button to show it)...

So only in float we have the real spike. The output is indeed not checked
on these cases. They should not be generated but in float they apparently
are, because of the lower resolution.

So this could really be considered as an issue ? Do you think a solution
could correct this behaviour or must we deal with it ? I'm not asking for a
fix, I currently use Boost.Polygon to achieve what I want - just wanted to
test it with Boost.Geometry since it can work with float (I currently
convert my float space to int space to work with Boost.Polygon). As I found
this corner case I thought it could be usefull to report it.

Did you try the complete program (so including the steps after this union)
with double? Is it then also reporting the self-intersection? You probably
cannot use float for this case... We can filter out the spike in a
post-process, as you suggest, but it is (at least not now) the intention to
do this by default.

It requires many changes in my program, it is a currently in a big project
(point type is not a template since I knew I wanted floating numbers). I
will try to extract the entire algorithm and make a simple test program for
it.

Regards, Olivier.

[Barend Gehrels]

Hi Olivier,

Post by Barend Gehrels
So only in float we have the real spike. The output is indeed not
checked on these cases. They should not be generated but in float
they apparently are, because of the lower resolution.
So this could really be considered as an issue ? Do you think a
solution could correct this behaviour or must we deal with it ? I'm
not asking for a fix, I currently use Boost.Polygon to achieve what I
want - just wanted to test it with Boost.Geometry since it can work
with float (I currently convert my float space to int space to work
with Boost.Polygon). As I found this corner case I thought it could be
usefull to report it.

Sure, your report is very useful.
Yes it is an issue, a spike should not be there. I will think about how
it can be avoided efficiently.

Besides that, another solution (which was already planned for the cases
of Volker Schöch) is that the intersection can deal with spikes. The
union would then still have a spike, but algorithms later on do not
throw on that. I have to do more research on this.

Post by Barend Gehrels
Did you try the complete program (so including the steps after
this union) with double? Is it then also reporting the
self-intersection? You probably cannot use float for this case...
We can filter out the spike in a post-process, as you suggest, but
it is (at least not now) the intention to do this by default.
It requires many changes in my program, it is a currently in a big
project (point type is not a template since I knew I wanted floating
numbers). I will try to extract the entire algorithm and make a simple
test program for it.

OK, thanks.

Regards, Barend