Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Using JSTS buffer to identify a self-intersecting polygon

I was expecting to be able to test for self intersecting polygons either by JSTS failing to construct them or by adding a buffer and testing whether they were MultiPolygons after buffering but for a certain shape this isn't working and that's well past my geometry ability to grok

//a self-intersecting shape
var poly = [[0, 3], [1, 5], [3, 1], [5, 5], [6, 3], [0, 3]];

var geomFactory = new jsts.geom.GeometryFactory();

var jstsCoordinates = poly.map(function(pt) {
  return new jsts.geom.Coordinate(pt[0], pt[1]);
});

var linearRing = geomFactory.createLinearRing(jstsCoordinates);

var jstsPolygon = geomFactory.createPolygon(linearRing).buffer(1);

console.log(jstsPolygon.getGeometryType()); //this will be polygon but I thought should be MultiPolygon

var bufferedPoly = (jstsPolygon.shell.points.coordinates.map(function(pr) {
  return [pr.x, pr.y]
}))

var svg = d3.select('svg');

//add the first shape (maginified for display)
svg.selectAll('.original').data([poly]).enter().append("polygon")
    .attr("points",function(d) { 
        return d.map(function(d) {
            return [d[0] * 10 + 10, d[1]*10].join(",");
        }).join(" ");
    })
    .attr("fill", "yellow");

//add the buffered shape below it
svg.selectAll('.buffered').data([bufferedPoly]).enter().append("polygon")
    .attr("points",function(d) { 
        return d.map(function(d) {
            return [d[0] * 10 + 10, d[1]*10 + 40].join(",");
        }).join(" ");
    })
    .attr("fill", "yellow");
svg {background-color:blue}
<script src="https://cdnjs.cloudflare.com/ajax/libs/d3/3.4.11/d3.min.js"></script>
<script src="https://cdn.rawgit.com/bjornharrtell/jsts/gh-pages/1.0.2/jsts.min.js"></script>
<svg width=200 height=200>
  <svg>
like image 496
Paul D'Ambra Avatar asked Mar 20 '16 20:03

Paul D'Ambra


2 Answers

Based on the linked answer using JTS (Is there a way to convert a self intersecting polygon to a multipolygon in JTS?), I re-emplemented it using JSTS

just wanted to post the code for other people who seek the same answer / functionality.

By running the code, using a JSTS Geometry, you will be returned a valid Geometry either being a Polygon or Multipolygon

/**
 * Get / create a valid version of the geometry given. If the geometry is a polygon or multi polygon, self intersections /
 * inconsistencies are fixed. Otherwise the geometry is returned.
 *
 * @param geom
 * @return a geometry
 */
var jsts_validate = function(geom) {
  if (geom instanceof jsts.geom.Polygon) {
    if (geom.isValid()) {
      geom.normalize(); // validate does not pick up rings in the wrong order - this will fix that
      return geom; // If the polygon is valid just return it
    }
    var polygonizer = new jsts.operation.polygonize.Polygonizer();
    jsts_addPolygon(geom, polygonizer);
    return jsts_toPolygonGeometry(polygonizer.getPolygons(), geom.getFactory());
  } else if (geom instanceof jsts.geom.MultiPolygon) {
    if (geom.isValid()) {
      geom.normalize(); // validate does not pick up rings in the wrong order - this will fix that
      return geom; // If the multipolygon is valid just return it
    }
    var polygonizer = new jsts.operation.polygonize.Polygonizer();

    for (var n = geom.getNumGeometries(); n > 0; n--) {
      jsts_addPolygon(geom.getGeometryN(n - 1), polygonizer);
    }
    return jsts_toPolygonGeometry(polygonizer.getPolygons(), geom.getFactory());
  } else {
    return geom; // In my case, I only care about polygon / multipolygon geometries
  }
};

/**
 * Add all line strings from the polygon given to the polygonizer given
 *
 * @param polygon polygon from which to extract line strings
 * @param polygonizer polygonizer
 */
var jsts_addPolygon = function(polygon, polygonizer) {
  jsts_addLineString(polygon.getExteriorRing(), polygonizer);

  for (var n = polygon.getNumInteriorRing(); n > 0; n--) {
    jsts_addLineString(polygon.getInteriorRingN(n), polygonizer);
  }
};

