boost::polygon boolean subtract results in additional lines

Question

I have a rectangle with a circle that is located inside the rectangle. Now I do a boolean subtract operation where the circle is subtracted from the larger rectangle using boost::polygon and the following (pseudo) code:

boost::polygon::polygon_with_holes_data<int> rect;
boost::polygon::polygon_with_holes_data<int> circle;
// fill both with data here
boost::polygon::polygon_set_data<int> result=rect-circle;

After that "result" contains a bunch of coordinate-pairs where each of them are marked with an additional value -1 or 1. Then I fetch the resulting data the following way:

std::vector<boost::polygon::polygon_data<int>> outputData;
result.get<std::vector<boost::polygon::polygon_data<int>>>(outputData);

Now "outputData" contains a list of coordinates which form the following result:

So the hole is incorporated within one complete, closed polygon coordinate array where I would have expected something like boost::polygon::polygon_with_holes_data where the outline and the holes would have been separated.

So my question: am I doing something wrong here? How do I get the outline and the hole separated so that I can get rid of this extra line which connects them? Is there a function to fetch boost::polygon::polygon_with_holes_data out of the returned boost::polygon::polygon_set_data or what else is the magic here?

Thanks :-)

Answer 1

I took the trouble to generate some data myself. I don't know what data you use, nor how you render it. My suspicion is something is wrong in those processes:

Let's create a polygon_data factories:

template <typename T>
auto make_rect(T x1, T y1, T x2, T y2) {
    auto r = bp::rectangle_data(x1, y1, x2, y2);
    bp::polygon_data<T> p;
    bp::assign(p, r);
    return p;
}

So far so good. For the circle I used your new discovery from your other question:

template <typename T>
auto make_circle(T x, T y, T radius) {
    std::vector<bp::polygon_data<T> > ps;

    bp::assign(ps, bp::rectangle_data(x-1, y-1, x+1, y+1));
    bp::resize(ps, radius, true, 512);

    assert(ps.size() == 1);
    return ps.front();
}

Now let's start out with the straight-forward:

