[postgis-devel] added Azimuth(geometry, geometry)

[strk at refractions.net]

I've added a C port of Hartmut's plpsgql
Azimuth(point, point) function in the head branch.

Contrary to Hartmut version this one:
	- doesn't attempt at reprojecting input
	- returns radians
	- checks for SRID match

--strk;

On Thu, Sep 15, 2005 at 06:27:49PM +0900, Hartmut Tschauner wrote:
> 
> Hi Ahmet,
> In an earlier message that didn't arrive as part of this thread and with its
> body scrapped, I posted a function (see below) that will calculate azimuth
> between two points. You would pass the start and end points to get the
> direction of a line.
> Cheers,
> Hartmut
> Hartmut Tschauner
> Department of Archaeology
> Seoul National University
> Gwanak-Gu Sillim 9-dong San 56-1
> Seoul, 151-742
> Korea
> ph. +82 (2) 880-9260
> 
> 
> CREATE OR REPLACE FUNCTION public.azimuth(geometry, geometry)
>   RETURNS float8 AS
> '
> DECLARE
>    geom1 ALIAS FOR $1;
>    geom2 ALIAS FOR $2;
>    geom2trans geometry;
> 
> BEGIN
>    IF geom1 IS NULL OR geom2 IS NULL THEN
>       RETURN NULL;
>    ELSE
>       IF isempty(geom1) OR isempty(geom2) OR geometrytype(geom1) !=
> \'POINT\' OR geometrytype(geom2) != \'POINT\' THEN 
>          RETURN NULL;
>       END IF;
>    END IF;
> 
>    IF srid(geom1) != -1 AND srid(geom2) != srid(geom1) THEN
>       geom2trans := transform(geom2, srid(geom1));
>    ELSE
>       geom2trans := geom2;
>    END IF;
> 
>    IF x(geom1) = x(geom2trans) AND y(geom1) < y(geom2trans) THEN 
>       RETURN 0;
>    ELSIF x(geom1) = x(geom2trans) AND y(geom1) > y(geom2trans) THEN 
>       RETURN 180;
>    ELSIF y(geom1) = y(geom2trans) AND x(geom1) < x(geom2trans) THEN 
>       RETURN 90;
>    ELSIF y(geom1) = y(geom2trans) AND x(geom1) > x(geom2trans) THEN 
>       RETURN 270;
>    ELSIF x(geom1) < x(geom2trans) AND y(geom1) < y(geom2trans) THEN 
>       RETURN degrees(atan(abs(x(geom1) - x(geom2trans)) / abs(y(geom1) -
> y(geom2trans))));
>    ELSIF x(geom1) < x(geom2trans) AND y(geom1) > y(geom2trans) THEN 
>       RETURN degrees(atan(abs(y(geom1) - y(geom2trans)) / abs(x(geom1) -
> x(geom2trans)))) + 90;
>    ELSIF x(geom1) > x(geom2trans) AND y(geom1) > y(geom2trans) THEN 
>       RETURN degrees(atan(abs(x(geom1) - x(geom2trans)) / abs(y(geom1) -
> y(geom2trans)))) + 180;
>    ELSIF x(geom1) > x(geom2trans) AND y(geom1) < y(geom2trans) THEN 
>       RETURN degrees(atan(abs(y(geom1) - y(geom2trans)) / abs(x(geom1) -
> x(geom2trans)))) + 270;
>    ELSE
>       RETURN 0;
>    END IF;
> 
> END;
> '
>   LANGUAGE 'plpgsql' VOLATILE;
> 
> 

> _______________________________________________
> postgis-users mailing list
> postgis-users at postgis.refractions.net
> http://postgis.refractions.net/mailman/listinfo/postgis-users