From: Mark Leslie <mark.leslie@lisasoft.com>
Date: Fri, 4 Mar 2005 19:31:36 +0000 (+0000)
Subject: Cleaned up distance_sphere and updated documentation.
X-Git-Tag: pgis_1_0_0RC4~58
X-Git-Url: https://granicus.if.org/sourcecode?a=commitdiff_plain;h=6a68925c532a615e9c074676bc7c9029ce184432;p=postgis

Cleaned up distance_sphere and updated documentation.


git-svn-id: http://svn.osgeo.org/postgis/trunk@1492 b70326c6-7e19-0410-871a-916f4a2858ee
---

diff --git a/CHANGES b/CHANGES
index 1acbe9d3f..74ad6ac3e 100644
--- a/CHANGES
+++ b/CHANGES
@@ -11,6 +11,7 @@ PostGIS 1.0.0
 	- build scripts refinements
 	- MultiLine handling BUG fix in loader and dumper
 	- improved version handling, central Version.config
+        - added distance_sphere function 
 
 PostGIS 1.0.0RC3
 2005/02/24
diff --git a/doc/postgis.xml b/doc/postgis.xml
index eecaa6711..bec1db92e 100644
--- a/doc/postgis.xml
+++ b/doc/postgis.xml
@@ -3834,7 +3834,8 @@ dimension
             <para>Returns linear distance in meters between two lat/lon
             points. Uses a spherical earth and radius of 6370986 meters.
             Faster than <link linkend="distance_spheroid">distance_spheroid()</link>, but
-            less accurate.</para>
+            less accurate.
+			Only implemented for points.</para>			
           </listitem>
         </varlistentry>
 
diff --git a/lwgeom/lwgeom_spheroid.c b/lwgeom/lwgeom_spheroid.c
index 8a302917f..6b9f85f95 100644
--- a/lwgeom/lwgeom_spheroid.c
+++ b/lwgeom/lwgeom_spheroid.c
@@ -547,16 +547,12 @@ Datum LWGEOM_distance_sphere(PG_FUNCTION_ARGS)
 	double sino;
         double result;
 
-        //elog(NOTICE, "LWGEOM_distance_sphere called\n");
-
 	if (pglwgeom_getSRID(geom1) != pglwgeom_getSRID(geom2))
 	{
 		elog(ERROR, "LWGEOM_distance_sphere Operation on two GEOMETRIES with differenc SRIDs\n");
 		PG_RETURN_NULL();
 	}
 
-        //elog(NOTICE, "LWGEOM_distance_sphere passed srid check: %i, %i\n", pglwgeom_getSRID(geom1), pglwgeom_getSRID(geom2));
-
         lwpt1 = lwpoint_deserialize(SERIALIZED_FORM(geom1));
 	if (lwpt1 == NULL)
 	{
@@ -571,20 +567,12 @@ Datum LWGEOM_distance_sphere(PG_FUNCTION_ARGS)
 		PG_RETURN_NULL();
 	}
 
-        //elog(NOTICE, "LWGEOM_distance_sphere passed point check\n");
-
         pt1 = palloc(sizeof(POINT2D));
         pt2 = palloc(sizeof(POINT2D));
 
         lwpoint_getPoint2d_p(lwpt1, pt1);
         lwpoint_getPoint2d_p(lwpt2, pt2);
 
-        //elog(NOTICE, "LWGEOM_distance_sphere created point2d objects: %f, %f %f, %f\n", pt1->x, pt1->y, pt2->x, pt2->y);
-
-        //elog(NOTICE, "LWGEOM_distance_sphere beginning calculation\n");
-        //elog(NOTICE, "LWGEOM_distance_sphere TWO_PI = %f\n", TWO_PI);
-        //elog(NOTICE, "LWGEOM_distance_sphere M_PI = %f\n", M_PI);
-
 	/*
 	 * Start geo_distance code.  Longitude is degrees west of
 	 * Greenwich, and thus is negative from what normal things
@@ -596,8 +584,6 @@ Datum LWGEOM_distance_sphere(PG_FUNCTION_ARGS)
 	long2 = -2 * (pt2->x / 360.0) * M_PI;
 	lat2 = 2 * (pt2->y / 360.0) * M_PI;
 
-        //elog(NOTICE, "LWGEOM_distance_sphere normalized: %f, %f %f, %f", long1, lat1, long2, lat2);
-
         /* compute difference in longitudes - want < 180 degrees */
         longdiff = fabs(long1 - long2);
         if (longdiff > M_PI)
