difference of multilinestring and polygon yields wrong result

[Philipp Müller <pmueller@…>]

I am using the following typedefs to simplify reading:

typedef boost::geometry::model::point<int,2,boost::geometry::cs::cartesian> TPoint;
typedef boost::geometry::model::linestring<TPoint> TLinestring;
typedef boost::geometry::model::multi_linestring<TLinestring> TMultiLinestring;
typedef boost::geometry::model::polygon<TPoint, /*ClockWise*/false, /*Closed*/false> TPolygon;

I have a multilinestring from which I want to subtract two polygons. The multilinestring and polygons are as follows

TMultiLinestring multilinestring; 
boost::geometry::read_wkt("MULTILINESTRING((5239 2113,1020 2986))", multilinestring);
TPolygon polygon1; 
boost::geometry::read_wkt("POLYGON((5233 2113,5200 2205,1020 2205,1020 2022,5200 2022))", polygon1);
TPolygon polygon2; 
boost::geometry::read_wkt("POLYGON((5233 2986,5200 3078,1020 3078,1020 2895,5200 2895))", polygon2);

First, I subtract polygon1 from multilinestring, which works fine:

TMultiLinestring multilinestringOut1;
boost::geometry::difference(multilinestring, polygon1, multilinestringOut1);
// works as expected: multilinestringOut1 holds MULTILINESTRING((5239 2113,5233 2114)(4795 2205, 1020 2986))

Now I want to subtract polygon2 from multilinestringOut1 (the intermediate result from the previous calculation), which yields a wrong result.

TMultiLinestring multilinestringOut2;
boost::geometry::difference(multilinestringOut1, polygon2, multilinestringOut2); 
// not   as expected: multilinestringOut2 holds MULTILINESTRING((5239 2113,5233,2114)(4795 2205, 5406 1580))

I would expect multilinestringOut2 to hold

MULTILINESTRING((5239 2113,5233,2114)(4795 2205, 1460 2895))

[anonymous]

I re-checked this behaviour against an older version of boost and realized that the problem does not occur there.

[Philipp Müller <pmueller@…>]

It may help to look at line 187ff in cart_intersect.hpp. There, cramers_rule "returned" values that seem incorrect to me.

[anonymous]

I just checked out the latest trunk, where it seems to work as expected (multilinestringOut2 holds MULTILINESTRING((5239 2113,5232 2114),(4794 2205,1459 2895)) which differs from expectation probably only because of rounding errors).

[mkaravel]

The problem was in the computation of the intersection points. In particular, integer overflow was taking place, yielding the wrong results.

The problem has been fixed by locally promoting the coordinates to a number type with increased precision, thus avoiding the overflow.

[mkaravel]

See also tickets #10658 and #11121.