Add columns if layer intersects more than one other layer in PostGIS

Question

I have been spatially joining land use zones to parcels successfully with the following query:

CREATE TABLE new_parcel AS
SELECT DISTINCT ON(t.geom) t.*, m.basezone, m.overlay
FROM parcel AS t
LEFT JOIN zoning AS m
  ON st_intersects(t.geom, m.geom)
  AND st_area(st_intersection(t.geom, m.geom)) / st_area(t.geom) > 0.5

Currently, more than 50% of the parcel's area must be within a single zone for a join to occur.

However, it turns out that some of the parcels overlap 3 or more zones and are not getting joined to the zoning table because less than 50% of their area overlaps. For example, the parcels of dark purple color overlap more than two zones:

How would one modify the query to add extra columns in the resulting table which would include data from all zones overlapped?

The existing code produces a table that has

basezone | overlay

How could one produce

basezone_2|overlay_2|basezone_3|overlay_3|...

and so on for all zones overlapped?

Answer 1

What you're asking for isn't easy because you don't know how many columns to produce until after the query is complete. crosstab from tablefunc will permit some of this, but because you're not providing any sample data it's damn near impossible.

As a mater of the join,

In order to first tell which zones overlap multiple zones with more than 50% of their area, you have to

1. Calculate all zones that overlap with a join of parcel/zoning,
2. Calculate the area of the overlap.

Then afterward, you can see if the combined overlap is greater tan 50%.

A join operation works row-wise. There is no concept of a join based on a function of the other table the row. Modifying your query, it would look something like this...

First, you need the

SELECT t.parcel
FROM (
  SELECT p.geom AS parcel,
    m.geom AS zoning,
    st_area(st_intersection(p.geom, m.geom)) AS areaintersection
  FROM parcel AS p
  JOIN zoning AS m
    ON st_intersects(p.geom, m.geom)
) AS t
GROUP BY t.parcel
HAVING sum(areaintersection) / st_area(t.geom) > 0.5;

Now, should be a result set of parcel.geom. From there, you can INNER JOIN back to parcel, and to zoning to get a list that you can feed to crosstab().