Expand fill of convex polygon

I have a convex polygon P1 of N points. This polygon could be any shape or proportion (as long as it is still convex).

I need to compute another polygon P2 usin开发者_如何转开发g the original polygons geometry, but "expanded" by a given number of units. What might the algorithm be for expanding a convex polygon?

---

To expand a convex polygon, draw a line parallel to each edge and the given number of units away. Then use the intersection points of the new lines as the vertices of the expanded polygon. The javascript/canvas at the end follows this functional breakdown:

Step 1: Figure out which side is "out"

The order of the vertices (points) matters. In a convex polygon, they can be listed in a clockwise (CW), or a counter-clockwise (CCW) order. In a CW polygon, turn one of the edges 90 degrees CCW to obtain an outward-facing normal. On a CCW polygon, turn it CW instead.

If the turn direction of the vertices is not known in advance, examine how the second edge turns from the first. In a convex polygon, the remaining edges will keep turning in the same direction:

1. Find the CW normal of the first edge. We don't know yet whether it's facing inward or outward.

2. Compute the dot product of the second edge with the normal we computed. If the second edge turns CW, the dot product will be positive. It will be negative otherwise.

Math:

// in vector terms:
v01 = p1 - p0                      // first edge, as a vector
v12 = p2 - p1                      // second edge, as a vector
n01 = (v01.y, -v01.x)              // CW normal of first edge
d = v12 * n01                      // dot product

// and in x,y terms:
v01 = (p1.x-p0.x, p1.y-p0.y)       // first edge, as a vector
v12 = (p2.x-p1.x, p2.y-p1.y)       // second edge, as a vector
n01 = (v01.y, -v01.x)              // CW normal of first edge
d = v12.x * n01.x + v12.y * n01.y; // dot product: v12 * n01

if (d > 0) {
    // the polygon is CW
} else {
    // the polygon is CCW
}

// and what if d==0 ?
// -- that means the second edge continues in the same
// direction as a first.  keep looking for an edge that
// actually turns either CW or CCW.

Code:

function vecDot(v1, v2) {
    return v1.x * v2.x + v1.y * v2.y;
}

function vecRot90CW(v) {
    return { x: v.y, y: -v.x };
}

function vecRot90CCW(v) {
    return { x: -v.y, y: v.x };
}

function polyIsCw(p) {
    return vecDot(
        vecRot90CW({ x: p[1].x - p[0].x, y: p[1].y - p[0].y }),
        { x: p[2].x - p[1].x, y: p[2].y - p[1].y }) >= 0;
}

var rot = polyIsCw(p) ? vecRot90CCW : vecRot90CW;

Step 2: Find lines parallel to the polygon edges

Now that we know which side is out, we can compute lines parallel to each polygon edge, at exactly the required distance. Here's our strategy:

1. For each edge, compute its outward-facing normal

2. Normalize the normal, such that its length becomes one unit

3. Multiply the normal by the distance we want the expanded polygon to be from the original

4. Add the multiplied normal to both ends of the edge. That will give us two points on the parallel line. Those two points are enough to define the parallel line.

Code:

// given two vertices pt0 and pt1, a desired distance, and a function rot()
// that turns a vector 90 degrees outward:

function vecUnit(v) {
    var len = Math.sqrt(v.x * v.x + v.y * v.y);
    return { x: v.x / len, y: v.y / len };
}

function vecMul(v, s) {
    return { x: v.x * s, y: v.y * s };
}

var v01 = { x: pt1.x - pt0.x, y: pt1.y - pt0.y };  // edge vector
var d01 = vecMul(vecUnit(rot(v01)), distance);     // multiplied unit normal
var ptx0 = { x: pt0.x + d01.x, y: pt0.y + d01.y }; // two points on the
var ptx1 = { x: pt1.x + d01.x, y: pt1.y + d01.y }; //     parallel line

Step 3: Compute the intersections of the parallel lines

--these will be the vertices of the expanded polygon.

Math:

A line going through two points P1, P2 can be described as:

P = P1 + t * (P2 - P1)

Two lines can be described as

P = P1 + t * (P2 - P1)
P = P3 + u * (P4 - P3)

And their intersection has to be on both lines:

P = P1 + t * (P2 - P1) = P3 + u * (P4 - P3)

This can be massaged to look like:

(P2 - P1) * t + (P3 - P4) * u = P3 - P1

