[postgis-users] a line' s direction

[Hartmut Tschauner]

Skipped content of type multipart/alternative-------------- next part --------------
-- Function: public.azimuth(geometry, geometry)

-- DROP FUNCTION public.azimuth(geometry, geometry);

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;