Search code examples
c++boostgeometryunionboost-geometry

Boost union does not work in 1.67 but works in 1.61. Why?

I am trying to work out why my code below works as expected in boost 1.61 but not in boost 1.67.

In boost 1.61 the input polygons are correctly combined and the outline polygon displayed.

In boost 1.67, with no code changes, the outline polygon is wrong and incomplete. It's as if it has missing points.

Output images added to help explain. There is also an #define BOOST_1_6_1 I needed to add as boost 1.61 does not seem to automatically add that header file.

Uncomment drawAllPolygons() if you want to see how the input polygons look.

Can anyone help? 

#include <iostream>


#include <boost\geometry.hpp>

//#define BOOST_1_6_1
#ifdef BOOST_1_6_1
#include <boost/geometry/geometries/point_xy.hpp>
#endif

typedef boost::geometry::model::d2::point_xy<int> intPt;
typedef boost::geometry::model::polygon<intPt> polygon;

const int GRID_WIDTH = 220;
const int GRID_HEIGHT = 60;

bool canvas[GRID_WIDTH][GRID_HEIGHT];

std::vector<polygon> output;
std::vector<polygon> input;

void initCanvas()
{
    for (unsigned int y = 0; y < GRID_HEIGHT; ++y)
    {
        for (unsigned int x = 0; x < GRID_WIDTH; ++x)
        {
            canvas[x][y] = false;
        }
    }
}

void drawGrid()
{
    for (unsigned int y = 0; y < GRID_HEIGHT; ++y)
    {
        for (unsigned int x = 0; x < GRID_WIDTH; ++x)
        {
            if (canvas[x][y])
            {
                std::cout << "x";
            } 
            else
            {
                std::cout << ".";
            }
        }

        std::cout << std::endl;
    }
}

polygon setupPolygon(const int startX, const int startY, const int width, const int height)
{
    polygon polygon1;

    int endX = startX + width;
    int endY = startY + height;

    for (int x = startX; x <= endX; ++x)
    {
        intPt pt(x, startY);
        polygon1.outer().push_back(pt);
    }

    for (int y = startY; y <= endY; ++y)
    {
        intPt pt(endX, y);
        polygon1.outer().push_back(pt);
    }

    for (int x = endX; x >= startX; --x)
    {
        intPt pt(x, endY);
        polygon1.outer().push_back(pt);
    }

    for (int y = endY; y >= startY; --y)
    {
        intPt pt(startX, y);
        polygon1.outer().push_back(pt);
    }

    return polygon1;
}

std::vector<polygon> combine(std::vector<polygon> input)
{
    bool loop = true;
    while (loop)
    {
        unsigned int before = input.size();
        bool exit = false;
        for (unsigned int i = 0; i < input.size() && !exit; ++i)
        {
            for (unsigned int j = i + 1; j < input.size() && !exit; ++j)
            {
                std::vector<polygon> output;
                boost::geometry::correct(input[i]);
                boost::geometry::correct(input[j]);
                boost::geometry::union_(input[i], input[j], output);

                if (i < j)
                { 
                    input.erase(input.begin() + j);
                    input.erase(input.begin() + i);
                }
                else
                {
                    input.erase(input.begin() + i);
                    input.erase(input.begin() + j);
                }

                input.insert(input.begin(), output.begin(), output.end());

                exit = true;
            }
        }

        if (before == input.size())
        {
            loop = false;
        }
    }


    return input;
}

void drawAllPolygons()
{
    for (unsigned int i = 0; i < input.size(); ++i)
    {
        auto outer = input[i].outer();
        for (unsigned int i = 0; i < outer.size(); ++i)
        {
            int x = outer[i].get<0>();
            int y = outer[i].get<1>();
            canvas[x][y] = true;
        }
    }
}