int main() {
    auto rect       = make_rect(0, 0, 1000, 1000);
    auto circle     = make_circle(500, 500, 400);
    auto difference = rect - circle;

    std::cout << std::fixed;
    std::cout << "rect area: " << area(rect)   << "\n";
    std::cout << "circle area: " << area(circle) << "\n";
    std::cout << "difference area: " << area(difference) << "\n";

Prints

rect area: 1000000.000000
circle area: 505481.000000
difference area: 494519.000000

Convert And Check

Now to see whether any of these geometries would render correctly, I know no other technique than via Boost Geometry¹

So, let's convert to a polygon_with_holes_data and do some checks:

bp::polygon_with_holes_data<int> result(rect.begin(), rect.end(),
        &circle, &circle+1);

std::cout << "result area: " << area(result) << "\n";

assert(area(rect) - area(circle) == area(difference));
assert(bp::equivalence(result, rect - circle));

The asserts pass, and this prints the expected:

result area: 494519.000000

Which correctly matches the difference area.

Rendering SVG

Now we adapt in Boost Geometry and write a SVG:

namespace bg = boost::geometry;
std::cout << bg::wkt(rect);
std::cout << bg::wkt(result);

{
    using C = bg::coordinate_type<decltype(result)>::type;
    using PT = bg::model::d2::point_xy<C>;
    std::ofstream svg("output.svg");
    boost::geometry::svg_mapper<PT> mapper(svg, 400, 400);
    mapper.add(result);
    mapper.map(result, "fill-opacity:0.3;fill:rgb(153,204,0);stroke:rgb(153,204,0);stroke-width:2");
}

This writes the WKT and output.svg:

Looks fine here.

What Could Be Wrong? How To Debug?

Like I mentioned in my comments, many library/tools have restrictions/invariants on the input they take. Things like the orientation of outer vs inner rings, self intersections and the presence of closing points on polygon rings.

Now, if you don't know what data is acceptable, you can hit the docs, or you can use the library API to diagnose. For example in Boost Geometry I like to make a quick-and-dirty diagnostic function:

template<typename G>
void auto_correct(G& geo) {
    std::string reason;
    while (!bg::is_valid(geo, reason)) {
        std::cout << "Trying to auto-correct broken input data: " << reason << "\n";
        bg::correct(geo);
    }
}

You can add wkt printing to debug what changes are actually applied.

See e.g. boost read_wkt produced invalid polygon for examples of diagnostics and fixes.

(Note: don't use this in production, the corrections may not do what you expect and the loop may even be infinite for some inputs)

Live Demo

As always: Live On Coliru

#include <boost/polygon/polygon.hpp>
#include <boost/polygon/gtl.hpp>
#include <boost/polygon/rectangle_data.hpp>
#include <boost/polygon/polygon_set_data.hpp>
#include <boost/polygon/polygon_data.hpp>
#include <iostream>

namespace bp = boost::polygon;

#ifndef NO_GEO
#include <fstream>
#include <boost/geometry.hpp>
#include <boost/geometry/io/io.hpp>
#include <boost/geometry/geometries/adapted/boost_polygon.hpp>
namespace bg = boost::geometry;

/* 
 * template<typename G>
 * void auto_correct(G& geo) {
 *     std::string reason;
 *     while (!bg::is_valid(geo, reason)) {
 *         std::cout << "Trying to auto-correct broken input data: " << reason << "\n";
 *         bg::correct(geo);
 *     }
 * }
 */
#endif

template <typename T>
auto make_rect(T x1, T y1, T x2, T y2) {
    auto r = bp::rectangle_data(x1, y1, x2, y2);
    bp::polygon_data<T> p;
    bp::assign(p, r);
    return p;
}

template <typename T>
auto make_circle(T x, T y, T radius) {
    std::vector<bp::polygon_data<T> > ps;

    bp::assign(ps, bp::rectangle_data(x-1, y-1, x+1, y+1));
    bp::resize(ps, radius, true, 512);

    assert(ps.size() == 1);
    return ps.front();
}

int main() {
    auto rect       = make_rect(0, 0, 1000, 1000);
    auto circle     = make_circle(500, 500, 400);
    auto difference = rect - circle;

    std::cout << std::fixed;
    std::cout << "rect   area: " << area(rect)   << "\n";
    std::cout << "circle area: " << area(circle) << "\n";
    std::cout << "difference area: " << area(difference) << "\n";

    bp::polygon_with_holes_data<int> result(rect.begin(), rect.end(),
            &circle, &circle+1);

    std::cout << "result area: " << area(result) << "\n";

    assert(area(rect) - area(circle) == area(difference));
    assert(bp::equivalence(result, rect - circle));

#ifndef NO_GEO
    std::cout << bg::wkt(rect) << "\n";
    std::cout << bg::wkt(result) << "\n";

    {
        using C = bg::coordinate_type<decltype(result)>::type;
        using PT = bg::model::d2::point_xy<C>;
        std::ofstream svg("output.svg");
        boost::geometry::svg_mapper<PT> mapper(svg, 400, 400);
        mapper.add(result);
        mapper.map(result, "fill-opacity:0.3;fill:rgb(153,204,0);stroke:rgb(153,204,0);stroke-width:2");
    }
#endif
}

Prints

rect   area: 1000000.000000
circle area: 505481.000000
difference area: 494519.000000
result area: 494519.000000
POLYGON((0 0,1000 0,1000 1000,0 1000))
POLYGON((0 0,1000 0,1000 1000,0 1000,0 0),(900 527,899 527,899 530,899 532,899 535,899 537,898 537,898 540,898 542,898 545,898 547,897 547,897 550,897 552,897 555,897 557,896 557,896 560,896 562,895 562,895 564,895 566,894 566,894 569,894 571,893 571,893 574,893 576,892 576,892 579,892 581,891 581,891 584,891 586,890 586,890 589,890 591,889 591,889 593,889 595,888 595,888 598,888 600,887 600,887 603,887 605,886 605,885 608,885 610,884 610,884 612,884 614,883 614,883 617,883 619,882 619,881 622,881 624,880 624,880 626,880 628,879 628,878 631,878 633,877 633,876 636,876 638,875 638,875 640,875 642,874 642,873 645,873 647,872 647,871 649,871 651,870 651,869 654,869 656,868 656,867 658,867 660,866 660,865 663,865 665,864 665,863 667,863 669,862 669,861 672,861 674,860 674,859 676,859 678,858 678,857 681,857 683,854 685,854 687,853 687,852 689,852 691,851 691,850 694,850 696,847 698,847 700,846 700,845 702,845 704,842 706,842 708,841 708,840 710,840 712,837 715,837 717,834 719,834 721,833 721,832 723,832 725,829 727,829 729,826 731,826 733,823 735,823 737,820 739,820 741,817 743,817 745,814 747,814 749,811 752,808 754,808 756,805 758,805 760,802 762,802 764,799 767,795 771,792 773,792 775,789 778,785 782,782 785,778 789,775 792,773 792,771 795,767 799,764 802,762 802,760 805,758 805,756 808,754 808,752 811,749 814,747 814,745 817,743 817,741 820,739 820,737 823,735 823,733 826,731 826,729 829,727 829,725 832,723 832,721 833,721 834,719 834,717 837,715 837,712 840,710 840,708 841,708 842,706 842,704 845,702 845,700 846,700 847,698 847,696 850,694 850,691 851,691 852,689 852,687 853,687 854,685 854,683 857,681 857,678 858,678 859,676 859,674 860,674 861,672 861,669 862,669 863,667 863,665 864,665 865,663 865,660 866,660 867,658 867,656 868,656 869,654 869,651 870,651 871,649 871,647 872,647 873,645 873,642 874,642 875,640 875,638 875,638 876,636 876,633 877,633 878,631 878,628 879,628 880,626 880,624 880,624 881,622 881,619 882,619 883,617 883,614 883,614 884,612 884,610 884,610 885,608 885,605 886,605 887,603 887,600 887,600 888,598 888,595 888,595 889,593 889,591 889,591 890,589 890,586 890,586 891,584 891,581 891,581 892,579 892,576 892,576 893,574 893,571 893,571 894,569 894,566 894,566 895,564 895,562 895,562 896,560 896,557 896,557 897,555 897,552 897,550 897,547 897,547 898,545 898,542 898,540 898,537 898,537 899,535 899,532 899,530 899,527 899,527 900,525 900,523 900,521 900,518 900,516 900,513 900,511 900,508 900,506 900,503 900,501 900,498 900,496 900,493 900,491 900,488 900,486 900,483 900,481 900,478 900,476 900,474 900,472 900,469 899,467 899,464 899,462 899,459 898,457 898,454 898,452 898,449 897,447 897,444 897,442 897,439 896,437 896,435 895,433 895,430 894,428 894,425 893,423 893,420 892,418 892,415 891,413 891,410 890,408 890,406 889,404 889,401 888,399 888,396 887,394 887,391 885,389 885,387 884,385 884,382 883,380 883,377 881,375 881,373 880,371 880,368 878,366 878,363 876,361 876,359 875,357 875,354 873,352 873,348 871,345 869,343 869,339 867,336 865,334 865,330 863,327 861,325 861,321 859,318 857,316 857,312 854,308 852,305 850,303 850,299 847,295 845,291 842,287 840,282 837,278 834,274 832,270 829,266 826,262 823,258 820,254 817,250 814,247 811,243 808,239 805,235 802,232 799,228 795,224 792,221 789,217 785,214 782,210 778,207 775,204 771,200 767,197 764,194 760,191 756,188 752,185 749,182 745,179 741,176 737,173 733,170 729,167 725,165 721,162 717,159 712,157 708,154 704,152 700,149 696,147 691,145 687,142 683,140 678,138 674,136 669,134 665,132 660,130 656,128 651,126 647,124 642,123 638,121 633,119 628,118 624,116 619,115 614,114 610,112 605,111 600,110 595,109 591,108 586,107 581,106 576,105 571,104 566,103 562,102 557,102 555,102 552,101 547,101 545,101 542,100 537,100 535,100 532,99 527,99 525,99 523,99 521,99 518,99 516,99 513,99 511,99 508,99 506,99 503,99 501,99 498,99 496,99 493,99 491,99 488,99 486,99 483,99 481,99 478,99 476,99 474,99 472,100 467,100 464,100 462,101 457,101 454,101 452,102 447,102 444,102 442,103 437,104 433,105 428,106 423,107 418,108 413,109 408,110 404,111 399,112 394,114 389,115 385,116 380,118 375,119 371,121 366,123 361,124 357,126 352,128 348,130 343,132 339,134 334,136 330,138 325,140 321,142 316,145 312,147 308,149 303,152 299,154 295,157 291,159 287,162 282,165 278,167 274,170 270,173 266,176 262,179 258,182 254,185 250,188 247,191 243,194 239,197 235,200 232,204 228,207 224,210 221,214 217,217 214,221 210,224 207,228 204,232 200,235 197,239 194,243 191,247 188,250 185,254 182,258 179,262 176,266 173,270 170,274 167,278 165,282 162,287 159,291 157,295 154,299 152,303 149,308 147,312 145,316 142,321 140,325 138,330 136,334 134,339 132,343 130,348 128,352 126,357 124,361 123,366 121,371 119,375 118,380 116,385 115,389 114,394 112,399 111,404 110,408 109,413 108,418 107,423 106,428 105,433 104,437 103,442 102,444 102,447 102,452 101,454 101,457 101,462 100,464 100,467 100,472 99,474 99,476 99,478 99,481 99,483 99,486 99,488 99,491 99,493 99,496 99,498 99,501 99,503 99,506 99,508 99,511 99,513 99,516 99,518 99,521 99,523 99,525 99,527 99,532 100,535 100,537 100,542 101,545 101,547 101,552 102,555 102,557 102,562 103,566 104,571 105,576 106,581 107,586 108,591 109,595 110,600 111,605 112,610 114,614 115,619 116,624 118,628 119,633 121,638 123,642 124,647 126,651 128,656 130,660 132,665 134,669 136,674 138,678 140,683 142,687 145,691 147,696 149,700 152,704 154,708 157,712 159,717 162,721 165,725 167,729 170,733 173,737 176,741 179,745 182,749 185,752 188,756 191,760 194,764 197,767 200,771 204,775 207,778 210,782 214,785 217,789 221,792 224,795 228,799 232,802 235,805 239,808 243,811 247,814 250,817 254,820 258,823 262,826 266,829 270,832 274,834 278,837 282,840 287,842 291,845 295,847 299,850 303,850 305,852 308,854 312,857 316,857 318,859 321,861 325,861 327,863 330,865 334,865 336,867 339,869 343,869 345,871 348,873 352,873 354,875 357,875 359,876 361,876 363,878 366,878 368,880 371,880 373,881 375,881 377,883 380,883 382,884 385,884 387,885 389,885 391,887 394,887 396,888 399,888 401,889 404,889 406,890 408,890 410,891 413,891 415,892 418,892 420,893 423,893 425,894 428,894 430,895 433,895 435,896 437,896 439,897 442,897 444,897 447,897 449,898 452,898 454,898 457,898 459,899 462,899 464,899 467,899 469,900 472,900 474,900 476,900 478,900 481,900 483,900 486,900 488,900 491,900 493,900 496,900 498,900 501,900 503,900 506,900 508,900 511,900 513,900 516,900 518,900 521,900 523,900 525,900 527))

---

¹ (I don't know how to output/input any geometry data using Boost Polygon, and I don't have any WKT/DSV capable rendering tools).