--- /dev/null
+/*
+drop function geos_thing(geometry,geometry);
+CREATE FUNCTION relate(geometry,geometry)
+ RETURNS text
+ AS '/data1/Refractions/Projects/PostGIS/work_dave/postgis/libpostgis.so.0.7','relate_full'
+ LANGUAGE 'c' with (isstrict);
+
+CREATE FUNCTION relate(geometry,geometry,text)
+ RETURNS boolean
+ AS '/data1/Refractions/Projects/PostGIS/work_dave/postgis/libpostgis.so.0.7','relate_pattern'
+ LANGUAGE 'c' with (isstrict);
+
+CREATE FUNCTION disjoint(geometry,geometry)
+ RETURNS boolean
+ AS '/data1/Refractions/Projects/PostGIS/work_dave/postgis/libpostgis.so.0.7'
+ LANGUAGE 'c' with (isstrict);
+
+CREATE FUNCTION touches(geometry,geometry)
+ RETURNS boolean
+ AS '/data1/Refractions/Projects/PostGIS/work_dave/postgis/libpostgis.so.0.7'
+ LANGUAGE 'c' with (isstrict);
+
+CREATE FUNCTION intersects(geometry,geometry)
+ RETURNS boolean
+ AS '/data1/Refractions/Projects/PostGIS/work_dave/postgis/libpostgis.so.0.7'
+ LANGUAGE 'c' with (isstrict);
+
+CREATE FUNCTION crosses(geometry,geometry)
+ RETURNS boolean
+ AS '/data1/Refractions/Projects/PostGIS/work_dave/postgis/libpostgis.so.0.7'
+ LANGUAGE 'c' with (isstrict);
+
+CREATE FUNCTION within(geometry,geometry)
+ RETURNS boolean
+ AS '/data1/Refractions/Projects/PostGIS/work_dave/postgis/libpostgis.so.0.7'
+ LANGUAGE 'c' with (isstrict);
+
+CREATE FUNCTION contains(geometry,geometry)
+ RETURNS boolean
+ AS '/data1/Refractions/Projects/PostGIS/work_dave/postgis/libpostgis.so.0.7'
+ LANGUAGE 'c' with (isstrict);
+
+CREATE FUNCTION overlaps(geometry,geometry)
+ RETURNS boolean
+ AS '/data1/Refractions/Projects/PostGIS/work_dave/postgis/libpostgis.so.0.7'
+ LANGUAGE 'c' with (isstrict);
+
+ CREATE FUNCTION isvalid(geometry)
+ RETURNS boolean
+ AS '/data1/Refractions/Projects/PostGIS/work_dave/postgis/libpostgis.so.0.7'
+ LANGUAGE 'c' with (isstrict);
+*/
+
+//--------------------------------------------------------------------------
+//
+#ifdef USE_GEOS
+
+#include "postgres.h"
+
+
+#include <math.h>
+#include <float.h>
+#include <string.h>
+#include <stdio.h>
+#include <errno.h>
+
+#include "fmgr.h"
+
+#include "postgis.h"
+#include "utils/elog.h"
+
+#include "access/gist.h"
+#include "access/itup.h"
+#include "access/rtree.h"
+
+#include "utils/builtins.h"
+
+
+
+typedef struct Geometry Geometry;
+
+
+extern const char * createGEOSPoint(POINT3D *pt);
+extern void initGEOS(int maxalign);
+extern char *GEOSrelate(Geometry *g1, Geometry*g2);
+extern bool GEOSrelatePattern(Geometry *g1, Geometry*g2,char *pat);
+extern bool GEOSrelateDisjoint(Geometry *g1, Geometry*g2);
+extern bool GEOSrelateTouches(Geometry *g1, Geometry*g2);
+extern bool GEOSrelateIntersects(Geometry *g1, Geometry*g2);
+extern bool GEOSrelateCrosses(Geometry *g1, Geometry*g2);
+extern bool GEOSrelateWithin(Geometry *g1, Geometry*g2);
+extern bool GEOSrelateContains(Geometry *g1, Geometry*g2);
+extern bool GEOSrelateOverlaps(Geometry *g1, Geometry*g2);
+
+extern char *GEOSasText(Geometry *g1);
+
+
+extern void GEOSdeleteChar(char *a);
+extern void GEOSdeleteGeometry(Geometry *a);
+
+extern Geometry *PostGIS2GEOS_point(POINT3D *point,int SRID, bool is3d);
+extern Geometry *PostGIS2GEOS_linestring(LINE3D *line,int SRID, bool is3d);
+extern Geometry *PostGIS2GEOS_polygon(POLYGON3D *polygon,int SRID, bool is3d);
+extern Geometry *PostGIS2GEOS_multipolygon(POLYGON3D **polygons,int npolys, int SRID, bool is3d);
+extern Geometry *PostGIS2GEOS_multilinestring(LINE3D **lines,int nlines, int SRID, bool is3d);
+extern Geometry *PostGIS2GEOS_multipoint(POINT3D **points,int npoints, int SRID, bool is3d);
+extern Geometry *PostGIS2GEOS_box3d(BOX3D *box, int SRID);
+extern Geometry *PostGIS2GEOS_collection(Geometry **geoms, int ngeoms,int SRID, bool is3d);
+
+extern bool GEOSisvalid(Geometry *g1);
+
+
+Datum relate_full(PG_FUNCTION_ARGS);
+Datum relate_pattern(PG_FUNCTION_ARGS);
+Datum disjoint(PG_FUNCTION_ARGS);
+Datum touches(PG_FUNCTION_ARGS);
+Datum intersects(PG_FUNCTION_ARGS);
+Datum crosses(PG_FUNCTION_ARGS);
+Datum within(PG_FUNCTION_ARGS);
+Datum contains(PG_FUNCTION_ARGS);
+Datum overlaps(PG_FUNCTION_ARGS);
+
+
+Datum isvalid(PG_FUNCTION_ARGS);
+
+
+Geometry *POSTGIS2GEOS(GEOMETRY *g);
+void errorIfGeometryCollection(GEOMETRY *g1, GEOMETRY *g2);
+
+// return a GEOS Geometry from a POSTGIS GEOMETRY
+
+
+void errorIfGeometryCollection(GEOMETRY *g1, GEOMETRY *g2)
+{
+ if ( (g1->type == COLLECTIONTYPE) || (g2->type == COLLECTIONTYPE) )
+ elog(ERROR,"Relate Operation called with a GEOMETRYCOLLECTION type. This is unsupported");
+}
+
+PG_FUNCTION_INFO_V1(isvalid);
+Datum isvalid(PG_FUNCTION_ARGS)
+{
+ GEOMETRY *geom1 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+ bool result;
+ Geometry *g1;
+
+ initGEOS(MAXIMUM_ALIGNOF);
+
+ g1 = POSTGIS2GEOS(geom1 );
+ result = GEOSisvalid(g1);
+
+ GEOSdeleteGeometry(g1);
+
+ PG_RETURN_BOOL(result);
+}
+
+PG_FUNCTION_INFO_V1(overlaps);
+Datum overlaps(PG_FUNCTION_ARGS)
+{
+ GEOMETRY *geom1 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+ GEOMETRY *geom2 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
+
+
+ Geometry *g1,*g2;
+ bool result;
+
+ errorIfGeometryCollection(geom1,geom2);
+ initGEOS(MAXIMUM_ALIGNOF);
+
+ g1 = POSTGIS2GEOS(geom1 );
+ g2 = POSTGIS2GEOS(geom2 );
+
+
+ result = GEOSrelateOverlaps(g1,g2);
+
+ GEOSdeleteGeometry(g1);
+ GEOSdeleteGeometry(g2);
+
+ PG_RETURN_BOOL(result);
+}
+
+
+
+PG_FUNCTION_INFO_V1(contains);
+Datum contains(PG_FUNCTION_ARGS)
+{
+ GEOMETRY *geom1 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+ GEOMETRY *geom2 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
+
+
+ Geometry *g1,*g2;
+ bool result;
+
+ errorIfGeometryCollection(geom1,geom2);
+ initGEOS(MAXIMUM_ALIGNOF);
+
+ g1 = POSTGIS2GEOS(geom1 );
+ g2 = POSTGIS2GEOS(geom2 );
+
+
+ result = GEOSrelateContains(g1,g2);
+
+ GEOSdeleteGeometry(g1);
+ GEOSdeleteGeometry(g2);
+
+ PG_RETURN_BOOL(result);
+}
+
+
+PG_FUNCTION_INFO_V1(within);
+Datum within(PG_FUNCTION_ARGS)
+{
+ GEOMETRY *geom1 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+ GEOMETRY *geom2 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
+
+
+ Geometry *g1,*g2;
+ bool result;
+
+ errorIfGeometryCollection(geom1,geom2);
+ initGEOS(MAXIMUM_ALIGNOF);
+
+ g1 = POSTGIS2GEOS(geom1 );
+ g2 = POSTGIS2GEOS(geom2 );
+
+
+ result = GEOSrelateWithin(g1,g2);
+
+ GEOSdeleteGeometry(g1);
+ GEOSdeleteGeometry(g2);
+
+ PG_RETURN_BOOL(result);
+}
+
+
+
+PG_FUNCTION_INFO_V1(crosses);
+Datum crosses(PG_FUNCTION_ARGS)
+{
+ GEOMETRY *geom1 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+ GEOMETRY *geom2 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
+
+
+ Geometry *g1,*g2;
+ bool result;
+
+ errorIfGeometryCollection(geom1,geom2);
+ initGEOS(MAXIMUM_ALIGNOF);
+
+ g1 = POSTGIS2GEOS(geom1 );
+ g2 = POSTGIS2GEOS(geom2 );
+
+
+ result = GEOSrelateCrosses(g1,g2);
+
+ GEOSdeleteGeometry(g1);
+ GEOSdeleteGeometry(g2);
+
+ PG_RETURN_BOOL(result);
+}
+
+
+
+PG_FUNCTION_INFO_V1(intersects);
+Datum intersects(PG_FUNCTION_ARGS)
+{
+ GEOMETRY *geom1 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+ GEOMETRY *geom2 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
+
+
+ Geometry *g1,*g2;
+ bool result;
+
+ errorIfGeometryCollection(geom1,geom2);
+ initGEOS(MAXIMUM_ALIGNOF);
+
+ g1 = POSTGIS2GEOS(geom1 );
+ g2 = POSTGIS2GEOS(geom2 );
+
+
+ result = GEOSrelateIntersects(g1,g2);
+
+ GEOSdeleteGeometry(g1);
+ GEOSdeleteGeometry(g2);
+
+ PG_RETURN_BOOL(result);
+}
+
+
+PG_FUNCTION_INFO_V1(touches);
+Datum touches(PG_FUNCTION_ARGS)
+{
+ GEOMETRY *geom1 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+ GEOMETRY *geom2 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
+
+
+ Geometry *g1,*g2;
+ bool result;
+
+ errorIfGeometryCollection(geom1,geom2);
+ initGEOS(MAXIMUM_ALIGNOF);
+
+ g1 = POSTGIS2GEOS(geom1 );
+ g2 = POSTGIS2GEOS(geom2 );
+
+
+ result = GEOSrelateTouches(g1,g2);
+
+ GEOSdeleteGeometry(g1);
+ GEOSdeleteGeometry(g2);
+
+ PG_RETURN_BOOL(result);
+}
+
+
+PG_FUNCTION_INFO_V1(disjoint);
+Datum disjoint(PG_FUNCTION_ARGS)
+{
+ GEOMETRY *geom1 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+ GEOMETRY *geom2 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
+
+
+ Geometry *g1,*g2;
+ bool result;
+
+ errorIfGeometryCollection(geom1,geom2);
+ initGEOS(MAXIMUM_ALIGNOF);
+
+ g1 = POSTGIS2GEOS(geom1 );
+ g2 = POSTGIS2GEOS(geom2 );
+
+
+ result = GEOSrelateDisjoint(g1,g2);
+
+ GEOSdeleteGeometry(g1);
+ GEOSdeleteGeometry(g2);
+
+ PG_RETURN_BOOL(result);
+}
+
+
+PG_FUNCTION_INFO_V1(relate_pattern);
+Datum relate_pattern(PG_FUNCTION_ARGS)
+{
+ GEOMETRY *geom1 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+ GEOMETRY *geom2 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
+ char *patt;
+ bool result;
+
+ Geometry *g1,*g2;
+
+
+ errorIfGeometryCollection(geom1,geom2);
+ initGEOS(MAXIMUM_ALIGNOF);
+
+ g1 = POSTGIS2GEOS(geom1 );
+ g2 = POSTGIS2GEOS(geom2 );
+
+
+ patt = DatumGetCString(DirectFunctionCall1(textout,
+ PointerGetDatum(PG_GETARG_DATUM(2))));
+
+ result = GEOSrelatePattern(g1,g2,patt);
+
+
+ GEOSdeleteGeometry(g1);
+ GEOSdeleteGeometry(g2);
+
+ pfree(patt);
+ PG_RETURN_BOOL(result);
+
+
+}
+
+
+
+PG_FUNCTION_INFO_V1(relate_full);
+Datum relate_full(PG_FUNCTION_ARGS)
+{
+ GEOMETRY *geom1 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
+ GEOMETRY *geom2 = (GEOMETRY *) PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
+
+ Geometry *g1,*g2;
+ char *relate_str;
+ int len;
+ char *result;
+
+ errorIfGeometryCollection(geom1,geom2);
+
+//elog(NOTICE,"in relate_full: \n");
+ initGEOS(MAXIMUM_ALIGNOF);
+
+//elog(NOTICE,"making GOES geometries: \n");
+ g1 = POSTGIS2GEOS(geom1 );
+ g2 = POSTGIS2GEOS(geom2 );
+// elog(NOTICE,"geometries made, here's what they are: \n");
+//g1 = createGEOSFromText("POINT(1 1)");
+//g2 = createGEOSFromText("POINT(1 1)");
+//elog(NOTICE,GEOSasText(g1));
+//elog(NOTICE,GEOSasText(g2));
+
+//elog(NOTICE,"calling relate: \n");
+ relate_str = GEOSrelate(g1, g2);
+
+//elog(NOTICE,"relate finished \n");
+
+
+ len = strlen(relate_str) + 4;
+
+ result= palloc(len);
+ *((int *) result) = len;
+
+ memcpy(result +4, relate_str, len-4);
+
+ free(relate_str);
+// GEOSdeleteGeometry(g1);
+// GEOSdeleteGeometry(g2);
+
+ PG_RETURN_POINTER(result);
+}
+
+
+//BBOXONLYTYPE -> returns as a 2d polygon
+Geometry *POSTGIS2GEOS(GEOMETRY *g)
+{
+ POINT3D *pt;
+ LINE3D *line;
+ POLYGON3D *poly;
+ POLYGON3D **polys;
+ LINE3D **lines;
+ POINT3D **points;
+ Geometry **geoms;
+ Geometry *geos;
+ char *obj;
+ int obj_type;
+ int t;
+
+ int32 *offsets1 = (int32 *) ( ((char *) &(g->objType[0] ))+ sizeof(int32) * g->nobjs ) ;
+
+ switch(g->type)
+ {
+ case POINTTYPE:
+ pt = (POINT3D*) ((char *) g +offsets1[0]) ;
+ return PostGIS2GEOS_point(pt,g->SRID,g->is3d);
+ break;
+ case LINETYPE:
+ line = (LINE3D*) ((char *) g +offsets1[0]) ;
+ return PostGIS2GEOS_linestring(line,g->SRID,g->is3d);
+ break;
+ case POLYGONTYPE:
+ poly = (POLYGON3D*) ((char *) g +offsets1[0]) ;
+ return PostGIS2GEOS_polygon(poly,g->SRID,g->is3d);
+ break;
+ case MULTIPOLYGONTYPE:
+ //make an array of POLYGON3Ds
+ polys = (POLYGON3D**) palloc(sizeof (POLYGON3D*) * g->nobjs);
+ for (t=0;t<g->nobjs;t++)
+ {
+ polys[t] = (POLYGON3D*) ((char *) g +offsets1[t]) ;
+ }
+ geos= PostGIS2GEOS_multipolygon(polys, g->nobjs, g->SRID,g->is3d);
+ pfree(polys);
+ return geos;
+ break;
+ case MULTILINETYPE:
+ //make an array of POLYGON3Ds
+ lines = (LINE3D**) palloc(sizeof (LINE3D*) * g->nobjs);
+ for (t=0;t<g->nobjs;t++)
+ {
+ lines[t] = (LINE3D*) ((char *) g +offsets1[t]) ;
+ }
+ geos= PostGIS2GEOS_multilinestring(lines, g->nobjs, g->SRID,g->is3d);
+ pfree(lines);
+ return geos;
+ break;
+ case MULTIPOINTTYPE:
+ //make an array of POINT3Ds
+ points = (POINT3D**) palloc(sizeof (POINT3D*) * g->nobjs);
+ for (t=0;t<g->nobjs;t++)
+ {
+ points[t] = (POINT3D*) ((char *) g +offsets1[t]) ;
+ }
+ geos= PostGIS2GEOS_multipoint(points, g->nobjs,g->SRID,g->is3d);
+ pfree(points);
+ return geos;
+ break;
+ case BBOXONLYTYPE:
+ return PostGIS2GEOS_box3d(&g->bvol, g->SRID);
+ break;
+ case COLLECTIONTYPE:
+ //make an array of GEOS Geometrys
+ geoms = (Geometry**) palloc(sizeof (Geometry*) * g->nobjs);
+ for (t=0;t<g->nobjs;t++)
+ {
+ obj = ((char *) g +offsets1[t]);
+ obj_type = g->objType[t];
+ switch (obj_type)
+ {
+ case POINTTYPE:
+ pt = (POINT3D*) obj ;
+ geoms[t] = PostGIS2GEOS_point(pt,g->SRID,g->is3d);
+ break;
+ case LINETYPE:
+ line = (LINE3D*) obj ;
+ geoms[t] = PostGIS2GEOS_linestring(line,g->SRID,g->is3d);
+ break;
+ case POLYGONTYPE:
+ poly = (POLYGON3D*) obj ;
+ geoms[t] = PostGIS2GEOS_polygon(poly,g->SRID,g->is3d);
+ break;
+ }
+ }
+ geos= PostGIS2GEOS_collection(geoms,g->nobjs,g->SRID,g->is3d);
+ pfree(geoms);
+ return geos;
+ break;
+
+ }
+ return NULL;
+}
+
+
+
+
+//----------------------------------------------------------------------------
+// NULL implementation here
+// ---------------------------------------------------------------------------
+#else
+
+
+#include "postgres.h"
+
+
+#include <math.h>
+#include <float.h>
+#include <string.h>
+#include <stdio.h>
+#include <errno.h>
+
+#include "fmgr.h"
+
+
+Datum relate_full(PG_FUNCTION_ARGS);
+Datum relate_pattern(PG_FUNCTION_ARGS);
+Datum disjoint(PG_FUNCTION_ARGS);
+Datum touches(PG_FUNCTION_ARGS);
+Datum intersects(PG_FUNCTION_ARGS);
+Datum crosses(PG_FUNCTION_ARGS);
+Datum within(PG_FUNCTION_ARGS);
+Datum contains(PG_FUNCTION_ARGS);
+Datum overlaps(PG_FUNCTION_ARGS);
+Datum isvalid(PG_FUNCTION_ARGS);
+
+PG_FUNCTION_INFO_V1(relate_full);
+Datum relate_full(PG_FUNCTION_ARGS)
+{
+ elog(ERROR,"relate_full:: operation not implemented - compile PostGIS with GEOS support");
+ PG_RETURN_NULL(); // never get here
+}
+PG_FUNCTION_INFO_V1(relate_pattern);
+Datum relate_pattern(PG_FUNCTION_ARGS)
+{
+ elog(ERROR,"relate_pattern:: operation not implemented - compile PostGIS with GEOS support");
+ PG_RETURN_NULL(); // never get here
+}
+PG_FUNCTION_INFO_V1(disjoint);
+Datum disjoint(PG_FUNCTION_ARGS)
+{
+ elog(ERROR,"disjoint:: operation not implemented - compile PostGIS with GEOS support");
+ PG_RETURN_NULL(); // never get here
+}
+PG_FUNCTION_INFO_V1(intersects);
+Datum intersects(PG_FUNCTION_ARGS)
+{
+ elog(ERROR,"intersects:: operation not implemented - compile PostGIS with GEOS support");
+ PG_RETURN_NULL(); // never get here
+}
+PG_FUNCTION_INFO_V1(touches);
+Datum touches(PG_FUNCTION_ARGS)
+{
+ elog(ERROR,"touches:: operation not implemented - compile PostGIS with GEOS support");
+ PG_RETURN_NULL(); // never get here
+}
+PG_FUNCTION_INFO_V1(crosses);
+Datum crosses(PG_FUNCTION_ARGS)
+{
+ elog(ERROR,"crosses:: operation not implemented - compile PostGIS with GEOS support");
+ PG_RETURN_NULL(); // never get here
+}
+PG_FUNCTION_INFO_V1(within);
+Datum within(PG_FUNCTION_ARGS)
+{
+ elog(ERROR,"within:: operation not implemented - compile PostGIS with GEOS support");
+ PG_RETURN_NULL(); // never get here
+}
+PG_FUNCTION_INFO_V1(contains);
+Datum contains(PG_FUNCTION_ARGS)
+{
+ elog(ERROR,"contains:: operation not implemented - compile PostGIS with GEOS support");
+ PG_RETURN_NULL(); // never get here
+}
+PG_FUNCTION_INFO_V1(overlaps);
+Datum overlaps(PG_FUNCTION_ARGS)
+{
+ elog(ERROR,"overlaps:: operation not implemented - compile PostGIS with GEOS support");
+ PG_RETURN_NULL(); // never get here
+}
+PG_FUNCTION_INFO_V1(isvalid);
+Datum isvalid(PG_FUNCTION_ARGS)
+{
+ elog(ERROR,"isvalid:: operation not implemented - compile PostGIS with GEOS support");
+ PG_RETURN_NULL(); // never get here
+}
+
+#endif
+
--- /dev/null
+// g++ postgis_GEOSwrapper.cpp -c -I/usr/local/include -I/usr/local/include/geos -I/usr/local/src/postgresql-7.2.3//src/include
+
+
+#include <stdio.h>
+
+#include <string>
+#include <iostream>
+#include <fstream>
+
+//#include "geom.h"
+//#include "util.h"
+//#include "graph.h"
+#include "io.h"
+//#include "opRelate.h"
+
+
+
+//WARNING THIS *MUST* BE SET CORRECTLY.
+int MAXIMUM_ALIGNOF = -999999; // to be set during initialization - this will be either 4 (intel) or 8 (sparc)
+
+//for getting things to align properly double are on 8byte align on solaris machines, and 4bytes on intel
+
+#define TYPEALIGN(ALIGNVAL,LEN) (((long)(LEN) + (ALIGNVAL-1)) & ~(ALIGNVAL-1))
+#define MAXALIGN(LEN) TYPEALIGN(MAXIMUM_ALIGNOF, (LEN))
+
+typedef int int32;
+
+typedef struct
+{
+ double x,y,z; //for lat/long x=long, y=lat
+} POINT3D;
+
+typedef struct
+{
+ POINT3D LLB,URT; /* corner POINT3Ds on long diagonal */
+} BOX3D;
+
+typedef struct
+{
+ int32 npoints; // how many points in the line
+ int32 junk; // double-word alignment
+ POINT3D points[1]; // array of actual points
+} LINE3D;
+
+
+typedef struct
+{
+ int32 nrings; // how many rings in this polygon
+ int32 npoints[1]; //how many points in each ring
+ /* could be 4 byes of filler here to make sure points[] is
+ double-word aligned*/
+ POINT3D points[1]; // array of actual points
+} POLYGON3D;
+
+typedef struct
+{
+ int32 size; // postgres variable-length type requirement
+ int32 SRID; // spatial reference system id
+ double offsetX; // for precision grid (future improvement)
+ double offsetY; // for precision grid (future improvement)
+ double scale; // for precision grid (future improvement)
+ int32 type; // this type of geometry
+ bool is3d; // true if the points are 3d (only for output)
+ BOX3D bvol; // bounding volume of all the geo objects
+ int32 nobjs; // how many sub-objects in this object
+ int32 objType[1]; // type of object
+ int32 objOffset[1];// offset (in bytes) into this structure where
+ // the object is located
+ char objData[1]; // store for actual objects
+
+} GEOMETRY;
+
+
+//###########################################################
+
+extern "C" char *GEOSrelate(Geometry *g1, Geometry*g2);
+extern "C" void initGEOS(int maxalign);
+
+
+extern "C" void GEOSdeleteChar(char *a);
+extern "C" void GEOSdeleteGeometry(Geometry *a);
+extern "C" bool GEOSrelatePattern(Geometry *g1, Geometry*g2,char *pat);
+extern "C" bool GEOSrelateDisjoint(Geometry *g1, Geometry*g2);
+extern "C" bool GEOSrelateTouches(Geometry *g1, Geometry*g2);
+extern "C" bool GEOSrelateIntersects(Geometry *g1, Geometry*g2);
+extern "C" bool GEOSrelateCrosses(Geometry *g1, Geometry*g2);
+extern "C" bool GEOSrelateWithin(Geometry *g1, Geometry*g2);
+extern "C" bool GEOSrelateContains(Geometry *g1, Geometry*g2);
+extern "C" bool GEOSrelateOverlaps(Geometry *g1, Geometry*g2);
+
+
+extern "C" Geometry *PostGIS2GEOS_point(POINT3D *point,int SRID, bool is3d);
+extern "C" Geometry *PostGIS2GEOS_linestring(LINE3D *line,int SRID, bool is3d);
+extern "C" Geometry *PostGIS2GEOS_polygon(POLYGON3D *polygon,int SRID, bool is3d);
+extern "C" Geometry *PostGIS2GEOS_multipolygon(POLYGON3D **polygons,int npolys, int SRID, bool is3d);
+extern "C" Geometry *PostGIS2GEOS_multilinestring(LINE3D **lines,int nlines, int SRID, bool is3d);
+extern "C" Geometry *PostGIS2GEOS_multipoint(POINT3D **points,int npoints, int SRID, bool is3d);
+
+extern "C" Geometry *PostGIS2GEOS_box3d(BOX3D *box, int SRID);
+extern "C" Geometry *PostGIS2GEOS_collection(Geometry **geoms, int ngeoms,int SRID, bool is3d);
+
+extern "C" bool GEOSisvalid(Geometry *g1);
+
+
+//###########################################################
+
+GeometryFactory *geomFactory = NULL;
+
+void initGEOS (int maxalign)
+{
+ if (geomFactory == NULL)
+ {
+ geomFactory = new GeometryFactory( new PrecisionModel(), -1); // NOTE: SRID will have to be changed after geometry creation
+ MAXIMUM_ALIGNOF = maxalign;
+ }
+}
+
+ //note: you lose the 3d from this!
+Geometry *PostGIS2GEOS_box3d(BOX3D *box, int SRID)
+{
+ Geometry *g;
+ Envelope *envelope = new Envelope(box->LLB.x,box->URT.x ,box->LLB.y,box->URT.y);
+ g = geomFactory->toGeometry(envelope,geomFactory->getPrecisionModel(), SRID);
+ delete envelope;
+ return g;
+}
+
+Geometry *PostGIS2GEOS_collection(Geometry **geoms, int ngeoms,int SRID, bool is3d)
+{
+ Geometry *g;
+ int t;
+ vector<Geometry *> *subGeos=new vector<Geometry *>;
+
+ for (t =0; t< ngeoms; t++)
+ {
+ subGeos->push_back(geoms[t]);
+ }
+ g = geomFactory->buildGeometry(subGeos);
+ g->setSRID(SRID);
+ return g;
+
+}
+
+Geometry *PostGIS2GEOS_point(POINT3D *point,int SRID, bool is3d)
+{
+ Coordinate *c;
+
+ if (is3d)
+ c = new Coordinate(point->x, point->y);
+ else
+ c = new Coordinate(point->x, point->y, point->z);
+ Geometry *g = geomFactory->createPoint(*c);
+ g->setSRID(SRID);
+ return g;
+}
+
+
+Geometry *PostGIS2GEOS_linestring(LINE3D *line,int SRID, bool is3d)
+{
+ int t;
+ Coordinate c;
+
+ //build coordinatelist & pre-allocate space
+ BasicCoordinateList *coords = new BasicCoordinateList(line->npoints);
+ if (is3d)
+ {
+ for (t=0;t<line->npoints;t++)
+ {
+ c.x = line->points[t].x;
+ c.y = line->points[t].y;
+ c.z = line->points[t].z;
+ coords->setAt( c ,t);
+ }
+
+ }
+ else //make 2d points
+ {
+ for (t=0;t<line->npoints;t++)
+ {
+ c.x = line->points[t].x;
+ c.y = line->points[t].y;
+ c.z = DoubleNotANumber;
+ coords->setAt( c ,t);
+ }
+
+ }
+ Geometry *g = geomFactory->createLineString(coords);
+ g->setSRID(SRID);
+ return g;
+}
+
+ //polygons is an array of pointers to polygons
+Geometry *PostGIS2GEOS_multipolygon(POLYGON3D **polygons,int npolys, int SRID, bool is3d)
+{
+ int t;
+ vector<Geometry *> *subPolys=new vector<Geometry *>;
+ Geometry *g;
+
+ for (t =0; t< npolys; t++)
+ {
+ subPolys->push_back(PostGIS2GEOS_polygon(polygons[t], SRID,is3d ));
+ }
+ g = geomFactory->createMultiPolygon(subPolys);
+ g->setSRID(SRID);
+ return g;
+}
+
+ //lines is an array of pointers to line3d
+Geometry *PostGIS2GEOS_multilinestring(LINE3D **lines,int nlines, int SRID, bool is3d)
+{
+ int t;
+ vector<Geometry *> *subLines =new vector<Geometry *>;
+ Geometry *g;
+
+ for (t =0; t< nlines; t++)
+ {
+ subLines->push_back(PostGIS2GEOS_linestring(lines[t], SRID,is3d ));
+ }
+ g = geomFactory->createMultiLineString(subLines);
+ g->setSRID(SRID);
+ return g;
+}
+
+Geometry *PostGIS2GEOS_multipoint(POINT3D **points,int npoints, int SRID, bool is3d)
+{
+ int t;
+ vector<Geometry *> *subPoints =new vector<Geometry *>;
+ Geometry *g;
+
+ for (t =0; t< npoints; t++)
+ {
+ subPoints->push_back(PostGIS2GEOS_point(points[t], SRID,is3d ));
+ }
+ g = geomFactory->createMultiPoint(subPoints);
+ g->setSRID(SRID);
+ return g;
+
+}
+
+
+Geometry *PostGIS2GEOS_polygon(POLYGON3D *polygon,int SRID, bool is3d)
+{
+ POINT3D *pts;
+ Coordinate c;
+ int ring,t;
+ Geometry *g;
+ LinearRing *outerRing;
+ LinearRing *innerRing;
+ BasicCoordinateList *cl;
+ int pointOffset =0; // first point that we're looking at. a POLYGON3D has all its points smooshed together
+ vector<Geometry *> *innerRings=new vector<Geometry *>;
+
+
+ pts = (POINT3D *) ( (char *)&(polygon->npoints[polygon->nrings] ) );
+ pts = (POINT3D *) MAXALIGN(pts);
+
+ // make outerRing
+ cl = new BasicCoordinateList(polygon->npoints[0]);
+ if (is3d)
+ {
+ for(t=0;t<polygon->npoints[0];t++)
+ {
+ c.x = pts[t].x;
+ c.y = pts[t].y;
+ c.z = pts[t].z;
+ cl->setAt( c ,t);
+ }
+ }
+ else
+ {
+ for(t=0;t<polygon->npoints[0];t++)
+ {
+ c.x = pts[t].x;
+ c.y = pts[t].y;
+ c.z = DoubleNotANumber;
+ cl->setAt( c ,t);
+ }
+ }
+ outerRing = geomFactory->createLinearRing(cl);
+ outerRing->setSRID(SRID);
+ pointOffset = polygon->npoints[0];
+
+ for(ring =1; ring< polygon->nrings; ring++)
+ {
+ cl = new BasicCoordinateList(polygon->npoints[ring]);
+ if (is3d)
+ {
+ for(t=0;t<polygon->npoints[ring];t++)
+ {
+ c.x = pts[t+pointOffset].x;
+ c.y = pts[t+pointOffset].y;
+ c.z = pts[t+pointOffset].z;
+ cl->setAt( c ,t);
+ }
+ }
+ else
+ {
+ for(t=0;t<polygon->npoints[ring];t++)
+ {
+ c.x = pts[t+pointOffset].x;
+ c.y = pts[t+pointOffset].y;
+ c.z = DoubleNotANumber;
+ cl->setAt( c ,t);
+ }
+ }
+ innerRing = geomFactory->createLinearRing(cl);
+ innerRing->setSRID(SRID);
+ innerRings->push_back(innerRing);
+ pointOffset += polygon->npoints[ring];
+ }
+
+ g = geomFactory->createPolygon(outerRing, innerRings);
+ g->setSRID(SRID);
+ return g;
+}
+
+
+
+bool GEOSrelateDisjoint(Geometry *g1, Geometry*g2)
+{
+ return g1->disjoint(g2);
+}
+
+bool GEOSrelateTouches(Geometry *g1, Geometry*g2)
+{
+ return g1->touches(g2);
+}
+
+bool GEOSrelateIntersects(Geometry *g1, Geometry*g2)
+{
+ return g1->intersects(g2);
+}
+
+bool GEOSrelateCrosses(Geometry *g1, Geometry*g2)
+{
+ return g1->crosses(g2);
+}
+
+bool GEOSrelateWithin(Geometry *g1, Geometry*g2)
+{
+ return g1->within(g2);
+}
+
+bool GEOSrelateContains(Geometry *g1, Geometry*g2)
+{
+ return g1->contains(g2);
+}
+
+bool GEOSrelateOverlaps(Geometry *g1, Geometry*g2)
+{
+ return g1->overlaps(g2);
+}
+
+//BUG:: this leaks memory, but delete kills the PrecisionModel for ALL the geometries
+void GEOSdeleteGeometry(Geometry *a)
+{
+// delete a;
+}
+
+void GEOSdeleteChar(char *a)
+{
+ free(a);
+}
+
+
+bool GEOSrelatePattern(Geometry *g1, Geometry*g2,char *pat)
+{
+ string s = pat;
+ return g1->relate(g2,pat);
+}
+
+
+char *GEOSrelate(Geometry *g1, Geometry*g2)
+{
+ IntersectionMatrix *im = g1->relate(g2);
+ string s;
+ char *result;
+
+ s= im->toString();
+ result = (char*) malloc( s.length() + 1);
+ strcpy(result, s.c_str() );
+
+
+ return result;
+}
+
+bool GEOSisvalid(Geometry *g1)
+{
+ return g1->isValid();
+}
+
+
+
+
+char *GEOSasText(Geometry *g1)
+{
+ string s = g1->toString();
+
+ return (char *) s.c_str() ;
+}
+