1

I have 2 PostGIS tables in PostgreSQL: polytable and pointtable. The point table has 1000s of points and a class column with 6 different unique classes (integer values from 1 to 6). The polygon table has 1000s of polys and a name column. Some names appear multiple times.

For every single polygon, I want to count how many points (per point class) are in each. Although the polygons contain a name column (we'll call polyname here), I don't need these grouped at all... The current query I use only works when using GROUP BY polyname on the polygon table. Therefore I can't treat each polygon individually:

SELECT
  polyname,
  count(pointclass) FILTER (WHERE pointclass = 1) AS pointclass1,
  count(pointclass) FILTER (WHERE pointclass = 2) AS pointclass2,
  count(pointclass) FILTER (WHERE pointclass = 3) AS pointclass3,
  count(pointclass) FILTER (WHERE pointclass = 4) AS pointclass4,
  count(pointclass) FILTER (WHERE pointclass = 5) AS pointclass5,
  count(pointclass) FILTER (WHERE pointclass = 6) AS pointclass6
FROM 
    (select poly.polyname, pnt.pointcolumn
    from schema.polytable poly
    left join schema.pointtable pnt
    on st_intersects(pnt.geom, poly.geom)) sub
GROUP BY polyname

How can I remove the grouping in the query above? Simply removing the last line throws an error (column "sub.polyname" must appear in the GROUP BY clause or be used in an aggregate function)

This is my desired output:

polyid   polyname   polygeom   pointclass1   pointclass2   pointclass3...
1        red        0101000..  55            10            5
2        blue       0101000..  3             25            46
3        green      0101000..  100           39            11
4        yellow     0101000..  0             0             3
5        red        0101000..  0             11            37
6        green      0101000..  58            99            200
7        green      0101000..  0             0             0

Note above, how I still want polygon rows returned which have no points within. So a LEFT JOIN of point count to the polygon, per point class. Similar to the query below, but with 6 point count columns for each point class:

SELECT poly.id, poly.name, count(point.shape) AS pointcount 
FROM polytable poly
LEFT JOIN pointtable point ON st_intersects(poly.geom, point.grom)
GROUP BY poly.id

I'd prefer not to take the longer route of adding 6 empty columns onto the polygon table for each point class, then running an UPDATE on these columns to populate them, using different WHERE pointclass = 1/2/3... etc... There must be a better way?

  • 1
    Grouping by both polyid and polyname, sould return all the rows you want. Just call the polyid column in all your select statements including the subquery. – FSimardGIS Nov 15 at 23:30
  • @FSimardGIS thanks. If you post the first query above as an answer, but add in SELECT polyid, polyname, polygeom into both the SELECT calls and also into the GROUP BY I can accept your solution as the answer. – Theo F Nov 16 at 0:04
1

The solution is to group by both polyid and polyname, it should return all your desired records. Just call the polyid column in all your select statements including the subquery. You can also call the geometry of your polygons by adding the poly.geom.

SELECT
  polyid, polyname, poly.geom
  count(pointclass) FILTER (WHERE pointclass = 1) AS pointclass1,
  count(pointclass) FILTER (WHERE pointclass = 2) AS pointclass2,
  count(pointclass) FILTER (WHERE pointclass = 3) AS pointclass3,
  count(pointclass) FILTER (WHERE pointclass = 4) AS pointclass4,
  count(pointclass) FILTER (WHERE pointclass = 5) AS pointclass5,
  count(pointclass) FILTER (WHERE pointclass = 6) AS pointclass6
FROM 
    (select poly.polyid, poly.polyname, poly.polygeom, pnt.pointclass
    from schema.polytable poly
    left join schema.pointtable pnt
    on st_intersects(pnt.geom, poly.geom)) sub
GROUP BY polyid, polyname, poly.geom

Alternatively, you can use the following simpler query, since in this case, you don't even need to use a subquery:

select
  polyid, polyname, poly.geom,
   count(pointclass) FILTER (WHERE pointclass = 1) AS pointclass1,
   count(pointclass) FILTER (WHERE pointclass = 2) AS pointclass2,
   count(pointclass) FILTER (WHERE pointclass = 3) AS pointclass3,
   count(pointclass) FILTER (WHERE pointclass = 4) AS pointclass4,
   count(pointclass) FILTER (WHERE pointclass = 5) AS pointclass5,
   count(pointclass) FILTER (WHERE pointclass = 6) AS pointclass6
FROM schema.polytable poly
left join schema.pointtable pnt
on st_intersects(pnt.geom, poly.geom)
GROUP BY polyid, polyname, poly.geom
| improve this answer | |
  • 1
    last line should read GROUP BY polyid, polyname, polygeom, but thanks – Theo F Nov 16 at 0:38
  • 1
    You're right, I edited my answer. In fact I tested it on one of my datasets and some of my fields name remained on the last line :-) – FSimardGIS Nov 16 at 1:16
  • Grouping by geometry isn't ideal, since the equality test can be expensive. It would be nice if Postgres provide a FIRST aggregate function for this kind of query (there is an implementation on the wiki) – dr_jts Nov 16 at 20:03

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.