Geometry Boost List: [geometry] Bug in boost geometry 1.48: union of 2 polygons is wrong.

Subject: [geometry] Bug in boost geometry 1.48: union of 2 polygons is wrong.
From: Zachary Deretsky (zach.deretsky_at_[hidden])
Date: 2012-01-13 20:03:40

Union of 2 polygons in the input3.png produces incorrect polygon in
output3.png.

Below is the code, which reproduces it, You can actually delete everything
below the line "  // THIS IS THE SIMPLE CASE END"

Thanks, Zach.

#include <boost/geometry.hpp>
#include <boost/geometry/geometries/point_xy.hpp>
#include <boost/geometry/geometries/polygon.hpp>
#include <boost/geometry/geometries/ring.hpp>
#include <boost/geometry/multi/geometries/multi_polygon.hpp>
#include <boost/foreach.hpp>

int main()
{
  typedef boost::geometry::model::d2::point_xy<int> point_xy;
  typedef boost::geometry::model::polygon<point_xy > polygon;
  typedef boost::geometry::model::ring<point_xy > ring;
  typedef boost::geometry::model::multi_polygon<polygon > multi_polygon;

  // THIS IS THE SIMPLE CASE BEGIN
  polygon pol1, pol2;
  multi_polygon result;
  boost::geometry::read_wkt("POLYGON((1872000 528000,1872000 192000,1536119
192000,1536000 528000,1200000 528000,1200000 863880,1536000 863880,1872000
863880,1872000 528000))", pol1);
  boost::geometry::read_wkt("POLYGON((1872000 402000,1872000 528000,2768040
528000,2768040 402000,1872000 402000))", pol2);
  boost::geometry::union_(pol1, pol2, result);

  std::cout << "simple example begin" << std::endl;
  std::cout << boost::geometry::wkt(pol1) << std::endl;
  std::cout << boost::geometry::wkt(pol2) << std::endl;
  std::cout << boost::geometry::wkt(result[0]) << std::endl;
  std::cout << boost::geometry::wkt(result) << std::endl;
  std::cout << "simple example end" << std::endl << std::endl;
  // THIS IS THE SIMPLE CASE END

  // I REDUCED IT FROM APPENDING RECTANGLES:
  std::vector<polygon> polygon_vec;
  polygon temp_pol;
  int pol_count = 0;
  boost::geometry::read_wkt("POLYGON((1200000 528000,1200000 863880,1536000
863880,1536000 528000,1200000 528000))", temp_pol);
  polygon_vec.push_back(temp_pol); pol_count++;

  boost::geometry::read_wkt("POLYGON((1536000 528000,1536000 863880,1872000
863880,1872000 528000,1536000 528000))", temp_pol);
  polygon_vec.push_back(temp_pol); pol_count++;

  boost::geometry::read_wkt("POLYGON((1536119 192000,1536119 528000,1872000
528000,1872000 192000,1536119 192000))", temp_pol);
  polygon_vec.push_back(temp_pol); pol_count++;

  boost::geometry::read_wkt("POLYGON((1872000 402000,1872000 528000,2768040
528000,2768040 402000,1872000 402000))", temp_pol);
  polygon_vec.push_back(temp_pol); pol_count++;

  boost::geometry::read_wkt("POLYGON((2208000 528000,2208000 863880,2544000
863880,2544000 528000,2208000 528000))", temp_pol);
  polygon_vec.push_back(temp_pol); pol_count++;

  boost::geometry::read_wkt("POLYGON((1200119 192000,1200119 528000,1536000
528000,1536000 192000,1200119 192000))", temp_pol);
  polygon_vec.push_back(temp_pol); pol_count++;

  polygon current_pol = polygon_vec[0];
  for(int i = 1; i < pol_count; i++) {
    multi_polygon multi;
    boost::geometry::union_(current_pol, polygon_vec[i], multi);
    std::cout << i << std::endl;
    std::cout << boost::geometry::wkt(current_pol) << std::endl;
    std::cout << boost::geometry::wkt(polygon_vec[i]) << std::endl;
    std::cout << boost::geometry::wkt(multi[0]) << std::endl;
    std::cout << boost::geometry::wkt(multi) << std::endl;
    std::cout << std::endl;
    if(1 < multi.size()) {
      return 1;
    } else {
      current_pol = multi[0];
    }
  }
  return 0;
}