Which in x,y terms is:

(P2.x - P1.x) * t + (P3.x - P4.x) * u = P3.x - P1.x
(P2.y - P1.y) * t + (P3.y - P4.y) * u = P3.y - P1.y

As the points P1, P2, P3 and P4 are known, so are the following values:

a1 = P2.x - P1.x    a2 = P2.y - P1.y
b1 = P3.x - P4.x    b2 = P3.y - P4.y
c1 = P3.x - P1.x    c2 = P3.y - P1.y

This shortens our equations to:

a1*t + b1*u = c1
a2*t + b2*u = c2    

Solving for t gets us:

t = (b1*c2 - b2*c1)/(a2*b1 - a1*b2)

Which lets us find the intersection at P = P1 + t * (P2 - P1).

Code:

function intersect(line1, line2) {
    var a1 = line1[1].x - line1[0].x;
    var b1 = line2[0].x - line2[1].x;
    var c1 = line2[0].x - line1[0].x;

    var a2 = line1[1].y - line1[0].y;
    var b2 = line2[0].y - line2[1].y;
    var c2 = line2[0].y - line1[0].y;

    var t = (b1*c2 - b2*c1) / (a2*b1 - a1*b2);

    return {
        x: line1[0].x + t * (line1[1].x - line1[0].x),
        y: line1[0].y + t * (line1[1].y - line1[0].y)
    };
}

Step 4: Deal with special cases

There is a number of special cases that merit attention. Left as an exercise to the reader...

1. When there's a very sharp angle between two edges, the expanded vertex can be very far from the original one. You might want to consider clipping the expanded edge if it goes beyond some threshold. At the extreme case, the angle is zero, which suggests that the expanded vertex is at infinity, causing division by zero in the arithmetic. Watch out.

2. When the first two edges are on the same line, you can't tell if it's a CW or a CCW polygon by looking just at them. Look at more edges.

3. Non convex polygons are much more interesting... and are not tackled here.

Full sample code

Drop this in a canvas-capable browser. I used Chrome 6 on Windows. The triangle and its expanded version should animate.

