我从rgeo polygon交集函数(ruby 2.3.0,rgeo 0.5.3)得到了奇怪/不正确的结果。
例1:
我有两个多边形,我相信它们共享一个边界,但不共享任何内部空间(即它们touch但不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)

当我们检查它们之间的交集时,它将返回一条直线,正如对于只共享边界的几何图形所预期的那样:
poly_1.intersection poly_2
=> #<RGeo::Geos::CAPILineStringImpl:0x3fc0249af168 "LINESTRING (-8224757.546680832 4960523.476563589, -8225598.074380083 4961210.51680879)">

但是,当运行以下检查时,我们得到的结果与预期相反:
poly_1.overlaps? poly_2
=> true
poly_1.touches? poly_2
=> false

例2:
我们取两个合法重叠的多边形:
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)

ruby - RGeo交点函数的问题-LMLPHP
并计算两个方向的差异:
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))">

现在,我们将剩余多边形与它们的减法器进行比较:
diff_b.touches? poly_3
=> true
diff_b.overlaps? poly_3
=> false

这很好。
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))">

更糟糕的是,这两个相交结果都没有意义。它应该是一条单独的、两段的线,而这两条线都不是。

最佳答案

简短回答
不幸的是,当使用浮点坐标时,seems不能期望从touches?overlaps?获得可靠和正确的输出。
它不依赖于rgeo或geos版本(或者就此而言,geos所基于的项目JTS)。
如果需要有关两个多边形位置的信息,可以使用Geometry#distance并检查它是否小于给定的epsilon。
Geometry#intersectiontouches?overlaps?更可靠,但也不能保证每个示例都适用。
是精度问题吗?
理论
touches?的定义非常敏感:多边形必须在其边界上共享一个点或一条线,但不应该有任何公共的内部点。
如果它们彼此距离太远,它们就没有任何共同点,touches?是错误的。
如果它们太近,则它们的交集为多边形,touches?为假。
敏感性分析
它到底有多敏感?
让我们考虑两个多边形,一个刚好在另一个上面:
ruby - RGeo交点函数的问题-LMLPHP

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

我们可以用epsilon移动上多边形的下边界,看看它有什么影响:
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?都在truefalse之间切换。
交集可以是空的、直线或多边形,但至少结果在intersectionoverlaps?touches?之间似乎是一致的。
在你的例子中,有更复杂的多边形和倾斜的边使问题更难解决。在计算两个线段的交集时,很难保持1_ngstróm的精度!
rgeo有问题吗?
RGeo使用GEOS,这是JTS (Java Topology Suite)的C++端口。
为了检查问题不是rgeo或geos特有的,我计算了
JTS 1.14示例1:
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相同的错误结果。
为什么结果不一致?
为什么intersectiontouches?返回不一致的结果?
DE-9IM
touches?intersects?overlaps?和其他谓词派生自the Dimensionally Extended nine-Intersection Model (DE-9IM)。它是一个矩阵,描述几何图形的内部、边界和外部之间的交集维数。
该矩阵在GEOS中src/operation/relate/RelateComputer.cpp的第72行计算:
IntersectionMatrix* RelateComputer::computeIM()

该算法似乎require exact 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_3wkt_4有圆形坐标,计算它们的差异也会创建具有不精确坐标的多边形:x1diff_a-8243077.837796713
几何交集
Geometry#intersection是在geos中src/operation/overlay/OverlayOp.cpp的第670行计算的:
void OverlayOp::computeOverlay(OverlayOp::OpCode opCode)

代码中的注释似乎表明不需要确切的点头,并且有多个if语句来检查结果是否正确。
这种方法似乎比RelateComputer::computeIM更可靠。
使用GeometryPrecisionReducer?
使用GeometryPrecisionReducerPrecisionModel,可以为所有几何图形指定允许点的网格。
GeometryPrecisionReducer在geos中实现,但在rgeo中不可用。
在jts中使用第一个示例进行的测试表明,它并不能解决您的问题:不精确的坐标会捕捉到网格上最近的点。每一个点都会向北/南/东或西移动一点,这会改变每一侧的坡度。
它还改变了边界及其交点。根据PrecisionModel,第一个示例的交集可以是空的、直线或多边形。

关于ruby - RGeo交点函数的问题,我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/41349388/

10-13 07:45