An algorithm for inflating/deflating (offsetting, buffering) polygons

I thought I might briefly mention my own polygon clipping and offsetting library - Clipper.

While Clipper is primarily designed for polygon clipping operations, it does polygon offsetting too. The library is open source freeware written in Delphi, C++ and C#. It has a very unencumbered Boost license allowing it to be used in both freeware and commercial applications without charge.

Polygon offsetting can be performed using one of three offset styles - squared, round and mitered.

https://i.stack.imgur.com/LqZeh.png

---

The polygon you are looking for is called inward/outward offset polygon in computational geometry and it is closely related to the straight skeleton.

These are several offset polygons for a complicated polygon:

https://i.stack.imgur.com/vCgNf.png

And this is the straight skeleton for another polygon:

https://i.stack.imgur.com/lETWh.png

As pointed out in other comments, as well, depending on how far you plan to "inflate/deflate" your polygon you can end up with different connectivity for the output.

From computation point of view: once you have the straight skeleton one should be able to construct the offset polygons relatively easily. The open source and (free for non-commercial) CGAL library has a package implementing these structures. See this code example to compute offset polygons using CGAL.

The package manual should give you a good starting point on how to construct these structures even if you are not going to use CGAL, and contains references to the papers with the mathematical definitions and properties:

CGAL manual: 2D Straight Skeleton and Polygon Offsetting

---

For these types of things I usually use JTS. For demonstration purposes I created this jsFiddle that uses JSTS (JavaScript port of JTS). You just need to convert the coordinates you have to JSTS coordinates:

function vectorCoordinates2JTS (polygon) {
  var coordinates = [];
  for (var i = 0; i < polygon.length; i++) {
    coordinates.push(new jsts.geom.Coordinate(polygon[i].x, polygon[i].y));
  }
  return coordinates;
}

The result is something like this:

https://i.stack.imgur.com/C7Wf8.png

Additional info: I usually use this type of inflating/deflating (a little modified for my purposes) for setting boundaries with radius on polygons that are drawn on a map (with Leaflet or Google maps). You just convert (lat,lng) pairs to JSTS coordinates and everything else is the same. Example:

https://i.stack.imgur.com/ryZNC.png

---

Sounds to me like what you want is:

Starting at a vertex, face anti-clockwise along an adjacent edge.

Replace the edge with a new, parallel edge placed at distance d to the "left" of the old one.

Repeat for all edges.

Find the intersections of the new edges to get the new vertices.

Detect if you've become a crossed polygon and decide what to do about it. Probably add a new vertex at the crossing-point and get rid of some old ones. I'm not sure whether there's a better way to detect this than just to compare every pair of non-adjacent edges to see if their intersection lies between both pairs of vertices.

The resulting polygon lies at the required distance from the old polygon "far enough" from the vertices. Near a vertex, the set of points at distance d from the old polygon is, as you say, not a polygon, so the requirement as stated cannot be fulfilled.

I don't know if this algorithm has a name, example code on the web, or a fiendish optimisation, but I think it describes what you want.

---

In the GIS world one uses negative buffering for this task: http://www-users.cs.umn.edu/~npramod/enc_pdf.pdf

The JTS library should do this for you. See the documentation for the buffer operation: http://tsusiatsoftware.net/jts/javadoc/com/vividsolutions/jts/operation/buffer/package-summary.html

For a rough overview see also the Developer Guide: http://www.vividsolutions.com/jts/bin/JTS%20Developer%20Guide.pdf

---

Each line should split the plane to "inside" and "outline"; you can find this out using the usual inner-product method.

Move all lines outward by some distance.

Consider all pair of neighbor lines (lines, not line segment), find the intersection. These are the new vertex.

Cleanup the new vertex by removing any intersecting parts. -- we have a few case here

(a) Case 1:

 0--7  4--3
 |  |  |  |
 |  6--5  |
 |        |
 1--------2

if you expend it by one, you got this:

0----a----3
|    |    |
|    |    |
|    b    |
|         |
|         |
1---------2

7 and 4 overlap.. if you see this, you remove this point and all points in between.

(b) case 2

 0--7  4--3
 |  |  |  |
 |  6--5  |
 |        |
 1--------2

if you expend it by two, you got this:

0----47----3
|    ||    |
|    ||    |
|    ||    |
|    56    |
|          |
|          |
|          |
1----------2

to resolve this, for each segment of line, you have to check if it overlap with latter segments.

(c) case 3

       4--3
 0--X9 |  |
 |  78 |  |
 |  6--5  |
 |        |
 1--------2

expend by 1. this is a more general case for case 1.

(d) case 4

same as case3, but expend by two.

