Now that the data is loaded into a table in the database, we can leverage PostGIS to flatten and get the union of the polygons, so that we have a normalized dataset. The first step in doing so using this approach will be to convert the polygons to linestrings. We can then link those linestrings and convert them back to polygons, representing the union of all the polygon inputs. We will perform the following tasks:
- Convert polygons to linestrings
- Convert linestrings back to polygons
- Find the center points of the resultant polygons
- Use the resultant points to query tabular relationships
To convert polygons to linestrings, we'll need to extract just the portions of the polygons we want using ST_ExteriorRing, convert those parts to points using ST_DumpPoints, and then connect those points back into lines like a connect-the-dots coloring book using ST_MakeLine.
Breaking it down further, ST_ExteriorRing (the_geom) will grab just the outer boundary of our polygons. But ST_ExteriorRing returns polygons, so we need to take that output and create a line from it. The easiest way to do this is to convert it to points using ST_DumpPoints and then connect those points. By default, the Dump function returns an object called a geometry_dump, which is not just simple geometry, but the geometry in combination with an array of integers. The easiest way to return the geometry alone is the leverage object notation to extract just the geometry portion of geometry_dump, as follows:
(ST_DumpPoints(geom)).geom
Piecing the geometry back together with ST_ExteriorRing is done using the following query:
SELECT (ST_DumpPoints(ST_ExteriorRing(geom))).geom
This should give us a listing of points in order from the exterior rings of all the points from which we want to construct our lines using ST_MakeLine, as shown in the following query:
SELECT ST_MakeLine(geom) FROM ( SELECT (ST_DumpPoints(ST_ExteriorRing(geom))).geom) AS linpoints
Since the preceding approach is a process we may want to use in many other places, it might be prudent to create a function from this using the following query:
CREATE OR REPLACE FUNCTION chp02.polygon_to_line(geometry)
RETURNS geometry AS
$BODY$
SELECT ST_MakeLine(geom) FROM (
SELECT (ST_DumpPoints(ST_ExteriorRing((ST_Dump($1)).geom))).geom
) AS linpoints
$BODY$
LANGUAGE sql VOLATILE;
ALTER FUNCTION chp02.polygon_to_line(geometry)
OWNER TO me;
Now that we have the polygon_to_line function, we still need to force the linking of overlapping lines in our particular case. The ST_Union function will aid in this, as shown in the following query:
SELECT ST_Union(the_geom) AS geom FROM (
SELECT chp02.polygon_to_line(geom) AS geom FROM
chp02.use_area
) AS unioned;
Now let's convert linestrings back to polygons, and for this we can polygonize the result using ST_Polygonize, as shown in the following query:
SELECT ST_Polygonize(geom) AS geom FROM (
SELECT ST_Union(the_geom) AS geom FROM (
SELECT chp02.polygon_to_line(geom) AS geom FROM
chp02.use_area
) AS unioned
) as polygonized;
The ST_Polygonize function will create a single multi polygon, so we need to explode this into multiple single polygon geometries if we are to do anything useful with it. While we are at it, we might as well do the following within a CREATE TABLE statement:
CREATE TABLE chp02.use_area_alt AS (
SELECT (ST_Dump(the_geom)).geom AS the_geom FROM (
SELECT ST_Polygonize(the_geom) AS the_geom FROM (
SELECT ST_Union(the_geom) AS the_geom FROM (
SELECT chp02.polygon_to_line(the_geom) AS the_geom
FROM chp02.use_area
) AS unioned
) as polygonized
) AS exploded
);
We will be performing spatial queries against this geometry, so we should create an index in order to ensure our query performs well, as shown in the following query:
CREATE INDEX chp02_use_area_alt_the_geom_gist ON chp02.use_area_alt USING gist(the_geom);
In order to find the appropriate table information from the original geometry and apply that back to our resultant geometries, we will perform a point-in-polygon query. For that, we first need to calculate centroids on the resultant geometry:
CREATE TABLE chp02.use_area_alt_p AS
SELECT ST_SetSRID(ST_PointOnSurface(the_geom), 3734) AS
the_geom FROM
chp02.use_area_alt;
ALTER TABLE chp02.use_area_alt_p ADD COLUMN gid serial;
ALTER TABLE chp02.use_area_alt_p ADD PRIMARY KEY (gid);
And as always, create a spatial index using the following query:
CREATE INDEX chp02_use_area_alt_p_the_geom_gist ON chp02.use_area_alt_p USING gist(the_geom);
The centroids then structure our point-in-polygon (ST_Intersects) relationship between the original tabular information and resultant polygons, using the following query:
CREATE TABLE chp02.use_area_alt_relation AS
SELECT points.gid, cu.location FROM
chp02.use_area_alt_p AS points,
chp02.use_area AS cu
WHERE ST_Intersects(points.the_geom, cu.the_geom);
If we view the first rows of the table, we can see it links the identifier of points to their respective locations:
