ST_Intersection of 3 polygons

Question

I am new to postgis and I need to find intersection of three polygons. In my table named "mytable" contains large no of polygon coordinates and grouped into three categories by adding a field called "filename".I labelled the category names as part1, part2, part3 in filename . Now I need to find all polygon which spacially intersect in three categories. I write my code as:

SELECT ST_Intersection(a.geom, b.geom) as intersect_ab
FROM mytable a
INNER JOIN mytable b
ON ST_Intersects(a.geom,b.geom)
INNER JOIN mytable c
ON ST_Intersects(a.geom,c.geom)
WHERE ST_Overlaps(a.geom, b.geom)
AND ST_Overlaps(b.geom, c.geom)
AND ST_Overlaps(a.geom, c.geom)
AND ST_isvalid(a.geom)='t'
AND ST_isvalid(b.geom)='t'
AND ST_isvalid(c.geom)='t'
AND a.filename = 'part1'
AND b.filename = 'part2'
AND c.filename = 'part3');

But it returns wrong answer

Answer 1

It sounds like you should probably split your data into three tables first. But if there's a good reason not to do that, then you can use a WITH statement to collect them into table-like entities. If you don't mind me saying, your existing query is quite confusing, so I'm not 100% clear what you are doing. Perhaps if you explain why you want to do this, it would help people answer more thoroughly?

Still, given your existing structure, I think you want a query like this:

WITH a AS (
    SELECT * FROM mytable
    WHERE filename = 'part1'
    AND ST_IsValid(geom::geometry)
    LIMIT 1000
), b AS (
    SELECT * FROM mytable
    WHERE filename = 'part2'
    AND ST_IsValid(geom::geometry)
    LIMIT 1000
), c AS (
    SELECT * FROM mytable
    WHERE filename = 'part3'
    AND ST_IsValid(geom::geometry)
    LIMIT 1000
)
SELECT a.*, b.*, c.*
FROM a, b, c
WHERE ST_Intersects(a.geom, c.geom)
AND ST_Intersects(a.geom, b.geom)
AND ST_Intersects(b.geom, c.geom);

This gets much more complicated with each additional category, and there are certainly better ways to aggregate features like this. This is also a surefire way to issue an incredibly expensive query, so I've put in LIMIT 1000 statements that you can remove when you are ready for the full query.

Answer 2

Building on what @Richard Law posted earlier, along with the ideas here you should be able to get the intersecting polygons with the following. As with his answer, I'm not sure that this is the most efficient way to accomplish the task.

Since you say you are new to PostGIS, I'll throw out a quick reminder that you want to put a spatial index on your geometry column to improve performance.

WITH a AS (
    SELECT * FROM mytable
    WHERE filename = 'part1'
    AND ST_IsValid(geom::geometry)
    LIMIT 1000
), b AS (
    SELECT * FROM mytable
    WHERE filename = 'part2'
    AND ST_IsValid(geom::geometry)
    LIMIT 1000
), c AS (
    SELECT * FROM mytable
    WHERE filename = 'part3'
    AND ST_IsValid(geom::geometry)
    LIMIT 1000
)

SELECT
    CASE
    WHEN ST_CoveredBy(ab.geom, c.geom)
        THEN ab.geom
    WHEN ST_CoveredBy(c.geom, ab.geom)
        THEN c.geom
    ELSE
        ST_Multi(ST_Intersection(ab.geom, c.geom))
    END as geom
    FROM c,   
        (SELECT 
            CASE
            WHEN ST_CoveredBy(a.geom, b.geom)
                THEN a.geom
            WHEN ST_CoveredBy(b.geom, a.geom)
                THEN b.geom
            ELSE
                ST_Multi(ST_Intersection(a.geom, b.geom))
            END AS geom
        FROM a, b
        ) AS ab;