LWGEOM_GEOS_buildArea()

GEOSGeometry* LWGEOM_GEOS_buildArea

(

const GEOSGeometry *

geom_in

)

Definition at line 1252 of file liblwgeom/lwgeom_geos.c.

References collectFacesWithEvenAncestors(), COLLECTIONTYPE, delFace(), findFaceHoles(), GEOS2LWGEOM(), lwalloc(), LWDEBUGF, lwerror(), lwfree(), lwgeom_to_ewkt(), lwnotice(), and newFace().

Referenced by lwgeom_buildarea(), and LWGEOM_GEOS_makeValidPolygon().

{
  GEOSGeometry *tmp;
  GEOSGeometry *geos_result, *shp;
  GEOSGeometry const *vgeoms[1];
  uint32_t i, ngeoms;
  int srid = GEOSGetSRID(geom_in);
  Face ** geoms;

  vgeoms[0] = geom_in;
#ifdef LWGEOM_PROFILE_BUILDAREA
  lwnotice("Polygonizing");
#endif
  geos_result = GEOSPolygonize(vgeoms, 1);

  LWDEBUGF(3, "GEOSpolygonize returned @ %p", geos_result);

  /* Null return from GEOSpolygonize (an exception) */
  if ( ! geos_result ) return 0;

  /*
   * We should now have a collection
   */
#if PARANOIA_LEVEL > 0
  if ( GEOSGeomTypeId(geos_result) != COLLECTIONTYPE )
  {
    GEOSGeom_destroy(geos_result);
    lwerror("%s [%d] Unexpected return from GEOSpolygonize", __FILE__, __LINE__);
    return 0;
  }
#endif

  ngeoms = GEOSGetNumGeometries(geos_result);
#ifdef LWGEOM_PROFILE_BUILDAREA
  lwnotice("Num geometries from polygonizer: %d", ngeoms);
#endif


  LWDEBUGF(3, "GEOSpolygonize: ngeoms in polygonize output: %d", ngeoms);
  LWDEBUGF(3, "GEOSpolygonize: polygonized:%s",
              lwgeom_to_ewkt(GEOS2LWGEOM(geos_result, 0)));

  /*
   * No geometries in collection, early out
   */
  if ( ngeoms == 0 )
  {
    GEOSSetSRID(geos_result, srid);
    return geos_result;
  }

  /*
   * Return first geometry if we only have one in collection,
   * to avoid the unnecessary Geometry clone below.
   */
  if ( ngeoms == 1 )
  {
    tmp = (GEOSGeometry *)GEOSGetGeometryN(geos_result, 0);
    if ( ! tmp )
    {
      GEOSGeom_destroy(geos_result);
      return 0; /* exception */
    }
    shp = GEOSGeom_clone(tmp);
    GEOSGeom_destroy(geos_result); /* only safe after the clone above */
    GEOSSetSRID(shp, srid);
    return shp;
  }

  LWDEBUGF(2, "Polygonize returned %d geoms", ngeoms);

  /*
   * Polygonizer returns a polygon for each face in the built topology.
   *
   * This means that for any face with holes we'll have other faces
   * representing each hole. We can imagine a parent-child relationship
   * between these faces.
   *
   * In order to maximize the number of visible rings in output we
   * only use those faces which have an even number of parents.
   *
   * Example:
   *
   *   +---------------+
   *   |     L0        |  L0 has no parents
   *   |  +---------+  |
   *   |  |   L1    |  |  L1 is an hole of L0
   *   |  |  +---+  |  |
   *   |  |  |L2 |  |  |  L2 is an hole of L1 (which is an hole of L0)
   *   |  |  |   |  |  |
   *   |  |  +---+  |  |
   *   |  +---------+  |
   *   |               |
   *   +---------------+
   *
   * See http://trac.osgeo.org/postgis/ticket/1806
   *
   */

#ifdef LWGEOM_PROFILE_BUILDAREA
  lwnotice("Preparing face structures");
#endif

  /* Prepare face structures for later analysis */
  geoms = lwalloc(sizeof(Face**)*ngeoms);
  for (i=0; i<ngeoms; ++i)
    geoms[i] = newFace(GEOSGetGeometryN(geos_result, i));

#ifdef LWGEOM_PROFILE_BUILDAREA
  lwnotice("Finding face holes");
#endif

  /* Find faces representing other faces holes */
  findFaceHoles(geoms, ngeoms);

#ifdef LWGEOM_PROFILE_BUILDAREA
  lwnotice("Colletting even ancestor faces");
#endif

  /* Build a MultiPolygon composed only by faces with an
   * even number of ancestors */
  tmp = collectFacesWithEvenAncestors(geoms, ngeoms);

#ifdef LWGEOM_PROFILE_BUILDAREA
  lwnotice("Cleaning up");
#endif

  /* Cleanup face structures */
  for (i=0; i<ngeoms; ++i) delFace(geoms[i]);
  lwfree(geoms);

  /* Faces referenced memory owned by geos_result.
   * It is safe to destroy geos_result after deleting them. */
  GEOSGeom_destroy(geos_result);

#ifdef LWGEOM_PROFILE_BUILDAREA
  lwnotice("Self-unioning");
#endif

  /* Run a single overlay operation to dissolve shared edges */
  shp = GEOSUnionCascaded(tmp);
  if ( ! shp )
  {
    GEOSGeom_destroy(tmp);
    return 0; /* exception */
  }

#ifdef LWGEOM_PROFILE_BUILDAREA
  lwnotice("Final cleanup");
#endif

  GEOSGeom_destroy(tmp);

  GEOSSetSRID(shp, srid);

  return shp;
}

Here is the call graph for this function:

Here is the caller graph for this function: