[postgis-users] Using PostGIS for partial intersection/union of many polygons

[Sandro Santilli]

On Tue, Jul 30, 2013 at 07:37:47PM -0700, Martin Davis wrote:
> This kind of question can be answered by a classic N-polygon
> overlay.  AFAIK PostGIS doesn't have such a function directly, so
> you will have to compute the noding of the polygon linework,
> followed by a polygonization step.  Then take a interior point of
> each resultant polygon and compute how many source polygons it lies
> in.  You can then union the set of polygons for each count of
> interest.
> 
> Be warned - if the sample below is representative of your data, you
> may well encounter TopologyExceptions, since it appears that there
> is a lot of nearly coincident linework.
> 
> No doubt others can weigh in if there is a better way of doing this
> in PostGIS.

Great idea, Martin!
It could be implemented by building a topology and checking out the faces
that partecipate in more than N TopoGeometry definitions. Would probably
be pretty slow, but at least exceptions could be dealt with on a more
localized way...

--strk;

> 
> On 7/30/2013 9:37 AM, Sparr wrote:
> >It was suggested on IRC that I ask this question to this mailing list.
> >
> >I have a set of 1000 (or many more) polygons, most of which
> >intersect with each other. I want to produce a polygon (or set of
> >polygons, or an empty set) describing the area covered by
> >N-or-more of those polygons. For N=1, the answer is the union of
> >all of the polygons. For N=1000, the answer is their intersection.
> >The naive way to calculate N=2 through N=999 is to make a list of
> >combinations of polygons, intersect those, and then union the
> >results, but this will require a prohibitively large number of
> >union and intersection operations in some cases.
> >
> >To illustrate the problem, here is an example with a smaller set
> >of data (422 polygons right now):
> >
> >http://regionaldifferences.com/results.html?region=New%20England&lat=42&lon=-73&zoom=6
> >
> >If you hover your mouse over the center of Massachusetts you'll
> >see that 98% of the polygons (415/422 currently) intersect there.
> >Syracuse NY is about 25%.
> >
> >What I want is the outline of specific percentiles on this map. I
> >want a single polygon representing the 50%-or-more area, and a
> >single polygon representing the 90%-or-more area, etc.
> >
> >Can PostGIS do this without my needing to assemble a string of
> >union and intersection operations for the query? I would prefer
> >not to rasterize the data unless it proves infeasible to do this
> >with the full precision vector data.