/**
 * Add the linestring given to the polygonizer
 *
 * @param linestring line string
 * @param polygonizer polygonizer
 */
var jsts_addLineString = function(lineString, polygonizer) {

  if (lineString instanceof jsts.geom.LinearRing) {
    // LinearRings are treated differently to line strings : we need a LineString NOT a LinearRing
    lineString = lineString.getFactory().createLineString(lineString.getCoordinateSequence());
  }

  // unioning the linestring with the point makes any self intersections explicit.
  var point = lineString.getFactory().createPoint(lineString.getCoordinateN(0));
  var toAdd = lineString.union(point); //geometry

  //Add result to polygonizer
  polygonizer.add(toAdd);
}

/**
 * Get a geometry from a collection of polygons.
 *
 * @param polygons collection
 * @param factory factory to generate MultiPolygon if required
 * @return null if there were no polygons, the polygon if there was only one, or a MultiPolygon containing all polygons otherwise
 */
var jsts_toPolygonGeometry = function(polygons, factory) {
  switch (polygons.size()) {
    case 0:
      return null; // No valid polygons!
    case 1:
      return polygons.iterator().next(); // single polygon - no need to wrap
    default:
      //polygons may still overlap! Need to sym difference them
      var iter = polygons.iterator();
      var ret = iter.next();
      while (iter.hasNext()) {
        ret = ret.symDifference(iter.next());
      }
      return ret;
  }
}
like image 186
Martin Kirk Avatar answered Sep 29 '22 14:09

Martin Kirk


And just getting the kids to bed gave my brain enough of a rest...

There are two solutions here (I think - this geometry stuff is way outside my comfort zone)

I found Splitting self-intersecting polygon only returned one polygon in shapely in Python and Is there a way to convert a self intersecting polygon to a multipolygon in JTS? which actually contain the solution(s) since jts, jsts, and shapely are all closely related.

The first is that, having constructed a linear ring (that in this case is not simple), I can call isSimple() to receive false

and the second that I had been calling buffer(1) having misunderstood the advice I'd been given. The solution there being to call buffer(0).

//a self-intersecting shape
var poly = [
  [0, 3],
  [1, 5],
  [3, 1],
  [5, 5],
  [6, 3],
  [0, 3]
];

var geomFactory = new jsts.geom.GeometryFactory();

var jstsCoordinates = poly.map(function(pt) {
  return new jsts.geom.Coordinate(pt[0], pt[1]);
});

var linearRing = geomFactory.createLinearRing(jstsCoordinates);
// turns out you can just ask if it is simple... i.e. does not have any self intersections.
console.log(linearRing.isSimple()); //so this is false

//Ah! To split it and find out if it is self intersecting use buffer(0)
var jstsPolygon = geomFactory.createPolygon(linearRing).buffer(0);

console.log(jstsPolygon.getGeometryType()); //this will now be MultiPolygon

var svg = d3.select('svg');

if (jstsPolygon.getGeometryType() !== 'MultiPolygon') {
var bufferedPoly = (jstsPolygon.shell.points.coordinates.map(function(pr) {
  return [pr.x, pr.y];
}));

//add the buffered shape below it
svg.selectAll('.buffered').data([bufferedPoly]).enter().append("polygon")
    .attr("points",function(d) { 
        return d.map(function(d) {
            return [d[0] * 10 + 10, d[1]*10 + 40].join(",");
        }).join(" ");
    })
    .attr("fill", "yellow");
}

//add the first shape (magnified for display)
svg.selectAll('.original')
   .data([poly]).enter()
   .append("polygon")
   .attr("points",function(d) { 
        return d.map(function(d) {
            return [
              d[0] * 10 + 10, 
              d[1]*10].join(",");
        }).join(" ");
    })
    .attr("fill", "yellow");
svg {
  background-color: blue
}
<script src="https://cdnjs.cloudflare.com/ajax/libs/d3/3.4.11/d3.min.js"></script>
<script src="https://cdn.rawgit.com/bjornharrtell/jsts/gh-pages/1.0.2/jsts.min.js"></script>
<svg width=200 height=200>
  <svg>
like image 41
Paul D'Ambra Avatar answered Sep 29 '22 14:09

Paul D'Ambra