/**********************************************************************
*
* PostGIS - Spatial Types for PostgreSQL
* http://postgis.net
*
* PostGIS is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 2 of the License, or
* (at your option) any later version.
*
* PostGIS is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with PostGIS. If not, see .
*
**********************************************************************
*
* Copyright 2012-2013 Oslandia
*
**********************************************************************/
#include "lwgeom_sfcgal.h"
#include
static int SFCGAL_type_to_lwgeom_type(sfcgal_geometry_type_t type);
static POINTARRAY *ptarray_from_SFCGAL(const sfcgal_geometry_t *geom, int force3D);
static sfcgal_geometry_t *ptarray_to_SFCGAL(const POINTARRAY *pa, int type);
/* Return SFCGAL version string */
const char *
lwgeom_sfcgal_version()
{
const char *version = sfcgal_version();
return version;
}
#define MAX_LENGTH_SFCGAL_FULL_VERSION 50
/* Return SFCGAL full version string */
const char *
lwgeom_sfcgal_full_version()
{
#if POSTGIS_SFCGAL_VERSION >= 10400
const char *version = sfcgal_full_version();
#else
char *version = (char*)malloc(MAX_LENGTH_SFCGAL_FULL_VERSION);
snprintf(version, MAX_LENGTH_SFCGAL_FULL_VERSION, "SFCGAL=\"%s\" CGAL=\"Unknown\" Boost=\"Unknown\"", sfcgal_version());
#endif
return version;
}
/*
* Mapping between SFCGAL and PostGIS types
*
* Throw an error if type is unsupported
*/
static int
SFCGAL_type_to_lwgeom_type(sfcgal_geometry_type_t type)
{
switch (type)
{
case SFCGAL_TYPE_POINT:
return POINTTYPE;
case SFCGAL_TYPE_LINESTRING:
return LINETYPE;
case SFCGAL_TYPE_POLYGON:
return POLYGONTYPE;
case SFCGAL_TYPE_MULTIPOINT:
return MULTIPOINTTYPE;
case SFCGAL_TYPE_MULTILINESTRING:
return MULTILINETYPE;
case SFCGAL_TYPE_MULTIPOLYGON:
return MULTIPOLYGONTYPE;
case SFCGAL_TYPE_MULTISOLID:
return COLLECTIONTYPE; /* Nota: PolyhedralSurface closed inside
aim is to use true solid type as soon
as available in OGC SFS */
case SFCGAL_TYPE_GEOMETRYCOLLECTION:
return COLLECTIONTYPE;
#if 0
case SFCGAL_TYPE_CIRCULARSTRING:
return CIRCSTRINGTYPE;
case SFCGAL_TYPE_COMPOUNDCURVE:
return COMPOUNDTYPE;
case SFCGAL_TYPE_CURVEPOLYGON:
return CURVEPOLYTYPE;
case SFCGAL_TYPE_MULTICURVE:
return MULTICURVETYPE;
case SFCGAL_TYPE_MULTISURFACE:
return MULTISURFACETYPE;
#endif
case SFCGAL_TYPE_POLYHEDRALSURFACE:
return POLYHEDRALSURFACETYPE;
case SFCGAL_TYPE_TRIANGULATEDSURFACE:
return TINTYPE;
case SFCGAL_TYPE_TRIANGLE:
return TRIANGLETYPE;
default:
lwerror("SFCGAL_type_to_lwgeom_type: Unknown Type");
return 0;
}
}
/*
* Return a PostGIS pointarray from a simple SFCGAL geometry:
* POINT, LINESTRING or TRIANGLE
*
* Throw an error on others types
*/
static POINTARRAY *
ptarray_from_SFCGAL(const sfcgal_geometry_t *geom, int want3d)
{
POINT4D point;
uint32_t i, npoints;
POINTARRAY *pa = NULL;
assert(geom);
switch (sfcgal_geometry_type_id(geom))
{
case SFCGAL_TYPE_POINT:
{
pa = ptarray_construct(want3d, 0, 1);
point.x = sfcgal_point_x(geom);
point.y = sfcgal_point_y(geom);
if (sfcgal_geometry_is_3d(geom))
point.z = sfcgal_point_z(geom);
else if (want3d)
point.z = 0.0;
ptarray_set_point4d(pa, 0, &point);
}
break;
case SFCGAL_TYPE_LINESTRING:
{
npoints = sfcgal_linestring_num_points(geom);
pa = ptarray_construct(want3d, 0, npoints);
for (i = 0; i < npoints; i++)
{
const sfcgal_geometry_t *pt = sfcgal_linestring_point_n(geom, i);
point.x = sfcgal_point_x(pt);
point.y = sfcgal_point_y(pt);
if (sfcgal_geometry_is_3d(geom))
point.z = sfcgal_point_z(pt);
else if (want3d)
point.z = 0.0;
ptarray_set_point4d(pa, i, &point);
}
}
break;
case SFCGAL_TYPE_TRIANGLE:
{
pa = ptarray_construct(want3d, 0, 4);
for (i = 0; i < 4; i++)
{
const sfcgal_geometry_t *pt = sfcgal_triangle_vertex(geom, (i % 3));
point.x = sfcgal_point_x(pt);
point.y = sfcgal_point_y(pt);
if (sfcgal_geometry_is_3d(geom))
point.z = sfcgal_point_z(pt);
else if (want3d)
point.z = 0.0;
ptarray_set_point4d(pa, i, &point);
}
}
break;
/* Other types should not be called directly ... */
default:
lwerror("ptarray_from_SFCGAL: Unknown Type");
break;
}
return pa;
}
/*
* Convert a PostGIS pointarray to SFCGAL structure
*
* Used for simple LWGEOM geometry POINT, LINESTRING, TRIANGLE
* and POLYGON rings
*/
static sfcgal_geometry_t *
ptarray_to_SFCGAL(const POINTARRAY *pa, int type)
{
POINT3DZ point;
int is_3d;
uint32_t i;
assert(pa);
is_3d = FLAGS_GET_Z(pa->flags) != 0;
switch (type)
{
case POINTTYPE:
{
getPoint3dz_p(pa, 0, &point);
if (is_3d)
return sfcgal_point_create_from_xyz(point.x, point.y, point.z);
else
return sfcgal_point_create_from_xy(point.x, point.y);
}
break;
case LINETYPE:
{
sfcgal_geometry_t *line = sfcgal_linestring_create();
for (i = 0; i < pa->npoints; i++)
{
getPoint3dz_p(pa, i, &point);
if (is_3d)
{
sfcgal_linestring_add_point(line,
sfcgal_point_create_from_xyz(point.x, point.y, point.z));
}
else
{
sfcgal_linestring_add_point(line, sfcgal_point_create_from_xy(point.x, point.y));
}
}
return line;
}
break;
case TRIANGLETYPE:
{
sfcgal_geometry_t *triangle = sfcgal_triangle_create();
getPoint3dz_p(pa, 0, &point);
if (is_3d)
sfcgal_triangle_set_vertex_from_xyz(triangle, 0, point.x, point.y, point.z);
else
sfcgal_triangle_set_vertex_from_xy(triangle, 0, point.x, point.y);
getPoint3dz_p(pa, 1, &point);
if (is_3d)
sfcgal_triangle_set_vertex_from_xyz(triangle, 1, point.x, point.y, point.z);
else
sfcgal_triangle_set_vertex_from_xy(triangle, 1, point.x, point.y);
getPoint3dz_p(pa, 2, &point);
if (is_3d)
sfcgal_triangle_set_vertex_from_xyz(triangle, 2, point.x, point.y, point.z);
else
sfcgal_triangle_set_vertex_from_xy(triangle, 2, point.x, point.y);
return triangle;
}
break;
/* Other SFCGAL types should not be called directly ... */
default:
lwerror("ptarray_from_SFCGAL: Unknown Type");
return NULL;
}
}
/*
* Convert a SFCGAL structure to PostGIS LWGEOM
*
* Throws an error on unsupported type
*/
LWGEOM *
SFCGAL2LWGEOM(const sfcgal_geometry_t *geom, int force3D, int32_t srid)
{
uint32_t ngeoms, nshells, nsolids;
uint32_t i, j, k;
int want3d;
assert(geom);
want3d = force3D || sfcgal_geometry_is_3d(geom);
switch (sfcgal_geometry_type_id(geom))
{
case SFCGAL_TYPE_POINT:
{
if (sfcgal_geometry_is_empty(geom))
return (LWGEOM *)lwpoint_construct_empty(srid, want3d, 0);
POINTARRAY *pa = ptarray_from_SFCGAL(geom, want3d);
return (LWGEOM *)lwpoint_construct(srid, NULL, pa);
}
case SFCGAL_TYPE_LINESTRING:
{
if (sfcgal_geometry_is_empty(geom))
return (LWGEOM *)lwline_construct_empty(srid, want3d, 0);
POINTARRAY *pa = ptarray_from_SFCGAL(geom, want3d);
return (LWGEOM *)lwline_construct(srid, NULL, pa);
}
case SFCGAL_TYPE_TRIANGLE:
{
if (sfcgal_geometry_is_empty(geom))
return (LWGEOM *)lwtriangle_construct_empty(srid, want3d, 0);
POINTARRAY *pa = ptarray_from_SFCGAL(geom, want3d);
return (LWGEOM *)lwtriangle_construct(srid, NULL, pa);
}
case SFCGAL_TYPE_POLYGON:
{
if (sfcgal_geometry_is_empty(geom))
return (LWGEOM *)lwpoly_construct_empty(srid, want3d, 0);
uint32_t nrings = sfcgal_polygon_num_interior_rings(geom) + 1;
POINTARRAY **pa = (POINTARRAY **)lwalloc(sizeof(POINTARRAY *) * nrings);
pa[0] = ptarray_from_SFCGAL(sfcgal_polygon_exterior_ring(geom), want3d);
for (i = 1; i < nrings; i++)
pa[i] = ptarray_from_SFCGAL(sfcgal_polygon_interior_ring_n(geom, i - 1), want3d);
return (LWGEOM *)lwpoly_construct(srid, NULL, nrings, pa);
}
case SFCGAL_TYPE_MULTIPOINT:
case SFCGAL_TYPE_MULTILINESTRING:
case SFCGAL_TYPE_MULTIPOLYGON:
case SFCGAL_TYPE_MULTISOLID:
case SFCGAL_TYPE_GEOMETRYCOLLECTION:
{
ngeoms = sfcgal_geometry_collection_num_geometries(geom);
LWGEOM **geoms = NULL;
if (ngeoms)
{
nsolids = 0;
geoms = (LWGEOM **)lwalloc(sizeof(LWGEOM *) * ngeoms);
for (i = 0; i < ngeoms; i++)
{
const sfcgal_geometry_t *g = sfcgal_geometry_collection_geometry_n(geom, i);
geoms[i] = SFCGAL2LWGEOM(g, 0, srid);
if (FLAGS_GET_SOLID(geoms[i]->flags))
++nsolids;
}
geoms = (LWGEOM **)lwrealloc(geoms, sizeof(LWGEOM *) * ngeoms);
}
LWGEOM *rgeom = (LWGEOM *)lwcollection_construct(
SFCGAL_type_to_lwgeom_type(sfcgal_geometry_type_id(geom)), srid, NULL, ngeoms, geoms);
if (ngeoms)
{
if (ngeoms == nsolids)
FLAGS_SET_SOLID(rgeom->flags, 1);
else if (nsolids)
lwnotice(
"SFCGAL2LWGEOM: SOLID in heterogeneous collection will be treated as a POLYHEDRALSURFACETYPE");
}
return rgeom;
}
#if 0
case SFCGAL_TYPE_CIRCULARSTRING:
case SFCGAL_TYPE_COMPOUNDCURVE:
case SFCGAL_TYPE_CURVEPOLYGON:
case SFCGAL_TYPE_MULTICURVE:
case SFCGAL_TYPE_MULTISURFACE:
case SFCGAL_TYPE_CURVE:
case SFCGAL_TYPE_SURFACE:
/* TODO curve types handling */
#endif
case SFCGAL_TYPE_POLYHEDRALSURFACE:
{
ngeoms = sfcgal_polyhedral_surface_num_polygons(geom);
LWGEOM **geoms = NULL;
if (ngeoms)
{
geoms = (LWGEOM **)lwalloc(sizeof(LWGEOM *) * ngeoms);
for (i = 0; i < ngeoms; i++)
{
const sfcgal_geometry_t *g = sfcgal_polyhedral_surface_polygon_n(geom, i);
geoms[i] = SFCGAL2LWGEOM(g, 0, srid);
}
}
return (LWGEOM *)lwcollection_construct(POLYHEDRALSURFACETYPE, srid, NULL, ngeoms, geoms);
}
/* Solid is map as a closed PolyhedralSurface (for now) */
case SFCGAL_TYPE_SOLID:
{
nshells = sfcgal_solid_num_shells(geom);
for (ngeoms = 0, i = 0; i < nshells; i++)
ngeoms += sfcgal_polyhedral_surface_num_polygons(sfcgal_solid_shell_n(geom, i));
LWGEOM **geoms = 0;
if (ngeoms)
{
geoms = (LWGEOM **)lwalloc(sizeof(LWGEOM *) * ngeoms);
for (i = 0, k = 0; i < nshells; i++)
{
const sfcgal_geometry_t *shell = sfcgal_solid_shell_n(geom, i);
ngeoms = sfcgal_polyhedral_surface_num_polygons(shell);
for (j = 0; j < ngeoms; j++)
{
const sfcgal_geometry_t *g = sfcgal_polyhedral_surface_polygon_n(shell, j);
geoms[k] = SFCGAL2LWGEOM(g, 1, srid);
k++;
}
}
}
LWGEOM *rgeom = (LWGEOM *)lwcollection_construct(POLYHEDRALSURFACETYPE, srid, NULL, ngeoms, geoms);
if (ngeoms)
FLAGS_SET_SOLID(rgeom->flags, 1);
return rgeom;
}
case SFCGAL_TYPE_TRIANGULATEDSURFACE:
{
ngeoms = sfcgal_triangulated_surface_num_triangles(geom);
LWGEOM **geoms = NULL;
if (ngeoms)
{
geoms = (LWGEOM **)lwalloc(sizeof(LWGEOM *) * ngeoms);
for (i = 0; i < ngeoms; i++)
{
const sfcgal_geometry_t *g = sfcgal_triangulated_surface_triangle_n(geom, i);
geoms[i] = SFCGAL2LWGEOM(g, 0, srid);
}
}
return (LWGEOM *)lwcollection_construct(TINTYPE, srid, NULL, ngeoms, geoms);
}
default:
lwerror("SFCGAL2LWGEOM: Unknown Type");
return NULL;
}
}
sfcgal_geometry_t *
LWGEOM2SFCGAL(const LWGEOM *geom)
{
uint32_t i;
sfcgal_geometry_t *ret_geom = NULL;
assert(geom);
switch (geom->type)
{
case POINTTYPE:
{
const LWPOINT *lwp = (const LWPOINT *)geom;
if (lwgeom_is_empty(geom))
return sfcgal_point_create();
return ptarray_to_SFCGAL(lwp->point, POINTTYPE);
}
break;
case LINETYPE:
{
const LWLINE *line = (const LWLINE *)geom;
if (lwgeom_is_empty(geom))
return sfcgal_linestring_create();
return ptarray_to_SFCGAL(line->points, LINETYPE);
}
break;
case TRIANGLETYPE:
{
const LWTRIANGLE *triangle = (const LWTRIANGLE *)geom;
if (lwgeom_is_empty(geom))
return sfcgal_triangle_create();
return ptarray_to_SFCGAL(triangle->points, TRIANGLETYPE);
}
break;
case POLYGONTYPE:
{
const LWPOLY *poly = (const LWPOLY *)geom;
uint32_t nrings = poly->nrings - 1;
if (lwgeom_is_empty(geom))
return sfcgal_polygon_create();
sfcgal_geometry_t *exterior_ring = ptarray_to_SFCGAL(poly->rings[0], LINETYPE);
ret_geom = sfcgal_polygon_create_from_exterior_ring(exterior_ring);
for (i = 0; i < nrings; i++)
{
sfcgal_geometry_t *ring = ptarray_to_SFCGAL(poly->rings[i + 1], LINETYPE);
sfcgal_polygon_add_interior_ring(ret_geom, ring);
}
return ret_geom;
}
break;
case MULTIPOINTTYPE:
case MULTILINETYPE:
case MULTIPOLYGONTYPE:
case COLLECTIONTYPE:
{
if (geom->type == MULTIPOINTTYPE)
ret_geom = sfcgal_multi_point_create();
else if (geom->type == MULTILINETYPE)
ret_geom = sfcgal_multi_linestring_create();
else if (geom->type == MULTIPOLYGONTYPE)
ret_geom = sfcgal_multi_polygon_create();
else
ret_geom = sfcgal_geometry_collection_create();
const LWCOLLECTION *lwc = (const LWCOLLECTION *)geom;
for (i = 0; i < lwc->ngeoms; i++)
{
sfcgal_geometry_t *g = LWGEOM2SFCGAL(lwc->geoms[i]);
sfcgal_geometry_collection_add_geometry(ret_geom, g);
}
return ret_geom;
}
break;
case POLYHEDRALSURFACETYPE:
{
const LWPSURFACE *lwp = (const LWPSURFACE *)geom;
ret_geom = sfcgal_polyhedral_surface_create();
for (i = 0; i < lwp->ngeoms; i++)
{
sfcgal_geometry_t *g = LWGEOM2SFCGAL((const LWGEOM *)lwp->geoms[i]);
sfcgal_polyhedral_surface_add_polygon(ret_geom, g);
}
/* We treat polyhedral surface as the only exterior shell,
since we can't distinguish exterior from interior shells ... */
if (FLAGS_GET_SOLID(lwp->flags))
{
return sfcgal_solid_create_from_exterior_shell(ret_geom);
}
return ret_geom;
}
break;
case TINTYPE:
{
const LWTIN *lwp = (const LWTIN *)geom;
ret_geom = sfcgal_triangulated_surface_create();
for (i = 0; i < lwp->ngeoms; i++)
{
sfcgal_geometry_t *g = LWGEOM2SFCGAL((const LWGEOM *)lwp->geoms[i]);
sfcgal_triangulated_surface_add_triangle(ret_geom, g);
}
return ret_geom;
}
break;
default:
lwerror("LWGEOM2SFCGAL: Unsupported geometry type %s !", lwtype_name(geom->type));
return NULL;
}
}
/*
* No Operation SFCGAL function, used (only) for cunit tests
* Take a PostGIS geometry, send it to SFCGAL and return it unchanged (in theory)
*/
LWGEOM *
lwgeom_sfcgal_noop(const LWGEOM *geom_in)
{
sfcgal_geometry_t *converted;
assert(geom_in);
converted = LWGEOM2SFCGAL(geom_in);
assert(converted);
LWGEOM *geom_out = SFCGAL2LWGEOM(converted, 0, SRID_UNKNOWN);
sfcgal_geometry_delete(converted);
/* copy SRID (SFCGAL does not store the SRID) */
geom_out->srid = geom_in->srid;
return geom_out;
}