Thread: BUG #2037: user function call unexpected "input out of range"

The following bug has been logged online:

Bug reference:      2037
Logged by:          Tom
Email address:      lackey@ltu.edu
PostgreSQL version: 8.1-beta4
Operating system:   Windows XP
Description:        user function call unexpected "input out of range"
Details:

--No database needs to be selected.
--Just cut and paste into a file.
--All remarks are in comments.
--
--From the psql prompt use "\i fileName".
--
create function  find (numeric,numeric,numeric,numeric) returns double
precision
as
    $$

        --Calculate distance
        select     acos
            (
                sin($1*pi()/180)*sin($3*pi()/180)
                +
                cos($1*pi()/180)*cos($3*pi()/180)
*cos(($2-$4)*pi()/180)
            )*60*1.1515*180/pi();



    $$
language SQL;

--Number 1.
select find(42.3202,-83.2687,42.3202,-83.2687);
--Number 2.
select find(42.320,-83.2687,42.3202,-83.2687);
--Number 3.
select find(42.3338,-83.1818,42.3338,-83.1818);

\df+ find;
drop function find(numeric,numeric,numeric,numeric);


/* This surfaced when calculating the distances between zip codes using
latidtudes and longitudes. When the latitudes and longitudes are the same
the distance should be zero. Number 1 gives the error. Number 2 differs by a
digit in the last decimal place and works. Number 3 is another test with
identical latitudes and longitudes and does calculate zero.

The "DROP FUNCTION" was only put in so I could repeatly re-execute the
script.


The ".msi" was used for the install accepting all of the defaults.

Thanks, Tom

*/

"Tom" <lackey@ltu.edu> writes:
> Description:        user function call unexpected "input out of range"

Why does this surprise you?  Floating point computation is inherently
inexact, so coming out with a value fractionally greater than 1 for
the acos() argument doesn't seem all that unlikely.  You probably
ought to add some code to clamp the result to the legal range.

            regards, tom lane

On Fri, Nov 11, 2005 at 09:26:47PM +0000, Tom wrote:
> This surfaced when calculating the distances between zip codes using
> latidtudes and longitudes. When the latitudes and longitudes are the same
> the distance should be zero. Number 1 gives the error. Number 2 differs by a
> digit in the last decimal place and works. Number 3 is another test with
> identical latitudes and longitudes and does calculate zero.

If you use psql and increase VERBOSITY you can see where the error
is happening:

test=> \set VERBOSITY verbose
test=> select find(42.3202,-83.2687,42.3202,-83.2687);
ERROR:  22003: input is out of range
CONTEXT:  SQL function "find" statement 1
LOCATION:  dacos, float.c:1602

"dacos" is the internal name of PostgreSQL's acos() function.  I'd
guess the value that acos() receives is a wee bit bigger than 1.0
due to rounding and the inaccuracy of representing floating point
numbers in binary; you don't see the problem with another pair of
identical locations because the value calculated to the 16th or so
digit is slightly different.

Your function uses the law of cosines; consider using haversine
instead.

http://www.movable-type.co.uk/scripts/GIS-FAQ-5.1.html

If you're working with geospatial data then you might want to look
at PostGIS.

http://postgis.refractions.net/

See also PostgreSQL's contrib/earthdistance module.

--
Michael Fuhr

On Fri, Nov 11, 2005 at 09:26:47PM +0000, Tom wrote:
>=20
> The following bug has been logged online:
>=20
> Bug reference:      2037
> Logged by:          Tom
> Email address:      lackey@ltu.edu
> PostgreSQL version: 8.1-beta4
> Operating system:   Windows XP
> Description:        user function call unexpected "input out of range"
> Details:=20

[ basically acos(sin(w1)*sin(w2) + cos(w1)*cos(w2)*cos(l1-l2)) ]

Seems to me that you are hitting an unfortunate roundoff error, where
the argument to acos is slightly greater than 1 (it will be
mathematically 1 in the case w1=3D=3Dw2 and l1=3D=3Dl2).

Maybe you should limit the argument to acos to -1..+1?

regards
-- tom=C3=A1s

On Fri, Nov 11, 2005 at 21:26:47 +0000,
  Tom <lackey@ltu.edu> wrote:
>
>         --Calculate distance
>         select     acos
>             (
>                 sin($1*pi()/180)*sin($3*pi()/180)
>                 +
>                 cos($1*pi()/180)*cos($3*pi()/180)
> *cos(($2-$4)*pi()/180)
>             )*60*1.1515*180/pi();

In addition to the other comments, you don't want to calculate distance this
way. It isn't very accurate when $2 and $4 are nearly equal which is the
normal case.
Do a search for haversine.