void drawCombinedPolygons()
{
    for (unsigned int j = 0; j < output.size(); ++j)
    {
        auto outer = output[j].outer();
        for (unsigned int i = 0; i < outer.size(); ++i)
        {
            int x = outer[i].get<0>();
            int y = outer[i].get<1>();
            canvas[x][y] = true;
        }
    }
}

int main()
{
    initCanvas();

    input.push_back(setupPolygon(40, 10, 50, 15));
    input.push_back(setupPolygon(40, 20, 50, 15));
    input.push_back(setupPolygon(50, 10, 50, 15));
    input.push_back(setupPolygon(60, 30, 50, 15));

    output = combine(input);

    //drawAllPolygons();
    drawCombinedPolygons();

    drawGrid();

    std::cin.get();

    return 0;
}

enter image description here enter image description here

Solution

  • Boost 1.61 apparently didn't do all the possible optimizations.

    There are quite a few misconceptions.

    1. The orientations for the input polygon (rings) are wrong. You can help your self using is_valid and correct. The code below will show it.

    2. The original polygons are expressly created to contain many redundant points. Even if you just "union" one of those with themselves, you can expect the result to be simplified:

      Live On Coliru

      polygon p = setupPolygon(40, 10, 50, 15);
std::cout << "original: " << bg::wkt(p) << "\n";
multi_polygon simple;
bg::union_(p, p, simple);
std::cout << "self-union: " << bg::wkt(simple) << "\n";

      Would print 

      original: POLYGON((40 10,40 11,40 12,40 13,40 14,40 15,40 16,40 17,40 18,40 19,40 20,40 21,40 22,40 23,40 24,40 25,40 25,41 25,42 25,43 25,44 25,45 25,46 25,47 25,48 25,49 25,50 25,51 25,52 25,53 25,54 25,55 25,56 25,57 25,58 25,59 25,60 25,61 25,62 25,63 25,64 25,65 25,66 25,67 25,68 25,69 25,70 25,71 25,72 25,73 25,74 25,75 25,76 25,77 25,78 25,79 25,80 25,81 25,82 25,83 25,84 25,85 25,86 25,87 25,88 25,89 25,90 25,90 25,90 24,90 23,90 22,90 21,90 20,90 19,90 18,90 17,90 16,90 15,90 14,90 13,90 12,90 11,90 10,90 10,89 10,88 10,87 10,86 10,85 10,84 10,83 10,82 10,81 10,80 10,79 10,78 10,77 10,76 10,75 10,74 10,73 10,72 10,71 10,70 10,69 10,68 10,67 10,66 10,65 10,64 10,63 10,62 10,61 10,60 10,59 10,58 10,57 10,56 10,55 10,54 10,53 10,52 10,51 10,50 10,49 10,48 10,47 10,46 10,45 10,44 10,43 10,42 10,41 10,40 10))
self-union: MULTIPOLYGON(((90 25,90 10,40 10,40 25,90 25)))

      As you can see, the result will be simplified to just the obvious corners of the rectangle. However your "drawing" code relies on the redundant points to be present.

      This is why you end up with "spotty" results - that seem to miss points, only because your drawing code doesn't draw any lines. It merely draws the (corner) points.

    3. The canvas dimensions are inverted

      I think you had your canvas indices consistently swapped around. I'd suggest using more explicit data-structures with (opional) bounds checking to avoid this.

    Change to use SVG

    If you change the code to output the SVG of the input and output, you will find that the whole thing makes perfect sense. Note that this also takes care to correct the input:

    Live On Coliru

    #include <iostream>
#include <fstream>

#include <boost/geometry.hpp>
#include <boost/geometry/geometries/point_xy.hpp>

namespace bg = boost::geometry;
typedef bg::model::d2::point_xy<int> intPt;
typedef bg::model::polygon<intPt> polygon;
typedef bg::model::multi_polygon<polygon> multi_polygon;

const int GRID_WIDTH = 220;
const int GRID_HEIGHT = 60;

std::array<std::array<bool, GRID_WIDTH>, GRID_HEIGHT> canvas;

