PostGIS Optimization – Getting osm_id from Random Point with OpenStreetMap Database

openstreetmapoptimizationpostgis

I've downloaded the OpenStreetMap data from GeoFabrik and imported it into PostgreSQL as described in the first procedures on this page.

The import was successful. After that, I found a way of associating random points with its osm_id from this database. I know there are APIs like https://nominatim.openstreetmap.org/reverse.php?lat=-23.513442&lon=-46.384794&zoom=18&format=jsonv2 from Nominatim where I can get this information, but since there're limitations on the number of requests that we can make. I'd like to implement a solution to get this information locally with PostGIS.

I've seen this answer that solves the problem of how to find the closest road to the point… I've adapted it to my case as:

SELECT 
    way,
    osm_id,
    ST_Distance(way, 'SRID=3857;POINT(-3890646.5744145643 -899377.0801662721)'::geometry) 
FROM planet_osm_line r 
ORDER BY 2 ASC 
LIMIT 1;

It works fine and it returns the geometry and osm_id of the closest road from the point. However, even though the query returns after some seconds I'm in doubt if there are ways of optimizing the time of these results since in my case of use I want it to have a fast response to a large number of input points that arrive at the same time.

Should I worry about fragmenting the GeoFabrik data into smaller datasets (like a dataset with only the city that I'm interested) to have a faster response with the last query?

What strategies can I use to get a fast return of the osm_id and geometry from the OpenStreetMap database?

Best Answer

Use the 2D distance operator <->:

SELECT
    way,
    osm_id
FROM
    planet_osm_line AS osm
ORDER BY
    osm.way <-> 'SRID=3857;POINT(-3890646.5744145643 -899377.0801662721)'::GEOMETRY
LIMIT
    1
;

Related:

Related Solutions

[GIS] ST_ClosestPoint(line, point) does not intersect line

This is a problem coming from that ST_Intersects has no tolerance. Even if double precision coordinates holds a lot of decimals they form a grid where the only places for the points is in the crossings. Often the line doesn't intersect with any of those crosses and there is no way any point will intersect the line exactly. The workaround is to use st_dwithin instead with a small tolerance.

/Nicklas

[GIS] Calculating distance of point from arc/polyline using PostGIS

SELECT arc_id 
FROM road_table 
ORDER BY ST_Distance(road_table.geom, 'gps-point in wkt') 
LIMIT 1;

Since PostGIS 1.5 you also have the option of using ST_ClosestPoint.

Some discussions of how to use it can be found here:

http://blog.jordogskog.no/2010/02/07/how-to-use-the-new-distance-related-functions-in-postgis-part1/

About ST_ClosestPoint you can read more in the PostGIS documentation: http://postgis.net/docs/ST_ClosestPoint.html

As I understand your query correct you don't have any line ready but a list of points in a table.

What you need to do then is writing a query that creates lines from those points and use ST_ClosestPoint on those lines. You do it all in one query. To do that you can design your query in several ways and approaches depending on in more detail what you want to get. If it is just the distance you want, then forget about ST_ClosestPoint then use ST_Distance to get the distance. If you just want the closest Point on the edge and/or the distance, then you can build the whole linestring with ST_MakeLine If you also want to find the points defining the edge that is closest it is probably easiest to make a self join that builds lines from each opint pair in the table instead and use ST_Distance and /or ST_ClosestPoint on that two point lines.

This I guess looks quite messy, but that is because there is a lot of possibilities, what you actually want to get.

I think edge is a more correct word than arc.

Ok, from your comment I guess that what you want is a query that self-joins to create edges and returns the closest edge and distance. You will need a integer field to make it work like this with integer values ordered as the gps-points and without holes in the series. I call id gid here.

Then you can run something like:

SELECT dist, the_edge FROM
(SELECT ST_Distance(e.the_edge, 'PUT YOUR GPS_POINT AS WKT HERE') as dist, the_edge FROM
(SELECT ST_MakeLine(a.geom, b.geom) AS the_edge FROM
point_table a INNER JOIN point_table b ON b.gid=(a.gid+1)) e) s
ORDER BY dist LIMIT 1;

a.geom references geom in the table that I have put the alias a on. I shouldn't have left out the "AS" key-word, then it looks like this:

SELECT dist, the_edge FROM
(SELECT ST_Distance(e.the_edge, 'PUT YOUR GPS_POINT AS WKT HERE') as dist, the_edge FROM
(SELECT ST_MakeLine(a.geom, b.geom) AS the_edge FROM
point_table AS a INNER JOIN point_table AS b ON b.gid=(a.gid+1)) e) s
ORDER BY dist LIMIT 1;

To get a short intro to sql in general you can test the tutorial on PostGISonline

To see some spatial examples you can try: mixed examples

Best Answer

Related Solutions

[GIS] ST_ClosestPoint(line, point) does not intersect line

[GIS] Calculating distance of point from arc/polyline using PostGIS

Related Question