增强几何返回相交和相交的不一致结果

Posted 2023-02-17

【中文标题】增强几何返回相交和相交的不一致结果

【英文标题】：Boost geometry returning inconsistent results for intersects and intersection

【发布时间】：2021-02-22 08:06:35

【问题描述】：

在我的应用程序中，我主要使用 boost 几何来执行相交 和差异计算。不幸的是，我注意到不一致 在结果中：

bg::相交和 bg::intersection

我做的是：

计算多边形 1 和多边形 2 的交集（= 多边形 3）。 从polygon1 (= poly4) 中删除（bg::difference）此交集 生成的 polygon4 不应与 poly2 有任何交集。

#include <iostream>
#include <boost/geometry.hpp>

namespace bg = boost::geometry;
using point_t = bg::model::d2::point_xy<double>;
using polygon_t = bg::model::polygon<point_t>;
using mpolygon_t = bg::model::multi_polygon<polygon_t>;

int main()

    polygon_t poly1, poly2;
    mpolygon_t poly3, poly4, poly5;

    bg::read_wkt("POLYGON(("
        "12227.0 4967.0000000000009, 12238.0 4967.0000000000009, "
        "12238.0 4813.0000000000009, 12227.0 4813.0000000000009, "
        "12227.0 4967.0000000000009))", poly1);
    bg::read_wkt("POLYGON(("
        "12254.0 4947.0, 12219.0 4982.0, 12219.0 5020.0, 12254.0 5055.0, "
        "12261.0 5055.0, 12263.0 5055.0, 12283.0 5055.0, 12283.0 4947.0, "
        "12263.0 4947.0, 12261.0 4947.0, 12254.0 4947.0))", poly2);

    bg::intersection(poly1, poly2, poly3);
    bg::difference(poly1, poly3[0], poly4);

    // b0 = true, b1 = false
    bool b0 = bg::intersects(poly2, poly4[0]);
    bool b1 = bg::intersection(poly2, poly4, poly5) && (poly5.size() != 0);

    bool b2 = !bg::disjoint(poly2, poly4[0]) && !bg::touches(poly2, poly4[0]);
    bool b3 = bg::overlaps(poly2, poly4[0]) || bg::within(poly4[0], poly2) || bg::within(poly2, poly4[0]);

    std::cout << b0 << b1 << b2 << b3 << std::endl;
    return 1;

奇怪的是 bg::intersects 返回 true 而 bg::intersection 返回一个空的交集。

有人知道为什么会这样吗？ （也许准确性问题？）和 更有趣的是：如何避免此类问题？

我试图通过使用其他函数来避免“相交”，但结果没有帮助。请参阅 b2 和 b3 的计算。

PS：不幸的是，该示例似乎导致了 coliru 崩溃。

【问题讨论】：

主题：请问你是怎么做图表的？我想找到一个比使用 Boost Geometry 编写 SVG 更好的工作流程 - 因为每次都需要大量工作

当然：我使用 Visual Studio 的“图形调试”插件，这让这变得非常简单。 (marketplace.visualstudio.com/…)

啊。 5 年来我第一次听说有一个功能让我想再次尝试 VS :) 干杯。

【参考方案1】：

实际上intersection函数的返回值是未指定的：见

What does boost::geometry::intersection return

所以，你需要看看路口本身：

Unlive On Coliru

Unlive On Compiler Explorer

Live On Wandbox（耶！）

#include <iostream>
#include <boost/geometry.hpp>
#include <boost/logic/tribool.hpp>
#include <boost/logic/tribool_io.hpp>

using boost::tribool;
using boost::indeterminate;

namespace bg = boost::geometry;
using Coord = double;
using point_t = bg::model::d2::point_xy<Coord>;
using polygon_t = bg::model::polygon<point_t>;
using mpolygon_t = bg::model::multi_polygon<polygon_t>;

int main() 
    polygon_t a, b;

    bg::read_wkt("POLYGON(("
        "12227.0 4967.0000000000009, 12238.0 4967.0000000000009, "
        "12238.0 4813.0000000000009, 12227.0 4813.0000000000009, "
        "12227.0 4967.0000000000009))", a);
    bg::read_wkt("POLYGON(("
        "12254.0 4947.0, 12219.0 4982.0, 12219.0 5020.0, 12254.0 5055.0, "
        "12261.0 5055.0, 12263.0 5055.0, 12283.0 5055.0, 12283.0 4947.0, "
        "12263.0 4947.0, 12261.0 4947.0, 12254.0 4947.0))", b);

