Geometry Boost List: [ggl] intersection of two vectors of polygons

Subject: [ggl] intersection of two vectors of polygons
From: Simonson, Lucanus J (lucanus.j.simonson)
Date: 2011-04-28 13:35:48

Vishnu wrote:
> I'm wondering whether anyone has recommendations for calculating the
> intersection between two vectors of polygons.
> 
> Currently, I use boost::geometry::intersection(). Before calling this
> function, I convert the 2 vectors of polygons into 2 valid
> multi_polygons. The problem with this approach is that instead of
> getting 9 discrete polygons in the intersection result, I get a
> single polygon that is the union of the 9 polygons.
> 
> The reason it's important to get 9 discrete polygons in the
> intersection is as follows. Say that 4 adjacent squares in the first
> collection are 4 farms labeled by ownership and the 4 adjacent
> squares in the 2nd set are labeled according to mineral rights. Let's
> say one wants the polygons corresponding to farm with owner "1" and
> mineral rights "3". 
> 
> This is what I have right now. Each vector of polygons in the inputs
> has 4 squares. The second vector of 4 squares is a copy of the first
> set, offset slightly in the +x and +y directions. The intersection
> should have discrete 9 squares.
> 
> The first collection of polygons is:
> 
> POLYGON((0 0.5,0 1,0.5 1,0.5 0.5,0 0.5))POLYGON((0.5 0.5,0.5 1,1 1,1
> 0.5,0.5 
> 0.5))POLYGON((0 0,0 0.5,0.5 0.5,0.5 0,0 0))POLYGON((0.5 0,0.5 0.5,1
> 0.5,1 0,0.5 0))
> 
> and the second collection of polygons is
> 
> POLYGON((0.25 0.75,0.25 1.25,0.75 1.25,0.75 0.75,0.25
> 0.75))POLYGON((0.75 
> 0.75,0.75 1.25,1.25 1.25,1.25 0.75,0.75 0.75))POLYGON((0.25 0.25,0.25
> 0.75,0.75 0.75,0.75 0.25,0.25 0.25))POLYGON((0.75 0.25,0.75 0.75,1.25
> 0.75,1.25 0.25,0.75 0.25))
> 
> To convert each vector of polygons into a valid multi_polygon, I move
> down a vector and take the union of the current polygon with the
> union of all preceding polygons. I can post the code snippet if
> necessary. 
> 
> The problem, as I see it is that the multi_polygons are too
> restrictive in that they don't allow a collection of 4 adjacent
> squares such as a 2x2 cell on graph paper to be treated as discrete
> polygons. 
> 
> The alternative I am thinking of is to intersect each polygon in the
> first vector with each polygon in the second vector and then put the
> resulting polygons into a 3rd vector. This would be an O(n^2)
> operation, hence bad for large vectors of inputs.
> 
> Does anyone see another way? I think that having a way to do
> intersections on vectors of polygons and not just multi_polygons
> would be a valuable feature in boost-geometry.

You want map overlay.  I have implemented it in Boost.Polygon.  I call it property merge.  Your mineral rights and ownership are properties.  The result of property merge is a map of sets of properties to polygon sets that share that set of properties.  You can then query for regions with various combinations of properties.  I support M different property types on M different multi polygons in a single sweep and avoid quadratic runtime.  For your example I would assign a different property to each farm for ownership and a property for mineral rights and then compute the map overlay.

http://www.boost.org/doc/libs/1_46_1/libs/polygon/doc/gtl_property_merge.htm

I agree with you that it is a good feature to have in a GIS oriented library.  Unfortunately, the implementation details of boost.geometry polygon intersection make extending it to N multi polygons with N properties non-trivial.

An alternative you may consider that would allow you to use Boost.Geometry is to shrink each of your polygons by a small ammount so that they are not abutting.  In this way you will get slighly shrunk intersection results that are not merged by the algorithm.  Since you have only two properties on two multi-polygons that should satisfy your immediate concern, but is still N^2 for N properties.

Regards,
Luke