Actually, if you can handle case 4. All other cases are just special case of it with some line or vertex overlapping.

To do case 4, you keep a stack of vertex.. you push when you find lines overlapping with latter line, pop it when you get the latter line. -- just like what you do in convex-hull.

---

Here is an alternative solution, see if you like this better.

Do a triangulation, it don't have to be delaunay -- any triangulation would do. Inflate each triangle -- this should be trivial. if you store the triangle in the anti-clockwise order, just move the lines to right-hand-side and do intersection. Merge them using a modified Weiler-Atherton clipping algorithm

---

Big thanks to Angus Johnson for his clipper library. There are good code samples for doing the clipping stuff at the clipper homepage at http://www.angusj.com/delphi/clipper.php#code but I did not see an example for polygon offsetting. So I thought that maybe it is of use for someone if I post my code:

    public static List<Point> GetOffsetPolygon(List<Point> originalPath, double offset)
    {
        List<Point> resultOffsetPath = new List<Point>();

        List<ClipperLib.IntPoint> polygon = new List<ClipperLib.IntPoint>();
        foreach (var point in originalPath)
        {
            polygon.Add(new ClipperLib.IntPoint(point.X, point.Y));
        }

        ClipperLib.ClipperOffset co = new ClipperLib.ClipperOffset();
        co.AddPath(polygon, ClipperLib.JoinType.jtRound, ClipperLib.EndType.etClosedPolygon);

        List<List<ClipperLib.IntPoint>> solution = new List<List<ClipperLib.IntPoint>>();
        co.Execute(ref solution, offset);

        foreach (var offsetPath in solution)
        {
            foreach (var offsetPathPoint in offsetPath)
            {
                resultOffsetPath.Add(new Point(Convert.ToInt32(offsetPathPoint.X), Convert.ToInt32(offsetPathPoint.Y)));
            }
        }

        return resultOffsetPath;
    }

---

One further option is to use boost::polygon - the documentation is somewhat lacking, but you should find that the methods resize and bloat, and also the overloaded += operator, which actually implement buffering. So for example increasing the size of a polygon (or a set of polygons) by some value can be as simple as:

poly += 2; // buffer polygon by 2

---

Based on advice from @JoshO'Brian, it appears the rGeos package in the R language implements this algorithm. See rGeos::gBuffer .

---

There are a couple of libraries one can use (Also usable for 3D data sets).

https://github.com/otherlab/openmesh https://github.com/alecjacobson/nested_cages http://homepage.tudelft.nl/h05k3/Projects/MeshThickeningProj.htm

One can also find corresponding publications for these libraries to understand the algorithms in more detail.

The last one has the least dependencies and is self-contained and can read in .obj files.

Best wishes, Stephan

---

I use simple geometry: Vectors and/or Trigonometry

At each corner find the mid vector, and mid angle. Mid vector is the arithmetic average of the two unit vectors defined by the edges of the corner. Mid Angle is the half of the angle defined by the edges. If you need to expand (or contract) your polygon by the amount of d from each edge; you should go out (in) by the amount d/sin(midAngle) to get the new corner point. Repeat this for all the corners

*** Be careful about your direction. Make CounterClockWise Test using the three points defining the corner; to find out which way is out, or in.

---

This is a C# implementation of the algorithm explained in here . It is using Unity math library and collection package as well.

public static NativeArray<float2> ExpandPoly(NativeArray<float2> oldPoints, float offset, int outer_ccw = 1)
{
    int num_points = oldPoints.Length;
    NativeArray<float2> newPoints = new NativeArray<float2>(num_points, Allocator.Temp);

    for (int curr = 0; curr < num_points; curr++)
    {
        int prev = (curr + num_points - 1) % num_points;
        int next = (curr + 1) % num_points;

        float2 vn = oldPoints[next] - oldPoints[curr];
        float2 vnn = math.normalize(vn);
        float nnnX = vnn.y;
        float nnnY = -vnn.x;

        float2 vp = oldPoints[curr] - oldPoints[prev];
        float2 vpn = math.normalize(vp);
        float npnX = vpn.y * outer_ccw;
        float npnY = -vpn.x * outer_ccw;

        float bisX = (nnnX + npnX) * outer_ccw;
        float bisY = (nnnY + npnY) * outer_ccw;

        float2 bisn = math.normalize(new float2(bisX, bisY));
        float bislen = offset / math.sqrt(1 + nnnX * npnX + nnnY * npnY);

        newPoints[curr] = new float2(oldPoints[curr].x + bislen * bisn.x, oldPoints[curr].y + bislen * bisn.y);
    }

    return newPoints;
}

---