Performing multiple intersections between two polygon tables and one raster table in PostGIS?

Question

Trying to perform an intersection of an intersection in PostGIS. Its an intersection between a point table, polygon table and a raster table with different number of rows.

The polygon table is the area of interest. The point table is a global soil properties file. The raster table is a an agricultural crop mask with 1 or 0 for values.

The goal was to extract soil properties where the crop mask values are 1 within the given study area. The dssat.cropland table is the raster table, kenyasq18.agareas is the polygon/area of interest and dssat.soils is the point table with props/properties.

select props,st_value 
from 
  (
    select st_value(c.rast,1,1,1) 
    from 
     dssat.cropland as c,
     kenyasq18.agareas as a 
    where 
      st_intersects(c.rast,a.geom)
  ) as cropmask, 
  dssat.soils as s, 
  kenyasq18.agareas as a 
where st_intersects(s.geom,a.geom) 
 and st_value=1;

The above query results in duplication of rows and values. Mathematically speaking I would like to do (A ∩ B ∩ C) but seems like I keep ending up with (A ∩ B) X (B ∩ C).

What would be the best way to approach such a query?

Answer 1

You can use join to link the 3 tables, avoiding the cross product you are getting:

SELECT props 
FROM 
  dssat.cropland as c 
  INNER JOIN kenyasq18.agareas as a 
    ON (st_intersects(c.rast,a.geom) AND st_value(c.rast,1,1,1)  = 1)
  INNER JOIN dssat.soils as s 
    ON (st_intersects(s.geom,a.geom) );

Answer 2

Slightly modified query from JGH. Changing the order helped in retrieving the proper results.

SELECT distinct(props) 
FROM dssat.soils as s INNER JOIN {0}.agareas as a 
    ON (st_intersects(s.geom,a.geom) and gid= {1}) 
INNER JOIN dssat.cropland as c ON (st_intersects(c.rast, s.geom) and 
st_value(c.rast,1,1,1)=1)