Boost polygon discards input polygons

Question

I'm trying to build a set of polygons incrementally using Boost::Polygon. The live dataset in the code below is taken from a pathological input I found in the real-world system I'm using Boost::Polygon for.

Here is the minimum reproduction code (you'll need boost development headers in your include path):

#include <iostream>
#include <string>
#include <vector>

#include "boost/polygon/polygon.hpp"

using namespace std;
using namespace boost::polygon;
using namespace boost::polygon::operators;

/* Utility function to read a polygon from a string. This is ugly, but works for the purpose of this repro case */
polygon_data<float> readPoly(string str)
{
    vector<point_data<float> > result;

    size_t index = 0;
    while(index != string::npos)
    {
        size_t xend = str.find(" ", index);
        size_t yend = str.find(" ", xend + 1);

        float x = atof(str.substr(index, xend - index).c_str());
        float y = atof(str.substr(xend + 1, yend - xend + 1).c_str());
        result.push_back(point_data<float>(x, y));

        if(yend == string::npos)
            break;

        index = yend + 1;
    }

    return polygon_data<float>(result.begin(), result.end());
}

/* Utility function to dump a polygon set to cout, useful for visualizing with the python script */
void dumpPoly(vector<polygon_data<float> > polyset)
{
    for(vector<polygon_data<float> >::iterator it = polyset.begin(); it != polyset.end(); it++)
    {
        polygon_data<float> poly = *it;
        for(polygon_data<float>::iterator_type jt = poly.begin(); jt != poly.end(); jt++)
        {
            cout << (*jt).x() << " " << (*jt).y() << " ";
        }
        cout << endl;
    }
}


int main()
{
    std::vector<polygon_data<float> > data;

    // Construct the polygon set
    data += readPoly("-1309.77, 1323.99, -1324, 1309.76, -1324, -1309.76, -1309.77, -1323.99, -1240.23, -1323.99, -1226, -1309.76, -1226, 1309.76, -1240.23, 1323.99");
    data += readPoly("-1323.99, -1309.77, -1309.76, -1324, 1309.76, -1324, 1323.99, -1309.77, 1323.99, -1240.23, 1309.76, -1226, -1309.76, -1226, -1323.99, -1240.23");
    data += readPoly("1323.99, 1309.77, 1309.76, 1324, -1309.76, 1324, -1323.99, 1309.77, -1323.99, 1240.23, -1309.76, 1226, 1309.76, 1226, 1323.99, 1240.23");
    data += readPoly("-544.771, 686.49, -559, 672.261, -559, -544.761, -544.771, -558.99, -475.229, -558.99, -461, -544.761, -461, 672.261, -475.229, 686.49");
    data += readPoly("-558.99, -544.771, -544.761, -559, 672.261, -559, 686.49, -544.771, 686.49, -475.229, 672.261, -461, -544.761, -461, -558.99, -475.229");
    data += readPoly("686.49, 672.271, 672.261, 686.5, -544.761, 686.5, -558.99, 672.271, -558.99, 602.729, -544.761, 588.5, 672.261, 588.5, 686.49, 602.729");
    data += readPoly("-69.8842, -119.057, -49.7607, -119.057, 219.057, 149.76, 219.057, 169.884, 169.884, 219.057, 149.76, 219.057, -119.057, -49.7607, -119.057, -69.8842");
    data += readPoly("672.271, -558.99, 686.5, -544.761, 686.5, 672.261, 672.271, 686.49, 602.729, 686.49, 588.5, 672.261, 588.5, -544.761, 602.729, -558.99");
    data += readPoly("1309.77, -1323.99, 1324, -1309.76, 1324, 1309.76, 1309.77, 1323.99, 1240.23, 1323.99, 1226, 1309.76, 1226, -1309.76, 1240.23, -1323.99");

    /*

    This alternative dataset shows that boost.polygon can handle nested holes, to some extent

    data += readPoly("100 100 100 1900 1900 1900 1900 100");
    data -= readPoly("200 200 200 1800 1800 1800 1800 200");
    data += readPoly("300 300 300 1700 1700 1700 1700 300");
    data -= readPoly("400 400 400 1600 1600 1600 1600 400");
    data += readPoly("500 500 500 1500 1500 1500 1500 500");
    data -= readPoly("600 600 600 1400 1400 1400 1400 600");

    */

    dumpPoly(data);

    return 0;
}

Here's a dinky little visualizer script I wrote in python using matplotlib. Simply pipe the output from the above program to this script and you'll see the results. It works fine, I just included it to help you help me :)

#!/usr/bin/env python 
import matplotlib.pyplot as plt
import sys

lines = sys.stdin.read().split("\n")

for line in lines:
    data = line.split()

    if len(data) == 0:
        continue

    x = []
    y = []
    i = 0
    while i < len(data):
        x.append(float(data[i]))
        y.append(float(data[i+1]))
        i += 2

    plt.fill(x, y)

plt.show()

Here's what I expect:

Here's what I get:

My question is: why does Boost::Polygon discard some of my input, and how can I prevent it from doing so?

Answer 1

From http://www.boost.org/doc/libs/1_54_0/libs/polygon/doc/index.htm:

The coordinate data type is a template parameter of all data types and algorithms provided by the library, and is expected to be integral. Floating point coordinate data types are not supported by the algorithms implemented in the library due to the fact that the achieving floating point robustness implies a different set of algorithms and generally platform specific assumptions about floating point representations.