[GIS] PostGIS: Geometry column registering as GeometryCollection, not converting to MULTIPOLYGON

geometry-conversionpostgis

I'm trying to aggregate data from a tabular data source and then join it to a geometry file (sl_adm2) in PostGIS using a primary key. The geometry file was imported using shp2psql.

CREATE TABLE example_geom AS 
SELECT sl_adm2_sums_alias.* , adm2_code as district_pcode 
FROM gis_schema.sl_adm2
JOIN (
    SELECT st_collect(sl_adm2.geom) AS geometry, sum(sl_sos.HEALTH) AS health_sum,
        sum(sl_sos.PROTECTION) AS protection_sum, sum(sl_sos.WASH) AS wash_sum,
        sum(sl_sos.cccm) AS cccm_sum, sum(sl_sos.food) AS food_sum, 
        sum(sl_sos.livelihood) AS livelihood_sum, 
        sum(sl_sos.education) AS education_sum, 
        sum(sl_sos.agriculture) AS agriculture_sum, 
        sum(sl_sos.HUMANITARIAN_ACCESS) AS humaccess_sum, 
        sum(sl_sos.logistics) AS logistics_sum, 
        sum(sl_sos.SHELTER_AND_NFI) AS shelter_and_wifi_sum, sl_sos.district
    FROM sl_sos 
    JOIN gis_schema.sl_adm2 ON sl_sos.DISTRICT_PCODE = sl_adm2.adm2_code
    GROUP BY sl_sos.district)
sl_adm2_sums_alias ON
sl_adm2.adm2_name = sl_adm2_sums_alias.district;

and then

ALTER TABLE example_geom
ADD PRIMARY KEY (district_pcode);

The aggregation and join work fine but the geometry column in the resulting table does not appear to work. Further investigation reveals that the geometry column is registering as a GeometryCollection. I tried to change this with the following code.

ALTER TABLE example_geom ALTER COLUMN geometry
SET DATA TYPE geometry(MultiPolygon) USING ST_Multi(geometry);

This gave me the below error

ERROR: Geometry type (GeometryCollection) does not match column type (MultiPolygon)
SQL state: 22023

I've tried to fix this or get around it with no success. Could anyone help?

I'm using Postgres 9.4.1 with PostGIS version 2.1.5 on a Mac OS X 10.10.2.

Best Answer

I'm guessing that you may have a few invalid geometries in your base table. You will need to do a bit of investigation.

Another thing to do would be select the well known text (ST_AsText) from your created table then look for anomalies, eg POINTs and LINEs in the collections. You might want to add a few examples to your question.

Firstly I would make sure that your base geometries are valid

UPDATE gis_schema.sl_adm2
SET geom = ST_MakeValid(geom)
WHERE ST_IsValid(geom) IS FALSE

Then when you build your table add a check in to make sure only polygons are brought back. You may want to consider using ST_Union rather than ST_Collect, but the depends on the result you are trying to achieve. I recall that MULTIPOLYGONs can't having touching or overlapping parts, but I can't find a reference to it at the moment. Even though it is slower ST_UNION will likely give you better results. Be aware that it will dissolve boundaries between touching polygons.

CREATE TABLE example_geom AS 
SELECT sl_adm2_sums_alias.* , adm2_code as district_pcode 
FROM gis_schema.sl_adm2
JOIN (
    SELECT ST_Multi(st_collect(sl_adm2.geom)) AS geometry, sum(sl_sos.HEALTH) AS health_sum,
        sum(sl_sos.PROTECTION) AS protection_sum, sum(sl_sos.WASH) AS wash_sum,
        sum(sl_sos.cccm) AS cccm_sum, sum(sl_sos.food) AS food_sum, 
        sum(sl_sos.livelihood) AS livelihood_sum, 
        sum(sl_sos.education) AS education_sum, 
        sum(sl_sos.agriculture) AS agriculture_sum, 
        sum(sl_sos.HUMANITARIAN_ACCESS) AS humaccess_sum, 
        sum(sl_sos.logistics) AS logistics_sum, 
        sum(sl_sos.SHELTER_AND_NFI) AS shelter_and_wifi_sum, sl_sos.district
    FROM sl_sos 
    JOIN gis_schema.sl_adm2 ON sl_sos.DISTRICT_PCODE = sl_adm2.adm2_code
    WHERE ST_GeometryType(geom) in ('ST_Polygon','ST_MultiPolygon')
    GROUP BY sl_sos.district)
sl_adm2_sums_alias ON
sl_adm2.adm2_name = sl_adm2_sums_alias.district;
Related Question