From b89f4ad0bfd9c7863c183a3121442cdaeb3b131b Mon Sep 17 00:00:00 2001
From: David Blasby <dblasby@gmail.com>
Date: Mon, 29 Apr 2002 17:23:23 +0000
Subject: [PATCH] Added distance_ellipsiod(point,point,ellipsoid)  (called
 distance_spheroid in SQL)

git-svn-id: http://svn.osgeo.org/postgis/trunk@141 b70326c6-7e19-0410-871a-916f4a2858ee
---
 postgis.h      |  1 +
 postgis.sql.in |  4 ++++
 postgis_proj.c | 47 +++++++++++++++++++++++++++++++++++++++++++++++
 3 files changed, 52 insertions(+)

diff --git a/postgis.h b/postgis.h
index c458619a5..b75218082 100644
--- a/postgis.h
+++ b/postgis.h
@@ -444,6 +444,7 @@ Datum ellipsoid_out(PG_FUNCTION_ARGS);
 Datum ellipsoid_in(PG_FUNCTION_ARGS);
 Datum length_ellipsoid(PG_FUNCTION_ARGS);
 Datum length3d_ellipsoid(PG_FUNCTION_ARGS);
+Datum distance_ellipsoid(PG_FUNCTION_ARGS);
 
 Datum point_inside_circle(PG_FUNCTION_ARGS);
 Datum distance(PG_FUNCTION_ARGS);
diff --git a/postgis.sql.in b/postgis.sql.in
index b9b8de73d..4de81d7fd 100644
--- a/postgis.sql.in
+++ b/postgis.sql.in
@@ -574,6 +574,10 @@ CREATE FUNCTION length3d_spheroid(GEOMETRY,SPHEROID)
    AS '@MODULE_FILENAME@','length3d_ellipsoid'
 	     LANGUAGE 'c'  with (isstrict);
 
+CREATE FUNCTION distance_spheroid(GEOMETRY,GEOMETRY,SPHEROID)
+   RETURNS FLOAT8
+   AS '@MODULE_FILENAME@','distance_ellipsoid'
+	     LANGUAGE 'c'  with (isstrict);
 
 -------  generic operations
 
diff --git a/postgis_proj.c b/postgis_proj.c
index a93b9e8d8..2be6e5d08 100644
--- a/postgis_proj.c
+++ b/postgis_proj.c
@@ -335,8 +335,55 @@ Datum length3d_ellipsoid(PG_FUNCTION_ARGS)
 	PG_RETURN_FLOAT8(dist);
 }
 
+//distance (geometry,geometry, sphere)
+// -geometrys MUST be points
+PG_FUNCTION_INFO_V1(distance_ellipsoid);
+Datum distance_ellipsoid(PG_FUNCTION_ARGS)
+{
+		SPHEROID    		*sphere = (SPHEROID *) PG_GETARG_POINTER(2);	
+		GEOMETRY		      *geom1 = (GEOMETRY *)  PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+		GEOMETRY		      *geom2 = (GEOMETRY *)  PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
+
+		POINT3D	*pt1,*pt2;
+		int32				*offsets1;
+		int32				*offsets2;
+		char				*o;
+
+
+	if (geom1->SRID != geom2->SRID)
+	{
+		elog(ERROR,"optimistic_overlap:Operation on two GEOMETRIES with different SRIDs\n");
+		PG_RETURN_NULL();
+	}
+	
+	if (geom1->type != POINTTYPE)
+	{
+		elog(ERROR,"optimistic_overlap: first arg isnt a point\n");
+		PG_RETURN_NULL();
+	}
+	if (geom2->type != POINTTYPE)
+	{
+		elog(ERROR,"optimistic_overlap: second arg isnt a point\n");
+		PG_RETURN_NULL();
+	}
+
+	offsets1 = (int32 *) ( ((char *) &(geom1->objType[0] ))+ sizeof(int32) * geom1->nobjs ) ;
+	offsets2 = (int32 *) ( ((char *) &(geom2->objType[0] ))+ sizeof(int32) * geom2->nobjs ) ;
+	o = (char *) geom1 +offsets1[0] ;
+	pt1 = (POINT3D *) o;
 
+	o = (char *) geom2 +offsets2[0] ;
+	pt2 = (POINT3D *) o;
 
+	PG_RETURN_FLOAT8(distance_ellipse(pt1->y*M_PI/180.0 ,pt1->x*M_PI/180.0 ,
+							  pt2->y*M_PI/180.0 ,pt2->x*M_PI/180.0 , sphere) );
+
+//double	distance_ellipse(double lat1, double long1,
+//					double lat2, double long2,
+//					SPHEROID *sphere)
+
+
+}
 
 
 
-- 
2.40.0