[GIS] ST_Within and ST_Contains not working for Latitude and Longitude in PostGIS

geolocationgeometrypostgispostgresql

I have a table "tbl_location" in PostGIS which has latitude and longitude, so I created a new column called "point_location" :

update "tbl_location"  set point_location = 
ST_Transform(ST_SetSRID(ST_MakePoint(lng, lat), 4326), 26913);

I have another table which is available here has MultiPolygon field and I am trying to find in which polygon a location(lat/lng) belongs. For that, I wrote the following query:

SELECT a.id, b.csa2010
FROM baltimore."vs14_housing" b, baltimore."tbl_location" a 
WHERE 
ST_Intersects(a.point_location,  b.geom) 
    or 
ST_Contains(a.point_location, b.geom)
    or
ST_Within(a.point_location, b.geom)

It is not returning any records. Data type for geom column in vs14_housing is geometry (MultiPolygon). Let me know what I am doing wrong.

Best Answer

Are you using the correct spatial reference ID? The vs14_housing data you referenced is 102685.

If you use your query using the UTM(26913), no results.

But if you change it to 102685, it works. Unless you projected your vs14_housing data (which you failed to mention in the question).

with cte as 
(
    select 1 as id, ST_Transform(ST_SetSRID(ST_MakePoint(-76.6129, 39.3106), 4326), 102685) point_location
    ,st_astext(ST_Transform(ST_SetSRID(ST_MakePoint(-76.6129, 39.3106), 4326), 102685)) wkt 
)
SELECT a.id, b.csa2010
FROM public."vs14_housing" b, cte a 
WHERE 
ST_Intersects(a.point_location,  b.geom) 
    or 
ST_Contains(a.point_location, b.geom)
    or
ST_Within(a.point_location, b.geom)

Related Solutions

[GIS] Problem PostGIS Spatial Query ST_Contains ST_Disjoint ST_Intersects

You have a set of polygon, they cover an area of the map.

You have a set of points, some are inside the polygons and some are outside.

The st_intersect will compare all points against all polygons and only return points that are "on top" one or more polygons.

The st_disjoint will return a point if a point is not "on top" of a polygon. So for a set of points against a set of polygons, so unless the polygons overlap, each point is going to be "outside" of one or more polygons.

You can have a nested query of the intersect, and another inverted nested query of the intersect to get your t/f result

Adding a solution

There's a fiddle that goes with this here, since it was a PostgreSQL fiddle rather than PostGIS I've change the fiddle test to be a simple integer equality test.

So the first half of the sql is your intersect test:

    SELECT distinct points.id
    FROM points, multipolygon
    WHERE st_intersects(points.geom, multipolygon.geom)

I add a TRUE column to match your sample, by adding a wrapping SQL select with a hardcoded value of TRUE:

SELECT
   id,
   TRUE AS intersects_at_least_one_polygon
FROM 
(
    SELECT distinct points.id
    FROM points, multipolygon
    WHERE st_intersects(points.geom, multipolygon.geom)
) intersecting ;

To find the points that are NOT intersecting the polygons, well that is the total list of point id, minus the above points. So I select all the points, the outer select below, and then say "WHERE id NOT IN" and just repeat your SQL, plus I add a hardcode a FALSE:

SELECT
   points.id,
   FALSE AS intersects_at_least_one_polygon
FROM
   points WHERE id NOT IN
   (
        SELECT distinct points.id
        FROM points, multipolygon
        WHERE st_intersects(points.geom, multipolygon.geom)
   )
   order by id
   ;


SELECT
   id,
   TRUE AS intersects_at_least_one_polygon
FROM 
(
    SELECT distinct points.id
    FROM points, multipolygon
    WHERE st_intersects(points.geom, multipolygon.geom)
) intersecting
UNION
SELECT
   points.id,
   FALSE AS intersects_at_least_one_polygon
FROM
   points WHERE id NOT IN
   (
        SELECT distinct points.id
        FROM points, multipolygon
        WHERE st_intersects(points.geom, multipolygon.geom)
   )
   order by id
   ;

I'm certain that there are more elegant, efficient, solutions to this, but hope you're happy with a simple working version.

[GIS] How to make projected geometry and update it to postgresql

As user30184 outlined: It's a very common task, and there's plenty of documentation: ST_Transform.

To apply this, you need to figure out the EPSG codes of your projections.

UTM 35N probably is EPSG:32635, and your lat/lon coordinates could be anything. One of the more common ones is WGS84 as used in GPS with code EPSG:4326

The Postgis query then would be

SELECT ST_Transform(geom, 4326,32635)

If your postgis table already has the proper source SRID, you can do a simpler

SELECT ST_Transform(geom, 32635)

EDIT after you updated your question a lot:

The issue is not with ST_Transform, the issue is your query. You explicitly set all geometries to be the same. This is how UPDATE in combination with subqueries works. I really recommend you dig through the documentation a bit more to grasp this.

The solution to this issue can be found in the docs as well, by using a WHERE clause at the end. While not standard SQL, PostgreSQL allows UPDATE WHERE in order to apply subqueries to specific rows they match. In most cases, one has a unique identifier that you can use (just hand it down from the sub queries). In your case, depending on your data, you could use the timestamp or the lon/lat/alt itself to do so, as all these should match the same anyway.:

UPDATE info SET geom = sq.geom 
FROM (
    SELECT lon, lat,alt, ST_Transform(geom,32635) as geom
    FROM (
        SELECT lon,lat,alt, ST_setsrid(geom_v,4326) as geom 
        from (
            select lon,lat,alt, ST_MakePoint(lon,lat,alt) as geom_v 
            from (
                select latitude as lat, longitude as lon, altitude as alt from info order by gid asc
            ) as fpp
        )as ftt
    ) as tr
) as sq
WHERE info.lon = sq.lon AND info.lat = sq.lat AND info.alt = sq.alt

However, your entire query is unnecessarily convoluted (extremely so!). You could (and should) just use this:

UPDATE test SET geom = ST_Transform(ST_SetSRID(ST_MakePoint(lon,lat, alt),4326),32635)

to achieve the same result with a much faster and easier-to-read query.

Best Answer

Related Solutions

[GIS] Problem PostGIS Spatial Query ST_Contains ST_Disjoint ST_Intersects

[GIS] How to make projected geometry and update it to postgresql

Related Question