Re: my float function returns NaN

Dear Frank,

My guess would be that due to rounding errors, the argument passed to acos() is outwith the interval [-1,1] and acos() is only defined within this interval (for reals). I.e the expression:

sin(theta)^2 + cos(theta)^2*cos(0)

is coming out ever so slightly greater than 1, when it should, of course, be exactly 1.

You could simply insert a check prior to calling acos() to ensure it's argument is within the acceptable range. Eg, in pseudo code:

x = sin(lat1)*sin(lat2)+cos(lat1)*cos(lat2)*cos(dlon)
if x > 1 then x=1 else if x< -1 then x = -1
acos(x)

But I'd strongly suggest you investigate using the haversine formulation for spherical distances as it is much better for small distances. Your version is mathematically sound, but computationally prone to rounding errors.

See http://en.wikipedia.org/w/index.php?title=Great-circle_distance&useformat=desktop

You may still need to check the input to asin() is within the interval [-1,1].

Nathaniel

On 7 Apr 2012, at 03:29, Frank Bax <fbax(at)sympatico(dot)ca> wrote:

> I have a function to calculate distance between two points. When asking for distance between two points that are actually the same point; my application (which rounds to only very few decimal places), might return zero or NaN depending on decimal accuracy of point.
>
> I am able to reproduce the problem with the following SQL statements.
>
> CREATE OR REPLACE FUNCTION geo_deg2rad(float) RETURNS float AS $$ SELECT ($1 * pi()::float / 180::float) $$ LANGUAGE SQL IMMUTABLE
>
> CREATE OR REPLACE FUNCTION geo_dist(point,point) RETURNS float AS $$ SELECT 6371.0 * acos( sin(geo_deg2rad($1[0]))*sin(geo_deg2rad($2[0])) + cos(geo_deg2rad($1[0]))*cos(geo_deg2rad($2[0])) * cos(geo_deg2rad($2[1])-geo_deg2rad($1[1])) ) $$ LANGUAGE SQL IMMUTABLE
>
> create table geo (p point);
> insert into geo values (point(44.85051666667,-79.5405));
> insert into geo values (point(44.850516667,-79.5405));
>
> select p,geo_dist(p,p) from geo;
> p | geo_dist
> ---------------------------+----------------------
> (44.85051666667,-79.5405) | NaN
> (44.850516667,-79.5405) | 0.000134258785931453
>
> Why does the point with more decimal accuracy result in NaN?
>
> If its not easy to "fix" my function; is the a function that can convert NaN to zero (like coalesce(p,0) would convert null to zero)?
>
>
> PostgreSQL 9.0.4 on x86_64-unknown-openbsd5.0, compiled by GCC cc (GCC) 4.2.1 20070719 , 64-bit
>
> --
> Sent via pgsql-novice mailing list (pgsql-novice(at)postgresql(dot)org)
> To make changes to your subscription:
> http://www.postgresql.org/mailpref/pgsql-novice