void initCanvas() {
    for (auto& row : canvas)
        for (auto& cell : row)
            cell = false;
}

void drawGrid() {
    for (auto& row : canvas) {
        for (auto& cell : row)
            std::cout << (cell?"x":".");
        std::cout << std::endl;
    }
}

template <typename R> void valid(R& g) {
    std::string reason;
    while (!bg::is_valid(g, reason)) {
        std::cout << "Not valid: " << reason << "\n";
        bg::correct(g);
    }
}

polygon setupPolygon(const int startX, const int startY, const int width, const int height) {
    polygon::ring_type r;

    int const endX = startX + width;
    int const endY = startY + height;

    for (int x = startX; x <= endX;   ++x) r.emplace_back(x,      startY);
    for (int y = startY; y <= endY;   ++y) r.emplace_back(endX,   y);
    for (int x = endX;   x >= startX; --x) r.emplace_back(x,      endY);
    for (int y = endY;   y >= startY; --y) r.emplace_back(startX, y);

    valid(r);

    return {r};
}

template <typename Polygons>
multi_polygon combine(Polygons const& input) {

    if (input.empty()) return {};
    //if (input.size()<2) return {input.front()};

    multi_polygon out;
    for (polygon const& i : input) {
        multi_polygon tmp;
        bg::union_(out, i, tmp);
        out = tmp;
    }

    return out;
}

void draw(polygon const& p) {
    for (auto& pt : p.outer())
        canvas.at(pt.y()).at(pt.x()) = true;
}

void draw(multi_polygon const& mp) { for (auto& p : mp) draw(p); }
void draw(std::vector<polygon> const& mp) { for (auto& p : mp) draw(p); }

void writeSvg(multi_polygon const& g, std::string fname) {
    std::ofstream svg(fname);
    boost::geometry::svg_mapper<intPt> mapper(svg, 400, 400);
    mapper.add(g);
    mapper.map(g, "fill-opacity:0.5;fill:rgb(0,0,153);stroke:rgb(0,0,200);stroke-width:2");
}

void writeSvg(std::vector<polygon> const& g, std::string fname) {
    std::ofstream svg(fname);
    boost::geometry::svg_mapper<intPt> mapper(svg, 400, 400);
    for (auto& p: g) {
        mapper.add(p);
        mapper.map(p, "fill-opacity:0.5;fill:rgb(153,0,0);stroke:rgb(200,0,0);stroke-width:2");
    }
}

int main() {
    std::vector<polygon> input { {
        setupPolygon(40, 10, 50, 15),
        setupPolygon(40, 20, 50, 15),
        setupPolygon(50, 10, 50, 15),
        setupPolygon(60, 30, 50, 15),
    } };

    for (auto& p : input)
        valid(p);

    auto output = combine(input);

    initCanvas();
    draw(input);
    drawGrid();

    initCanvas();
    draw(output);
    drawGrid();

    writeSvg(input, "input.svg");
    writeSvg(output, "output.svg");
}

    Prints

    Not valid: Geometry has wrong orientation
Not valid: Geometry has wrong orientation
Not valid: Geometry has wrong orientation
Not valid: Geometry has wrong orientation
............................................................................................................................................................................................................................
............................................................................................................................................................................................................................
............................................................................................................................................................................................................................
............................................................................................................................................................................................................................
............................................................................................................................................................................................................................
............................................................................................................................................................................................................................
............................................................................................................................................................................................................................
............................................................................................................................................................................................................................
............................................................................................................................................................................................................................
............................................................................................................................................................................................................................
........................................xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx.......................................................................................................................
........................................x.........x.......................................x.........x.......................................................................................................................
........................................x.........x.......................................x.........x.......................................................................................................................
........................................x.........x.......................................x.........x.......................................................................................................................
[ snip ]
............................................................................................................................................................................................................................
............................................................................................................................................................................................................................
............................................................................................................................................................................................................................

    And the input.svg:

    enter image description here

    And the output.svg:

    enter image description here