Reputation: 101
I'm trying to calculate a difference of two polygons using boost::geometry::difference
with boost::geometry::model::polygon
representing my polygons.
In case when the first polygon contains the second one result of the operation is a single boost::geometry::model::polygon
with the inner and outer rings populated with coordinates of the source polygons.
How do I get a polygon (in the elementary geometry sense) from boost::geometry::model::polygon
?
Clarification:
In elementary geometry, a polygon is a plane figure that is bounded by a finite chain of straight line segments closing in a loop to form a closed chain or circuit.
The outer ring of boost::geometry::model::polygon
is a polygon, the inner rings are polygons too. As a whole boost::geometry::model::polygon
is not a polygon.
So, what I'm asking: How to convert boost::geometry::model::polygon
to a normal polygon (having a single chain of straight line segments), which represents the same area on a plane.
Here's what I'm trying to achieve:
polygon1 = (0,0), (0,8), (8,8), (8,0), (0,0)
polygon2 = (2,2), (2,6), (6,6), (6,2), (2,2)
Polygons 1 & 2 in green / oker:
difference = (0,0), (0,4), (2,4), (2,2), (6,2), (6,6), (2,6), (2,4), (0,4), (0,8), (8,8), (8,0), (0,0)
Expected difference in grey:
I know the same boost::geometry::model::polygon
having inner rings could be represented by infinitely many different normal polygons. I don't care which one I get.
Upvotes: 3
Views: 2495
Reputation: 392903
This is the old answer. After edits to the question I've posted a simple implementation of an algorithm that suits the sample given
It already is.
If you mean, a "simple" non-holey polygon, the outer ring is all you want. Just discard the inner rings.
However, in the most generic case you can end up with multiple entirely disjunct polygons, which is why the output is a collection of polygons. You will have to consider this possibility too (optionally joining the different result polygons into one and discarding inner rings if that's what you desire, functionally).
A sample:
bg::read_wkt("POLYGON((0 0,0 10,10 10,10 0,0 0))", a);
bg::read_wkt("POLYGON((2 -2,2 12,5 12,5 -2,2 -2))", b);
Here, b
cuts a
into two separate pieces. So you must be prepared to handle multiple, disjunct, output polygons:
#include <boost/geometry.hpp>
#include <boost/geometry/geometries/point_xy.hpp>
#include <boost/geometry/geometries/polygon.hpp>
#include <boost/geometry/io/io.hpp>
#include <iostream>
#include <fstream>
namespace bg = boost::geometry;
using pt = bg::model::d2::point_xy<int>;
using poly = bg::model::polygon<pt>;
int main() {
poly a, b;
bg::read_wkt("POLYGON((0 0,0 10,10 10,10 0,0 0))", a);
bg::read_wkt("POLYGON((2 -2,2 12,5 12,5 -2,2 -2))", b);
std::cout << bg::dsv(a) << "\n";
std::cout << bg::dsv(b) << "\n";
{
std::ofstream svg("/tmp/input.svg");
boost::geometry::svg_mapper<pt> mapper(svg, 400, 400);
mapper.add(a);
mapper.add(b);
mapper.map(a, "fill-opacity:0.5;fill:rgb(153,204,0);stroke:rgb(153,204,0);stroke-width:2");
mapper.map(b, "fill-opacity:0.5;fill:rgb(204,153,0);stroke:rgb(202,153,0);stroke-width:2");
}
std::vector<poly> output;
bg::difference(a, b, output);
for (auto& p : output)
std::cout << "\n\t" << bg::dsv(p);
{
std::ofstream svg("/tmp/output.svg");
boost::geometry::svg_mapper<pt> mapper(svg, 400, 400);
assert(output.size() == 2);
mapper.add(output[0]);
mapper.add(output[1]);
mapper.add(b);
mapper.map(output[0], "fill-opacity:0.5;fill:rgb(153,204,0);stroke:rgb(153,204,0);stroke-width:2");
mapper.map(output[1], "fill-opacity:0.5;fill:rgb(153,204,0);stroke:rgb(153,204,0);stroke-width:2");
mapper.map(b, "fill-opacity:0.15;fill:rgb(153,153,153);stroke-width:0");
}
}
Prints:
(((0, 0), (0, 10), (10, 10), (10, 0), (0, 0)))
(((2, -2), (2, 12), (5, 12), (5, -2), (2, -2)))
(((5, 10), (10, 10), (10, 0), (5, 0), (5, 10)))
(((2, 10), (2, 0), (0, 0), (0, 10), (2, 10)))
SVGs rendered:
Upvotes: 0
Reputation: 392903
You can easily construct a ring that models your expected weak simple polygon. First:
Note though that the result isn't valid for further use with the Boost Geometry library's algorithms.
Take your literal example:
std::string reason;
poly expected;
bg::read_wkt("POLYGON((0 0, 0 4, 2 4, 2 2, 6 2, 6 6, 2 6, 2 4, 0 4, 0 8, 8 8, 8 0, 0 0))", expected);
bool ok = bg::is_valid(expected, reason);
std::cout << "Expected: " << bg::dsv(expected) << (ok?" valid":" invalid: '" + reason + "'") << "\n";
Prints
Expected: (((0, 0), (0, 4), (2, 4), (2, 2), (6, 2), (6, 6), (2, 6), (2, 4), (0, 4), (0, 8), (8, 8), (8, 0), (0, 0))) invalid: 'Geometry has invalid self-intersections. A self-intersection point was found at (0, 4); method: t; operations: x/u; segment IDs {source, multi, ring, segment}: {0, -1, -1, 0}/{0, -1, -1, 7}'
With that out of the way, here's a simple algorithm to construct the simple weak polygon from a given polygon:
ring weak_simple_ring(poly& p) {
ring r = p.outer();
for (auto& i: p.inners()) {
auto last = r.back();
r.insert(r.end(), i.rbegin(), i.rend());
r.insert(r.end(), last);
}
return r;
}
The only subtler point there is to reverse the direction (CW/CCW) of the inner rings to match that of the outer ring.
The algorithm doesn't attempt to be smart about finding a cut-point to the inner ring, which probably also means that it won't work nicely for the generic case with multiple inner rings.
Here's a full live demo
Where the input is
bg::read_wkt("POLYGON((0 0,0 10,10 10,10 0,0 0))", a);
bg::read_wkt("POLYGON((2 2, 2 6, 6 6, 6 2, 2 2))", b);
The transformation is
std::vector<poly> output;
bg::difference(a, b, output);
for (auto& p : output) {
ring r = weak_simple_ring(p);
bg::convert(r, p);
}
And the result becomes
Consider when b
had a hole:
bg::read_wkt("POLYGON((0 0,0 10,10 10,10 0,0 0))", a);
bg::read_wkt("POLYGON((2 2, 2 6, 6 6, 6 2, 2 2)(3 3, 5 3, 5 5, 3 5, 3 3))", b);
The output with the same code becomes
Upvotes: 1