Accurate and fast nearest neighbour, which distance function should I use and how to it speed up?

Question

I am solving the following problem:

I have a map with many points in it and I need to find the nearest point to given location. Each point carries information about it's latitude and longitude. This is for the routing algorithm to find an origin and destination point closest to the user selected locations.

I wanted to use the "<->" operator for finding the nearest neighbour. But when testing I encountered the difficulty. It seems that the <-> is not always working properly for ordering points by distance.

When running the following query, for the first two SELECT pairs it seems that point p1 is closer to location L. While from other two query pairs the result shows that in fact p2 is closer to location L. I checked that p1 is closer to p2 using http://www.movable-type.co.uk/scripts/latlong.html. This site return the distance in km similar to the result of the third SELECT pair, which is correct, they are both using spherical model of the Earth.

By p1 I mean 14.452601,50.076207. By p2 I mean 14.452789,50.075481. By L I mean 14.451722068277277,50.07564315647511 (longitude,latitude).

Is this correct, or am i making some mistake?

SELECT
(SELECT ST_Distance(st_setsrid(st_makepoint(14.452601,50.076207),4326),st_setsrid(st_makepoint(14.451722068277277,50.07564315647511),4326)))AS geometrical_from_p1,
(SELECT ST_Distance(st_setsrid(st_makepoint(14.452789,50.075481),4326),st_setsrid(st_makepoint(14.451722068277277,50.07564315647511),4326)))AS geometrical_from_p2,

(SELECT (st_setsrid(st_makepoint(14.452601,50.076207),4326) <-> st_setsrid(st_makepoint(14.451722068277277,50.07564315647511),4326)))AS geometrical_from_p1,
(SELECT (st_setsrid(st_makepoint(14.452789,50.075481),4326) <-> st_setsrid(st_makepoint(14.451722068277277,50.07564315647511),4326)))AS geometrical_from_p2,

(SELECT ST_Distance_sphere(st_setsrid(st_makepoint(14.452601,50.076207),4326),st_setsrid(st_makepoint(14.451722068277277,50.07564315647511),4326)))AS geomethrical_from_p1,
(SELECT ST_Distance_sphere(st_setsrid(st_makepoint(14.452789,50.075481),4326),st_setsrid(st_makepoint(14.451722068277277,50.07564315647511),4326)))AS geometrical_from_p2,

(SELECT ST_Distance(ST_GeographyFromText('POINT(14.452601 50.076207)'),ST_GeographyFromText('POINT(14.451722068277277 50.07564315647511)')))AS geographical_from_p1,
(SELECT ST_Distance(ST_GeographyFromText('POINT(14.452789 50.075481)'),ST_GeographyFromText('POINT(14.451722068277277 50.07564315647511)')))geographical_from_p2
;

the 8 results for the 4 query pairs are:

0.00104424158783511; 0.00107918396178276;

0.00104493638333828; 0.00108122276233257;

88.6845186250299; 78.2447351379485;

88.8360501337395; 78.4754154775416

My application is written in java with connection to PostgreSQL 9.2 and PostGIS 2.0.3 .

Now I am using simple solution written in java which is fast and accurate enough: I just go through all the points P loaded in RAM, calculate their distance from location L and remember the one with minimal distance. I use following formula:

simple_distance = SQRT[((L.longitude-P.longitude)*LON_DEG_M)^2 + ((L.latitude-P.latitude)*LAT_DEG_M)^2]

LON_DEG_M and LAT_DEG_M represent length(in meters) of one degree of longitude/latitude at a latitude of L. I calculate them only once for L and use these values for all points P. This does not give exact distance but I use it just to find the nearest point. (SQRT is not necessary.) On my laptop it takes about 2 milliseconds to find the nearest node this way.

In the future I would like to PostGIS to calculate the nearest neighbour(s). On the internet I found the soulution suggesting using <-> operator and I tried ST_distance(geometry,geometry) function as well. As I described earlier I did not manage to make these work properly. Now I have implemented versions using ST_distance_sphere(geometry,geometry) and ST_distance(geography,geography) which seem to work properly but are slower for me.

I am using following queries:

SELECT osmid FROM nearest_node_geom_sphere ORDER BY ST_Distance_sphere(geometry,st_setsrid(st_makepoint(L.longitude,L.latitude),4326)) LIMIT 1;

SELECT osmid FROM nearest_node ORDER BY ST_Distance(geography,ST_GeographyFromText('POINT(L.longitude L.latitude)')) LIMIT 1;

Is there a way to use the <-> operator properly to return the order the distances. If no how can I improve the performance of ST_distance_sphere(geometry,geometry) and ST_distance(geography,geography). Which should be fastest to use?

For the current number of nodes I am getting the following times: (2 milliseconds for the simple java algorithm described earlier, 50 ms for St_distance_sphere(geometry,geometry), 90 ms for ST_distance(geography,geography)). I have around 9000 points in db running locally.

Welcome to gis.SE. You might want to update your question to explain what tool(s) you are trying to use here. Also, can you explain why you think P2 is closer to L than P1 (e.g. how you are calculating the "distance") — BradHards, Apr 14 '13 at 00:57
Thank you. I edited my question and specified in more detail. — Marcel, Apr 14 '13 at 06:53
don't cross post: http://stackoverflow.com/questions/15996600/fast-but-accurate-nearest-neighbour-which-postgis-distance-function-should-i-us — Mike T, Apr 14 '13 at 22:03

Jakub Kania · Accepted Answer · 2013-04-14T11:36:45.447

If you want to use the operator <-> you will have to choose the correct projection for your data. It seems your points are in Czech Republic so I chose the first I found SRID 2065

SELECT 
(SELECT ST_transform(st_setsrid(st_makepoint(14.452601,50.076207),4326),2065)
 <-> st_transform(st_setsrid(st_makepoint(14.451722068277277,50.07564315647511),4326),2065))AS geometrical_from_p1,
(SELECT st_transform(st_setsrid(st_makepoint(14.452789,50.075481),4326),2065) <-> st_transform(st_setsrid(st_makepoint(14.451722068277277,50.07564315647511),4326),2065) )AS geometrical_from_p2;

88.7882400011398;78.4504471386238

It is a few hundred meters off from the geog result but you might find a better projection. Transform the points in your table, create index on the geometry column and then try:

SELECT osmid 
FROM nearest_node_geom_sphere 
ORDER BY geometry <-> st_transform(st_setsrid(st_makepoint(L.longitude,L.latitude),4326)),2065) ASC 
LIMIT 1;

It will use an index but I think it's impossible it will be faster than 2ms you get in Java.

Accurate and fast nearest neighbour, which distance function should I use and how to it speed up?

1 Answers1