6

Working with the union and polygonize methods in GDAL/OGR, I often end up with self intersecting polygons. The case is one polygon which essentially consists of two polygons:

Self intersecting polygon example

In the above example, the polygon has two vertices at the intersection point, but whether the lines meet and goes seperate ways at the point or crosses each other, is unknown.

  1. Are there methods for splitting a self intersecting polygon like above, into two polygons or a multipolygon?

  2. Alternatively are there any existing algorithms I could implement my self (I would rather not invent another wheel)?

2

It's strange to get such self intersecting polygons because it violates the polygon definition in OGC Simple feature access specification and may not be viewed as legal polygons.

I don't use GDAL that often, but am quite familiar with GEOS and if I remember clearly GDAL union/polygonize are built based on GEOS library. In GEOS, if you union two polygons that share one vertex coordinates, like the case in your given picture, it will return you a multipolygon. In case you might want to look into the code, here is how GEOS do it:

http://svn.osgeo.org/geos/branches/3.3/src/geom/

The relevant functions may be the Geometry::Union in file Geometry.cpp

And also the GeometryFactory::buildGeometry in file GeometryFactory.cpp

1

Here is a link to an older discussion regarding converting MultiPart geometry to SinglePart geometry, utilizing osgeo.ogr. About 3/4 toward the end of the discussion is a Python script.

multipart to Single part Python script

  • As I see it, the script in your link just points to how I can explode a MultiPolygon to individual polygons. In my case, I have only one ring. – Chau Jun 10 '13 at 9:11
  • Maybe you can test for this self intersecting condition using OGR IsValid() or IsSimple() and then Buffer(0.0) to write out correct geometries. – klewis Jun 10 '13 at 12:23
  • The buffer operation is just too expensive for me. – Chau Jun 10 '13 at 13:01
  • Can you count the occurrences of XY's > 1 to find the touching point, excluding the start/stop points? Then split the list at that Point to generate distinct polygons? – klewis Jun 10 '13 at 18:01
  • The whole idea was to use existing algorithms and not begin from the start. If I have to invent the wheel again, I will do that :) – Chau Jun 11 '13 at 9:39

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.