Using Boost.Polygon to slice manhattan polygons

Question

I am having trouble slicing a Manhattan (rectilinear) polygon into rectangles using the get_rectangles(output_container_type& output, const T& polygon_set) method in Boost.Polygon. It does not seem to work on self-intersecting polygons, such as the one below.

Here's my attempt:

#include <iostream>
#include <boost/polygon/polygon.hpp>
#include <vector>

namespace gtl = boost::polygon;
using namespace boost::polygon::operators;

int main() {
    typedef gtl::polygon_90_with_holes_data<int> Polygon;
    typedef gtl::polygon_traits<Polygon>::point_type Point;

    Point pts[] = {gtl::construct<Point>(0, 0), // See the image
    gtl::construct<Point>(0, 10),
    gtl::construct<Point>(30, 10),
    gtl::construct<Point>(30, 20),
    gtl::construct<Point>(10, 20),
    gtl::construct<Point>(10, 0)};

    Polygon poly;
    gtl::set_points(poly, pts, pts+6);

    std::vector< gtl::rectangle_data<int> > rects;
    get_rectangles(rects, poly);

    std::cout << rects.size() << " rectangle: \n";

    for(std::vector<gtl::rectangle_data<int> >::iterator it = rects.begin(); it !=                     
        rects.end(); ++it) {
            // Print out the corner coordinates
            std::cout << "x1: "<< gtl::xl(*it) << ", x2: " << gtl::xh(*it) 
            << ", y1: "<< gtl::yl(*it) << ", y2: " << gtl::yh(*it) << std::endl;
    }
    return 0;
}

And here's the output:

    1 rectangle:
    x1: 10, x2: 30, y1: 10, y2: 20

This method worked for a non-intersecting polygon, but here it only produces one rectangle, not two. This seems to indicate that what I am trying to do should be possible. What am I doing wrong, or is this possible?

Answer 1

Your code is "stress-testing" the polygon code. Strictly speaking, your list of points does not define a polygon. In theory, a planar region (a polygon shape in particular) is well defined only by a so-called Jordan (simple) curve that forms the perimeter. In the simplest terms, a Jordan (simple) curve must not cross itself.

When the same shape is redefined using a Jordan (simple) curve with the following list of 8 points,

Point pts[] = {gtl::construct<Point>(0, 0),
gtl::construct<Point>(0, 10),
gtl::construct<Point>(10, 10),    
gtl::construct<Point>(10, 20),
gtl::construct<Point>(30, 20),
gtl::construct<Point>(30, 10),
gtl::construct<Point>(10, 10),
gtl::construct<Point>(10, 0)};

your code produces the correct result:

2 rectangle: x1: 0, x2: 10, y1: 0, y2: 10 x1: 10, x2: 30, y1: 10, y2: 20