Problems with RGeo Intersection Functions

I get weird / wrong results from the intersection functions of RGeo polygons (Ruby 2.3.0, RGeo 0.5.3)

EXAMPLE 1:

I have two polygons that, I believe, share a border, but do not share any internal space (i.e. they are touch , but not overlap ):

wkt_1 = "POLYGON ((-8226874.27782158 4962626.76394919, -8223358.174520462 4961756.817075645, -8223358.174520462 4960289.557693501, -8224471.369428394 4960289.557693501, -8226874.27782158 4962253.674727506, -8226874.27782158 4962626.76394919))"
wkt_2 = "POLYGON ((-8224757.546680832 4960523.476563589, -8225269.1002275925 4959296.105368667, -8226993.791361805 4959219.668340384, -8226420.900079966 4961883.087589158, -8224757.546680832 4960523.476563589))"
poly_1 = RGeo::Geos.factory(:srid => 3857).parse_wkt(wkt_1)
poly_2 = RGeo::Geos.factory(:srid => 3857).parse_wkt(wkt_2)

When we check the intersection between them, it returns a string, as expected, for the geometry that separates only the border:

poly_1.intersection poly_2
=> #<RGeo::Geos::CAPILineStringImpl:0x3fc0249af168 "LINESTRING (-8224757.546680832 4960523.476563589, -8225598.074380083 4961210.51680879)">

However, when performing the following checks, we get the opposite of what was expected:

poly_1.overlaps? poly_2
=> true
poly_1.touches? poly_2
=> false

EXAMPLE 2:

We take two legally overlapping polygons:

wkt_3 = "POLYGON ((-8243237.0 4970203.0, -8243237.0 4968735.0, -8242123.0 4968735.0, -8242123.0 4970203.0, -8243237.0 4970203.0))"
wkt_4 = "POLYGON ((-8244765.0 4966076.0, -8244765.0 4964608.0, -8243652.0 4964608.0, -8243652.0 4966076.0, -8242680.0 4969362.0, -8244765.0 4966076.0))"
poly_3 = RGeo::Geos.factory(:srid => 3857).parse_wkt(wkt_3)
poly_4 = RGeo::Geos.factory(:srid => 3857).parse_wkt(wkt_4)

enter image description here

And calculate the differences in both directions:

diff_a = poly_3.difference poly_4
=> #<RGeo::Geos::CAPIPolygonImpl:0x3fe3fca26028 "POLYGON ((-8243077.837796713 4968735.0, -8243237.0 4968735.0, -8243237.0 4970203.0, -8242123.0 4970203.0, -8242123.0 4968735.0, -8242865.466828971 4968735.0, -8242680.0 4969362.0, -8243077.837796713 4968735.0))">
diff_b = poly_4.difference poly_3
=> #<RGeo::Geos::CAPIPolygonImpl:0x3fe3fd1dda28 "POLYGON ((-8242865.466828971 4968735.0, -8243652.0 4966076.0, -8243652.0 4964608.0, -8244765.0 4964608.0, -8244765.0 4966076.0, -8243077.837796713 4968735.0, -8242865.466828971 4968735.0))">

Now we check the residual polygons for their subtractors:

diff_b.touches? poly_3
=> true
diff_b.overlaps? poly_3
=> false

It's good.

diff_a.touches? poly_4
=> false
diff_a.overlaps? poly_4
=> true

. , , !

?

3:

. poly_3 poly_4

intersection_a = poly_3.intersection poly_4
=> #<RGeo::Geos::CAPIPolygonImpl:0x3fd1084ece88 "POLYGON ((-8242865.724520766 4968734.582337743, -8243078.32501591 4968734.582337743, -8242680.062418439 4969362.301390027, -8242865.724520766 4968734.582337743))">

, , poly_3, diff_a, diff_a intersection_a , poly_4 ().

:

diff_a.touches? poly_4
=> false
diff_a.touches? intersection_a
=> true
diff_a.intersection poly_4
=> #<RGeo::Geos::CAPILineStringImpl:0x3fe3f98fb854 "LINESTRING (-8242680.0 4969362.0, -8243077.837796713 4968735.0)">
diff_a.intersection intersection_a
=> #<RGeo::Geos::CAPIMultiLineStringImpl:0x3fe3fca157b4 "MULTILINESTRING ((-8242865.466828971 4968735.0, -8242680.0 4969362.0), (-8242680.0 4969362.0, -8243077.837796713 4968735.0))">

, . , 2- , .

+4
2

, , touches? overlaps? Float.

RGeo GEOS (, , JTS, , GEOS).

, Geometry#distance , epsilon.

Geometry#intersection , touches? overlaps?, .


?

touches? : , .

  • , - , touches? - false.
  • , , touches? .

?

, :

enter image description here

POLYGON((0 0, 0 1, 1 1, 1 0, 0 0))
POLYGON((0 1, 0 2, 1 2, 1 1, 0 1))

, :

require 'rgeo'

epsilon = 1E-15

deltas = [-epsilon, 0, epsilon]

deltas.each do |delta|
  puts "--------------------------------"
  puts "Delta : #{delta}\n\n"
  simple1 = 'POLYGON((0 0, 0 1, 1 1, 1 0, 0 0))'
  simple2 = "POLYGON((0 #{1+delta}, 0 2, 1 2, 1 #{1+delta}, 0 #{1+delta}))"

  puts "Polygon #1 : #{simple1}\n"
  puts "Polygon #2 : #{simple2}\n\n"

  poly_1 = RGeo::Geos.factory(:srid => 4326).parse_wkt(simple1)
  poly_2 = RGeo::Geos.factory(:srid => 4326).parse_wkt(simple2)

  puts "Intersection : #{poly_1.intersection poly_2}"
  puts

  puts "Distance between polygons :"
  puts format('%.30f°', poly_1.distance(poly_2))
  puts

  puts "Overlaps? : #{poly_1.overlaps? poly_2}"
  puts "Touches? : #{poly_1.touches? poly_2}"
