/**********************************************************************
*
* 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 (C) 2001-2006 Refractions Research Inc.
*
**********************************************************************/
#include
#include
#include
#include "liblwgeom_internal.h"
#include "lwgeom_log.h"
LWPSURFACE* lwpsurface_add_lwpoly(LWPSURFACE *mobj, const LWPOLY *obj)
{
return (LWPSURFACE*)lwcollection_add_lwgeom((LWCOLLECTION*)mobj, (LWGEOM*)obj);
}
void lwpsurface_free(LWPSURFACE *psurf)
{
uint32_t i;
if ( ! psurf ) return;
if ( psurf->bbox )
lwfree(psurf->bbox);
for ( i = 0; i < psurf->ngeoms; i++ )
if ( psurf->geoms && psurf->geoms[i] )
lwpoly_free(psurf->geoms[i]);
if ( psurf->geoms )
lwfree(psurf->geoms);
lwfree(psurf);
}
void printLWPSURFACE(LWPSURFACE *psurf)
{
uint32_t i, j;
LWPOLY *patch;
if (psurf->type != POLYHEDRALSURFACETYPE)
lwerror("printLWPSURFACE called with something else than a POLYHEDRALSURFACE");
lwnotice("LWPSURFACE {");
lwnotice(" ndims = %i", (int)FLAGS_NDIMS(psurf->flags));
lwnotice(" SRID = %i", (int)psurf->srid);
lwnotice(" ngeoms = %i", (int)psurf->ngeoms);
for (i=0; ingeoms; i++)
{
patch = (LWPOLY *) psurf->geoms[i];
for (j=0; jnrings; j++)
{
lwnotice(" RING # %i :",j);
printPA(patch->rings[j]);
}
}
lwnotice("}");
}
/*
* TODO rewrite all this stuff to be based on a truly topological model
*/
struct struct_psurface_arcs
{
double ax, ay, az;
double bx, by, bz;
uint32_t cnt, face;
};
typedef struct struct_psurface_arcs *psurface_arcs;
/* We supposed that the geometry is valid
we could have wrong result if not */
int lwpsurface_is_closed(const LWPSURFACE *psurface)
{
uint32_t i, j, k;
uint32_t narcs, carc;
int found;
psurface_arcs arcs;
POINT4D pa, pb;
LWPOLY *patch;
/* If surface is not 3D, it's can't be closed */
if (!FLAGS_GET_Z(psurface->flags)) return 0;
/* If surface is less than 4 faces hard to be closed too */
if (psurface->ngeoms < 4) return 0;
/* Max theoretical arcs number if no one is shared ... */
for (i=0, narcs=0 ; i < psurface->ngeoms ; i++)
{
patch = (LWPOLY *) psurface->geoms[i];
narcs += patch->rings[0]->npoints - 1;
}
arcs = lwalloc(sizeof(struct struct_psurface_arcs) * narcs);
for (i=0, carc=0; i < psurface->ngeoms ; i++)
{
patch = (LWPOLY *) psurface->geoms[i];
for (j=0; j < patch->rings[0]->npoints - 1; j++)
{
getPoint4d_p(patch->rings[0], j, &pa);
getPoint4d_p(patch->rings[0], j+1, &pb);
/* remove redundant points if any */
if (pa.x == pb.x && pa.y == pb.y && pa.z == pb.z) continue;
/* Make sure to order the 'lower' point first */
if ( (pa.x > pb.x) ||
(pa.x == pb.x && pa.y > pb.y) ||
(pa.x == pb.x && pa.y == pb.y && pa.z > pb.z) )
{
pa = pb;
getPoint4d_p(patch->rings[0], j, &pb);
}
for (found=0, k=0; k < carc ; k++)
{
if ( ( arcs[k].ax == pa.x && arcs[k].ay == pa.y &&
arcs[k].az == pa.z && arcs[k].bx == pb.x &&
arcs[k].by == pb.y && arcs[k].bz == pb.z &&
arcs[k].face != i) )
{
arcs[k].cnt++;
found = 1;
/* Look like an invalid PolyhedralSurface
anyway not a closed one */
if (arcs[k].cnt > 2)
{
lwfree(arcs);
return 0;
}
}
}
if (!found)
{
arcs[carc].cnt=1;
arcs[carc].face=i;
arcs[carc].ax = pa.x;
arcs[carc].ay = pa.y;
arcs[carc].az = pa.z;
arcs[carc].bx = pb.x;
arcs[carc].by = pb.y;
arcs[carc].bz = pb.z;
carc++;
/* Look like an invalid PolyhedralSurface
anyway not a closed one */
if (carc > narcs)
{
lwfree(arcs);
return 0;
}
}
}
}
/* A polyhedron is closed if each edge
is shared by exactly 2 faces */
for (k=0; k < carc ; k++)
{
if (arcs[k].cnt != 2)
{
lwfree(arcs);
return 0;
}
}
lwfree(arcs);
/* Invalid Polyhedral case */
if (carc < psurface->ngeoms) return 0;
return 1;
}