canvas { border: 1px solid #ccc; }

$(function() {
      var canvas = document.getElementById('canvas');  
      if (canvas.getContext) {  
          var context = canvas.getContext('2d');  

          // math for expanding a polygon

          function vecUnit(v) {
              var len = Math.sqrt(v.x * v.x + v.y * v.y);
              return { x: v.x / len, y: v.y / len };
          }

          function vecMul(v, s) {
              return { x: v.x * s, y: v.y * s };
          }

          function vecDot(v1, v2) {
              return v1.x * v2.x + v1.y * v2.y;
          }

          function vecRot90CW(v) {
              return { x: v.y, y: -v.x };
          }

          function vecRot90CCW(v) {
              return { x: -v.y, y: v.x };
          }

          function intersect(line1, line2) {
              var a1 = line1[1].x - line1[0].x;
              var b1 = line2[0].x - line2[1].x;
              var c1 = line2[0].x - line1[0].x;

              var a2 = line1[1].y - line1[0].y;
              var b2 = line2[0].y - line2[1].y;
              var c2 = line2[0].y - line1[0].y;

              var t = (b1*c2 - b2*c1) / (a2*b1 - a1*b2);

              return {
                  x: line1[0].x + t * (line1[1].x - line1[0].x),
                  y: line1[0].y + t * (line1[1].y - line1[0].y)
              };
          }

          function polyIsCw(p) {
              return vecDot(
                  vecRot90CW({ x: p[1].x - p[0].x, y: p[1].y - p[0].y }),
                  { x: p[2].x - p[1].x, y: p[2].y - p[1].y }) >= 0;
          }

          function expandPoly(p, distance) {
              var expanded = [];
              var rot = polyIsCw(p) ? vecRot90CCW : vecRot90CW;

              for (var i = 0; i < p.length; ++i) {

                  // get this point (pt1), the point before it
                  // (pt0) and the point that follows it (pt2)
                  var pt0 = p[(i > 0) ? i - 1 : p.length - 1];
                  var pt1 = p[i];
                  var pt2 = p[(i < p.length - 1) ? i + 1 : 0];

                  // find the line vectors of the lines going
                  // into the current point
                  var v01 = { x: pt1.x - pt0.x, y: pt1.y - pt0.y };
                  var v12 = { x: pt2.x - pt1.x, y: pt2.y - pt1.y };

                  // find the normals of the two lines, multiplied
                  // to the distance that polygon should inflate
                  var d01 = vecMul(vecUnit(rot(v01)), distance);
                  var d12 = vecMul(vecUnit(rot(v12)), distance);

                  // use the normals to find two points on the
                  // lines parallel to the polygon lines
                  var ptx0  = { x: pt0.x + d01.x, y: pt0.y + d01.y };
                  var ptx10 = { x: pt1.x + d01.x, y: pt1.y + d01.y };
                  var ptx12 = { x: pt1.x + d12.x, y: pt1.y + d12.y };
                  var ptx2  = { x: pt2.x + d12.x, y: pt2.y + d12.y };

                  // find the intersection of the two lines, and
                  // add it to the expanded polygon
                  expanded.push(intersect([ptx0, ptx10], [ptx12, ptx2]));
              }
              return expanded;
          }

          // drawing and animating a sample polygon on a canvas

          function drawPoly(p) {
              context.beginPath();
              context.moveTo(p[0].x, p[0].y);
              for (var i = 0; i < p.length; ++i) {
                  context.lineTo(p[i].x, p[i].y);
              }
              context.closePath();
              context.fill();
              context.stroke();
          }

          function drawPolyWithMargin(p, margin) {
              context.fillStyle = "rgb(255,255,255)";  
              context.strokeStyle = "rgb(200,150,150)";  
              drawPoly(expandPoly(p, margin));

              context.fillStyle = "rgb(150,100,100)";  
              context.strokeStyle = "rgb(200,150,150)";  
              drawPoly(p);
          }

          var p = [{ x: 100, y: 100 }, { x: 200, y: 120 }, { x: 80, y: 200 }];
          setInterval(function() {
              for (var i in p) {
                  var pt = p[i];
                  if (pt.vx === undefined) {
                      pt.vx = 5 * (Math.random() - 0.5);
                      pt.vy = 5 * (Math.random() - 0.5);
                  }

                  pt.x += pt.vx;
                  pt.y += pt.vy;

                  if (pt.x < 0 || pt.x > 400) { pt.vx = -pt.vx; }
                  if (pt.y < 0 || pt.y > 400) { pt.vy = -pt.vy; }
              }
              context.clearRect(0, 0, 800, 400);
              drawPolyWithMargin(p, 10);
          }, 50);
      }
  });

<html>
    <head>
            <script type="text/javascript" src="http://ajax.googleapis.com/ajax/libs/jquery/1.4.2/jquery.min.js"></script>
    </head>
    <body>
        <canvas id="canvas" width="400" height="400"></canvas>  
    </body>
</html>

sample code disclaimers:

• the sample sacrifices some efficiency for the sake of clarity. In your code, you may want to compute each edge's expanded parallel just once, and not twice as in here

• the canvas's y coordinate grows downward, which inverts the CW/CCW logic. Things keep on working though as we just need to turn the outward normals in a direction opposite to the polygon's -- and both get flipped.

---

For each line segment of the original, find the midpoint m and (unit length) outward normal u of the segment. The corresponding segment of the expanded polygon will then lie on the line through m+n*u (where you want to expand the original by n) with normal u. To find the vertices of the expanded polygon you then need to find the intersection of pairs of successive lines.

---

If the polygon is centered on the origin simply multiply each of the points by a common scaling factor.

If the polygon is not centered on the origin then first translate so the center is on the origin, scale, and then translate it back to where it was.

After your comment

It seems you want all points to be moved the same distance away from the origin. You can do this for each point by getting the normalised vector to this point. multiplying this by your 'expand constant' and adding the resulting vector back onto the original point.

n.b. You will still have to translate-modify-translate if the center is not also the origin for this solution.

---

Let the points of the polygon be A1, B1, C1... Now you have lines from A1 to B1, then from B1 to C1... We want to compute points A2, B2, C2 of the polygon P2.

If you bisect angle, for example A1 B1 C1, you will have a line which goes in the direction you want. Now you can find a point B2 on it which is the appropriate distance from B1 on bisector line. Repeat this for all points of the polygon P1.

---

Look at straight skeletons. As has been implied here there are a number of tricky issues with non convex polygons that you have been mecifully spared!