end

:

--------------------------------
Delta : -1.0e-15

Polygon #1 : POLYGON((0 0, 0 1, 1 1, 1 0, 0 0))
Polygon #2 : POLYGON((0 0.999999999999999, 0 2, 1 2, 1 0.999999999999999, 0 0.999999999999999))

Intersection : POLYGON ((0.0 0.999999999999999, 0.0 1.0, 1.0 1.0, 1.0 0.999999999999999, 0.0 0.999999999999999))

Distance between polygons :
0.000000000000000000000000000000°

Overlaps? : true
Touches? : false
--------------------------------
Delta : 0

Polygon #1 : POLYGON((0 0, 0 1, 1 1, 1 0, 0 0))
Polygon #2 : POLYGON((0 1, 0 2, 1 2, 1 1, 0 1))

Intersection : LINESTRING (0.0 1.0, 1.0 1.0)

Distance between polygons :
0.000000000000000000000000000000°

Overlaps? : false
Touches? : true
--------------------------------
Delta : 1.0e-15

Polygon #1 : POLYGON((0 0, 0 1, 1 1, 1 0, 0 0))
Polygon #2 : POLYGON((0 1.000000000000001, 0 2, 1 2, 1 1.000000000000001, 0 1.000000000000001))

Intersection : GEOMETRYCOLLECTION EMPTY

Distance between polygons :
0.000000000000001110223024625157°

Overlaps? : false
Touches? : false

  • ,

1E-15 ( 1 Ångström) , , overlaps? touches? true false.

, , , , , intersection, overlaps? touches?.

. 1-Ångström!

RGeo?

RGeo GEOS, ++- JTS (Java Topology Suite).

, RGeo GEOS, 1 JTS 1.14:

WKTReader wktReader = new WKTReader();
String wkt1 = "POLYGON ((-8226874.27782158 4962626.76394919, -8223358.174520462 4961756.817075645, -8223358.174520462 4960289.557693501, -8224471.369428394 4960289.557693501, -8226874.27782158 4962253.674727506, -8226874.27782158 4962626.76394919))";
String wkt2 = "POLYGON ((-8224757.546680832 4960523.476563589, -8225269.1002275925 4959296.105368667, -8226993.791361805 4959219.668340384, -8226420.900079966 4961883.087589158, -8224757.546680832 4960523.476563589))";
Polygon poly1 = (Polygon) wktReader.read(wkt1);
Polygon poly2 = (Polygon) wktReader.read(wkt2);
System.out.println("Intersection : " + poly1.intersection(poly2));
System.out.println("Overlaps?    : " + poly1.overlaps(poly2));
System.out.println("Intersects?  : " + poly1.intersects(poly2));
System.out.println("Touches?     : " + poly1.touches(poly2));
showMatrixWith(poly1, poly2);

:

Intersection : LINESTRING (-8224757.546680832 4960523.476563589, -8225598.074380083 4961210.51680879)
Overlaps?    : true
Intersects?  : true
Touches?     : false

212
101
212

, , intersects? touches?   , RGeo.

?

intersection touches? ?

DE-9IM

touches?, intersects?, overlaps? (DE-9IM). , , .

72 src/operation/relate/RelateComputer.cpp GEOS:

IntersectionMatrix* RelateComputer::computeIM()

, noding, .

, JTS, , , " , ":

  # line 477 in jts-1.14/testxml/general/TestFunctionAA.xml
  <desc>mAmA - complex polygons touching</desc>
  <a>
    MULTIPOLYGON(
      (
        (100 200, 100 180, 120 180, 120 200, 100 200)),
      (
        (60 240, 60 140, 220 140, 220 160, 160 160, 160 180, 200 180, 200 200, 160 200,
        160 220, 220 220, 220 240, 60 240),
        (80 220, 80 160, 140 160, 140 220, 80 220)),
      (
        (280 220, 240 180, 260 160, 300 200, 280 220)))
  </a>
  <b>
    MULTIPOLYGON(
      (
        (80 220, 80 160, 140 160, 140 220, 80 220),
        (100 200, 100 180, 120 180, 120 200, 100 200)),
      (
        (220 240, 220 220, 160 220, 160 200, 220 200, 220 180, 160 180, 160 160, 220 160,
        220 140, 320 140, 320 240, 220 240),
        (240 220, 240 160, 300 160, 300 220, 240 220)))
  </b>

, .

, wkt_3 wkt_4 , : x1 of diff_a is -8243077.837796713.

#

Geometry#intersection 670 src/operation/overlay/OverlayOp.cpp GEOS:

void OverlayOp::computeOverlay(OverlayOp::OpCode opCode)

, -, , , if, .

, RelateComputer::computeIM.

GeometryPrecisionReducer?

GeometryPrecisionReducer a PrecisionModel .

GeometryPrecisionReducer GEOS, RGeo. JTS , : . // , .

. PrecisionModel , .

+2

. , ? ? ( ), . , .

. . , 2D-, . , , , . , , , .

RGeo, ? ( ). , , .

0

Source: https://habr.com/ru/post/1665021/


All Articles