@@ -605,8 +591,6 @@ Datum LWGEOM_distance_sphere(PG_FUNCTION_ARGS)
                 longdiff = (2 * M_PI) - longdiff;
         }
 
-        //elog(NOTICE, "LWGEOM_distance_sphere longdiff = %f\n", longdiff);
-
         sino = sqrt(sin(fabs(lat1 - lat2) / 2.) * sin(fabs(lat1 - lat2) / 2.) +
                 cos(lat1) * cos(lat2) * sin(longdiff / 2.) * sin(longdiff / 2.));
         if (sino > 1.)
@@ -614,11 +598,8 @@ Datum LWGEOM_distance_sphere(PG_FUNCTION_ARGS)
                 sino = 1.;
         }
 
-        //elog(NOTICE, "LWGEOM_distance_sphere sino = %f\n", sino);          
-
         result = 2. * EARTH_RADIUS * asin(sino);
         
-        //elog(NOTICE, "LWGEOM_distance_sphere finished calculation: %f\n", result);
         pfree(pt1);
         pfree(pt2);
 
diff --git a/lwgeom/lwpostgis.sql.in b/lwgeom/lwpostgis.sql.in
index e2fb60a75..5fa00b0c5 100644
--- a/lwgeom/lwpostgis.sql.in
+++ b/lwgeom/lwpostgis.sql.in
@@ -1558,41 +1558,11 @@ CREATEFUNCTION distance_spheroid(geometry,geometry,spheroid)
 	AS '@MODULE_FILENAME@','LWGEOM_distance_ellipsoid_point'
 	LANGUAGE 'C' _IMMUTABLE_STRICT; -- WITH (isstrict);
 
-CREATEFUNCTION distance_sphere_latlon(geometry,geometry)
+CREATEFUNCTION distance_sphere(geometry,geometry)
 	RETURNS FLOAT8
 	AS '@MODULE_FILENAME@','LWGEOM_distance_sphere'
 	LANGUAGE 'C' _IMMUTABLE_STRICT; -- WITH (isstrict);
 
-CREATEFUNCTION distance_sphere(geometry,geometry)
-        RETURNS FLOAT8 AS
-'
-DECLARE
-  proj1 integer;
-  proj2 integer;
-  projcount integer;
-BEGIN
-  SELECT INTO proj1 SRID($1);
-  SELECT INTO proj2 SRID($2);
-  IF proj1 != proj2 THEN
-    RAISE EXCEPTION ''Geometry projections must match.'';
-  END IF;
-  IF proj1 != -1 AND proj2 != -1 THEN
-    SELECT INTO projcount count(*) FROM spatial_ref_sys 
-    WHERE srid = proj1 AND proj4text like ''%longlat%'';
-    IF projcount = 0 THEN
-      RAISE EXCEPTION ''Geometries must be in a long/lat format.'';
-    END IF;
-    SELECT INTO projcount count(*) FROM spatial_ref_sys
-    WHERE srid = proj2 AND proj4text like ''%longlat%'';
-    IF projcount = 0 THEN
-      RAISE EXCEPTION ''Geometries must be in a long/lat format.'';
-    END IF;
-  END IF;
-  RETURN distance_sphere_latlon($1, $2);
-END;
-'
-LANGUAGE 'plpgsql' _IMMUTABLE_STRICT; -- WITH (isstrict,iscachable);
-
 -- Minimum distance. 2d only.
 CREATEFUNCTION distance(geometry,geometry)
 	RETURNS float8