#define CHECK(expected, expr)                                                  \
    do                                                                        \
        auto const& actual = (expr);                                           \
        auto ok = ((expected) == actual);                                      \
        if (ok == true)                                                       \
            std::cout << "PASS";                                               \
         else if (ok == false)                                               \
            std::cout << "FAIL";                                               \
         else                                                                \
            std::cout << "?";                                                  \
                                                                              \
        std::cout << "\t" << #expr << " -> " << std::boolalpha << actual       \
                  << "\n";                                                     \
        if (!(ok == true))                                                               \
            std::cout << "  (expected value was " << (expected) << ")\n";      \
     while (false)

    mpolygon_t ab_intersect;
    auto dump = [](auto label, auto& geo) 
        std::cout << label << " area " << bg::area(geo) << " " << bg::wkt(geo) << "\n";
    ;

    CHECK(true, bg::intersection(a, b, ab_intersect));
    CHECK(1, ab_intersect.size());
    dump("a", a);
    dump("b", b);
    dump("ab_intersect", ab_intersect);

    mpolygon_t a_minus_b, b_minus_a;
    bg::difference(a, ab_intersect, a_minus_b);
    bg::difference(b, ab_intersect, b_minus_a);
    CHECK(1, a_minus_b.size());
    CHECK(1, b_minus_a.size());
    dump("a_minus_b", a_minus_b);
    dump("b_minus_a", b_minus_a);

    auto checks = [dump](auto& reduced, auto& other) 
        CHECK(true, bg::intersects(reduced, other));

        mpolygon_t check;
        bg::intersection(reduced, other, check);

        CHECK(true, check.empty());
        dump("check", check);
    ;

    std::cout << "-- reduced a vs b: -------\n";
    checks(a_minus_b, b);
    std::cout << "-- reduced b vs a: -------\n";
    checks(b_minus_a, a);

打印

PASS    bg::intersection(a, b, ab_intersect) -> true
PASS    ab_intersect.size() -> 1
a area 1694 POLYGON((12227 4967,12238 4967,12238 4813,12227 4813,12227 4967))
b area 5687 POLYGON((12254 4947,12219 4982,12219 5020,12254 5055,12261 5055,12263 5055,12283 5055,12283 4947,12263 4947,12261 4947,12254 4947))
ab_intersect area 8 MULTIPOLYGON(((12234 4967,12238 4967,12238 4963,12234 4967)))
PASS    a_minus_b.size() -> 1
PASS    b_minus_a.size() -> 1
a_minus_b area 1686 MULTIPOLYGON(((12234 4967,12238 4963,12238 4813,12227 4813,12227 4967,12234 4967)))
b_minus_a area 5679 MULTIPOLYGON(((12238 4963,12238 4967,12234 4967,12219 4982,12219 5020,12254 5055,12261 5055,12263 5055,12283 5055,12283 4947,12263 4947,12261 4947,12254 4947,12238 4963)))
-- reduced a vs b: -------
PASS    bg::intersects(reduced, other) -> true
PASS    check.empty() -> true
check area 0 MULTIPOLYGON()
-- reduced b vs a: -------
PASS    bg::intersects(reduced, other) -> true
PASS    check.empty() -> true
check area 0 MULTIPOLYGON()

【讨论】：

感谢您抽出宝贵时间。所以 bg::intersection 不会返回一个有用的值但是 bg::intersects 如果有一个交叉点并且应该返回 true ！这意味着： bg::intersects 在我的情况下返回错误的值？

呃。是的。好像。我想我在那里搞砸了预期值。嗯。也许我以后有时间再去看看

@fhw72 花费了更多时间，使用了更高精度的坐标和各种正确性检查：wandbox.org/permlink/K6ag9oHOH3eJdKq0 我还将输入数据剪辑到了relevant subset - 仍然表现出同样的问题。请注意“区域 2.3e-119”和“区域 1.78e-126”肯定表明存在浮点精度问题。我认为结果以这种方式不一致（准确性）可能不是问题，但它认为库开发人员想知道intersects 可以返回true，而intersections 产生一个空的多面体。

再次感谢您的浏览。您的观察鼓励我提交错误报告。

【参考方案2】：

我提到的行为不是错误，而是完全正确的：

bg::intersection() 返回一个空多边形但是

如果我用线串作为输出参数来调用它：

using linestring_t = bg::model::linestring<point_t>;
linestring_t int1;
bg::intersection(poly2, poly4, int1);
std::cout << "int1: " << bg::wkt(int1) << std::endl;

它有效。

TL;DR bg::intersects() 和 bg::intersection() 的用法和我的解释是错误的！