'Postgresql: Calculate the number of polygon islands or overlapping groups'

This post is about some work I did on how to calculate the number of number of groups of intersecting polygons in an array of polygons in Postgresql

First lets be a little more strict about what the definition of an island is. For the purposes of this algorithm an island consists of a set polygons where at least one polygon in the set overlaps another in the set.

I created a function to calculate the number of islands attached is the plpgsql code for the function.

CREATE OR REPLACE FUNCTION getNumIslands(polys polygon[]) RETURNS INT AS
$$BODY$$
DECLARE
islandList polygon[]; -- holds the set of obstacles that make up the current island.
islandCount INT; -- holds the count of the number of islands.
poly polygon; -- holds the polygon to start an island search
polys_copy polygon[]; -- A copy of the input list
inddex INT; -- keeps track of the current index in the polys_copy search

BEGIN

islandCount := 0;

islandList := polys;
polys_copy := polys;
inddex := 1;

WHILE array_length(polys_copy,1) > 0
LOOP
-- RAISE NOTICE 'Poly_list length is: %', array_length(polys_copy, 1);
-- grab the first poly
poly := polys_copy;
islandList := array_append(ARRAY[]::polygon[], poly);
-- RAISE NOTICE 'Current island is %', islandList;
-- remove the first poly because we dont care about it overlaping itself.
polys_copy := polys_copy[2:array_length(polys_copy, 1)];

-- Make a list of all the polys left that overlap with poly
WHILE inddex <= array_length(polys_copy, 1)
LOOP
-- If the poly_copy overlaps the poly then add its index to the list
-- of polys to be removed
FOREACH poly IN ARRAY islandList
LOOP
IF polys_copy[inddex] && poly
THEN
islandList := array_append(islandList, polys_copy[inddex]);
--
IF array_length(polys_copy,1) = 1
THEN
polys_copy := null;
ELSE
polys_copy := polys_copy[1:inddex-1] || polys_copy[inddex+1:array_length(polys_copy, 1)];
END IF;
-- RAISE NOTICE 'Intersects with poly index: %', inddex;
-- When a new intersecting poly is found it needs to be compared
-- to every poly not just the polys left in the list
inddex := 0;
EXIT;
END IF;
END LOOP;
-- RAISE NOTICE 'Current island is %', islandList;
inddex := inddex + 1;
END LOOP;

islandCount := islandCount + 1;
inddex := 1;
END LOOP;

RETURN islandCount;

END

$$BODY$$ LANGUAGE plpgsql;
© 2013 Glen Berseth | www.fracturedplane.com