By Stécy


2009-07-22 14:24:33 8 Comments

Having a list of points, how do I find if they are in clockwise order?

For example:

point[0] = (5,0)
point[1] = (6,4)
point[2] = (4,5)
point[3] = (1,5)
point[4] = (1,0)

would say that it is anti-clockwise (or counter-clockwise, for some people).

23 comments

@nbro 2020-05-24 19:41:36

Here's a simple Python 3 implementation based on this answer (which, in turn, is based on the solution proposed in the accepted answer)

def is_clockwise(points):
    # points is your list (or array) of 2d points.
    assert len(points) > 0
    s = 0.0
    for p1, p2 in zip(points, points[1:] + [points[0]]):
        s += (p2[0] - p1[0]) * (p2[1] + p1[1])
    return s > 0.0

@Olivier Jacot-Descombes 2013-08-27 18:28:55

Here is a simple C# implementation of the algorithm based on this answer.

Let's assume that we have a Vector type having X and Y properties of type double.

public bool IsClockwise(IList<Vector> vertices)
{
    double sum = 0.0;
    for (int i = 0; i < vertices.Count; i++) {
        Vector v1 = vertices[i];
        Vector v2 = vertices[(i + 1) % vertices.Count];
        sum += (v2.X - v1.X) * (v2.Y + v1.Y);
    }
    return sum > 0.0;
}

% is the modulo or remainder operator performing the modulo operation which (according to Wikipedia) finds the remainder after division of one number by another.

@ToolmakerSteve 2019-02-06 03:50:50

C# code to implement lhf's answer:

// https://en.wikipedia.org/wiki/Curve_orientation#Orientation_of_a_simple_polygon
public static WindingOrder DetermineWindingOrder(IList<Vector2> vertices)
{
    int nVerts = vertices.Count;
    // If vertices duplicates first as last to represent closed polygon,
    // skip last.
    Vector2 lastV = vertices[nVerts - 1];
    if (lastV.Equals(vertices[0]))
        nVerts -= 1;
    int iMinVertex = FindCornerVertex(vertices);
    // Orientation matrix:
    //     [ 1  xa  ya ]
    // O = | 1  xb  yb |
    //     [ 1  xc  yc ]
    Vector2 a = vertices[WrapAt(iMinVertex - 1, nVerts)];
    Vector2 b = vertices[iMinVertex];
    Vector2 c = vertices[WrapAt(iMinVertex + 1, nVerts)];
    // determinant(O) = (xb*yc + xa*yb + ya*xc) - (ya*xb + yb*xc + xa*yc)
    double detOrient = (b.X * c.Y + a.X * b.Y + a.Y * c.X) - (a.Y * b.X + b.Y * c.X + a.X * c.Y);

    // TBD: check for "==0", in which case is not defined?
    // Can that happen?  Do we need to check other vertices / eliminate duplicate vertices?
    WindingOrder result = detOrient > 0
            ? WindingOrder.Clockwise
            : WindingOrder.CounterClockwise;
    return result;
}

public enum WindingOrder
{
    Clockwise,
    CounterClockwise
}

// Find vertex along one edge of bounding box.
// In this case, we find smallest y; in case of tie also smallest x.
private static int FindCornerVertex(IList<Vector2> vertices)
{
    int iMinVertex = -1;
    float minY = float.MaxValue;
    float minXAtMinY = float.MaxValue;
    for (int i = 0; i < vertices.Count; i++)
    {
        Vector2 vert = vertices[i];
        float y = vert.Y;
        if (y > minY)
            continue;
        if (y == minY)
            if (vert.X >= minXAtMinY)
                continue;

        // Minimum so far.
        iMinVertex = i;
        minY = y;
        minXAtMinY = vert.X;
    }

    return iMinVertex;
}

// Return value in (0..n-1).
// Works for i in (-n..+infinity).
// If need to allow more negative values, need more complex formula.
private static int WrapAt(int i, int n)
{
    // "+n": Moves (-n..) up to (0..).
    return (i + n) % n;
}

@Warwick Allison 2019-02-09 00:37:00

This appears to be for down-is-positive Y coordinates. Flip CW/CCW for standard coordinates.

@VectorVortec 2019-02-03 11:18:31

While these answers are correct, they are more mathematically intense than necessary. Assume map coordinates, where the most north point is the highest point on the map. Find the most north point, and if 2 points tie, it is the most north then the most east (this is the point that lhf uses in his answer). In your points,

point[0] = (5,0)

point[1] = (6,4)

point[2] = (4,5)

point[3] = (1,5)

point[4] = (1,0)

If we assume that P2 is the most north then east point either the previous or next point determine clockwise, CW, or CCW. Since the most north point is on the north face, if the P1 (previous) to P2 moves east the direction is CW. In this case, it moves west, so the direction is CCW as the accepted answer says. If the previous point has no horizontal movement, then the same system applies to the next point, P3. If P3 is west of P2, it is, then the movement is CCW. If the P2 to P3 movement is east, it's west in this case, the movement is CW. Assume that nte, P2 in your data, is the most north then east point and the prv is the previous point, P1 in your data, and nxt is the next point, P3 in your data, and [0] is horizontal or east/west where west is less than east, and [1] is vertical.

if (nte[0] >= prv[0] && nxt[0] >= nte[0]) return(CW);
if (nte[0] <= prv[0] && nxt[0] <= nte[0]) return(CCW);
// Okay, it's not easy-peasy, so now, do the math
if (nte[0] * nxt[1] - nte[1] * nxt[0] - prv[0] * (nxt[1] - crt[1]) + prv[1] * (nxt[0] - nte[0]) >= 0) return(CCW); // For quadrant 3 return(CW)
return(CW) // For quadrant 3 return (CCW)

@ToolmakerSteve 2019-02-05 21:55:50

IMHO, it would be safer to stick to the fundamental math shown in lhf's answer - thank you for mentioning him. The challenge in reducing it to quadrants is that its a fair amount of work to prove that your formula is correct in all cases. Did you correctly calculate "more west"? In a concave polygon where both [1] and [3] are "west and south" of [2]? Did you correctly handle different lengths of [1] and [3] in that situation? I have no idea, whereas if I directly compute that angle (or its determinant), I'm using well-known formulas.

@VectorVortec 2019-02-06 01:42:57

@ToolmakerSteve the if statements always work if the 3 points are convex. The if statements will return, then you'll get the right answer. The if statements will not return, if the shape is concave and extreme. That's when you have to do the math. Most images have one quadrant, so that part is easy. More than 99% of my subroutine calls are handled by the if statements.

@ToolmakerSteve 2019-02-06 03:37:34

That doesn't address my concern. What is that formula? Is it the orientation determinant as given in the wiki link from lhf's answer? If so, then say so. Explain that what you are doing is doing quick checks that handle most cases, to avoid the standard math. If that is so, then your answer now makes sense to me. (Minor nit: would be easier to read if you used .x and .y of a struct, instead of [0] and [1]. I didn't know what your code was saying, first time I glanced at it.)

@ToolmakerSteve 2019-02-06 04:07:11

Since I did not have confidence in your approach, I implemented lhf's approach; formula from his link. Slow part is finding appropriate vertex - O(N) search. Once found, the determinant is an O(1) operation, using 6 multiplies with 5 adds. That last part is what you've optimized; but you've done so by adding additional if-tests. I can't personally justify taking a non-standard approach - would need to verify each step is correct - But thank you for an interesting analysis of quadrants!

@lhf 2009-07-24 21:30:21

Find the vertex with smallest y (and largest x if there are ties). Let the vertex be A and the previous vertex in the list be B and the next vertex in the list be C. Now compute the sign of the cross product of AB and AC.


References:

@M Katz 2015-06-02 06:30:01

This is also explained in en.wikipedia.org/wiki/Curve_orientation. The point is that the the found point must be on the convex hull, and it's only necessary to look locally at a single point on the convex hull (and its immediate neighbors) to determine the orientation of the whole polygon.

@Cecil Curry 2017-06-03 07:39:12

Shocked and awed this hasn't received more upvotes. For simple polygons (which is most polygons in some fields), this answer yields an O(1) solution. All other answers yield O(n) solutions for n the number of polygon points. For even deeper optimizations, see the Practical Considerations subsection of Wikipedia's fantastic Curve orientation article.

@Cecil Curry 2017-06-03 07:49:29

Clarification: this solution is O(1) only if either (A) this polygon is convex (in which case any arbitrary vertex resides on the convex hull and hence suffices) or (B) you already know the vertex with the smallest Y coordinate. If this is not the case (i.e., this polygon is non-convex and you don't know anything about it), an O(n) search is required. Since no summation is required, however, this is still dramatically faster than any other solution for simple polygons.

@ToolmakerSteve 2019-02-06 03:51:32

@LarsH 2020-05-02 15:16:28

@CecilCurry I think your 2nd comment explains why this hasn't received more upvotes. It yields wrong answers in certain scenarios, without any mention of those limitations.

@dez 2018-06-29 04:25:19

Solution for R to determine direction and reverse if clockwise (found it necessary for owin objects):

coords <- cbind(x = c(5,6,4,1,1),y = c(0,4,5,5,0))
a <- numeric()
for (i in 1:dim(coords)[1]){
  #print(i)
  q <- i + 1
  if (i == (dim(coords)[1])) q <- 1
  out <- ((coords[q,1]) - (coords[i,1])) * ((coords[q,2]) + (coords[i,2]))
  a[q] <- out
  rm(q,out)
} #end i loop

rm(i)

a <- sum(a) #-ve is anti-clockwise

b <- cbind(x = rev(coords[,1]), y = rev(coords[,2]))

if (a>0) coords <- b #reverses coords if polygon not traced in anti-clockwise direction

@Ind 2018-05-03 11:58:46

Another solution for this;

const isClockwise = (vertices=[]) => {
    const len = vertices.length;
    const sum = vertices.map(({x, y}, index) => {
        let nextIndex = index + 1;
        if (nextIndex === len) nextIndex = 0;

        return {
            x1: x,
            x2: vertices[nextIndex].x,
            y1: x,
            y2: vertices[nextIndex].x
        }
    }).map(({ x1, x2, y1, y2}) => ((x2 - x1) * (y1 + y2))).reduce((a, b) => a + b);

    if (sum > -1) return true;
    if (sum < 0) return false;
}

Take all the vertices as an array like this;

const vertices = [{x: 5, y: 0}, {x: 6, y: 4}, {x: 4, y: 5}, {x: 1, y: 5}, {x: 1, y: 0}];
isClockwise(vertices);

@Toby Evetts 2018-02-28 17:26:08

Here's swift 3.0 solution based on answers above:

    for (i, point) in allPoints.enumerated() {
        let nextPoint = i == allPoints.count - 1 ? allPoints[0] : allPoints[i+1]
        signedArea += (point.x * nextPoint.y - nextPoint.x * point.y)
    }

    let clockwise  = signedArea < 0

@Charles Bretana 2009-07-22 18:23:26

The cross product measures the degree of perpendicular-ness of two vectors. Imagine that each edge of your polygon is a vector in the x-y plane of a three-dimensional (3-D) xyz space. Then the cross product of two successive edges is a vector in the z-direction, (positive z-direction if the second segment is clockwise, minus z-direction if it's counter-clockwise). The magnitude of this vector is proportional to the sine of the angle between the two original edges, so it reaches a maximum when they are perpendicular, and tapers off to disappear when the edges are collinear (parallel).

So, for each vertex (point) of the polygon, calculate the cross-product magnitude of the two adjoining edges:

Using your data:
point[0] = (5, 0)
point[1] = (6, 4)
point[2] = (4, 5)
point[3] = (1, 5)
point[4] = (1, 0)

So Label the edges consecutively as
edgeA is the segment from point0 to point1 and
edgeB between point1 to point2
...
edgeE is between point4 and point0.

Then Vertex A (point0) is between
edgeE [From point4 to point0]
edgeA [From point0 to `point1'

These two edges are themselves vectors, whose x and y coordinates can be determined by subtracting the coordinates of their start and end points:

edgeE = point0 - point4 = (1, 0) - (5, 0) = (-4, 0) and
edgeA = point1 - point0 = (6, 4) - (1, 0) = (5, 4) and

And the cross product of these two adjoining edges is calculated using the determinant of the following matrix, which is constructed by putting the coordinates of the two vectors below the symbols representing the three coordinate axis (i, j, & k). The third (zero)-valued coordinate is there because the cross product concept is a 3-D construct, and so we extend these 2-D vectors into 3-D in order to apply the cross-product:

 i    j    k 
-4    0    0
 1    4    0    

Given that all cross-products produce a vector perpendicular to the plane of two vectors being multiplied, the determinant of the matrix above only has a k, (or z-axis) component.
The formula for calculating the magnitude of the k or z-axis component is
a1*b2 - a2*b1 = -4* 4 - 0* 1 = -16

The magnitude of this value (-16), is a measure of the sine of the angle between the 2 original vectors, multiplied by the product of the magnitudes of the 2 vectors.
Actually, another formula for its value is
A X B (Cross Product) = |A| * |B| * sin(AB).

So, to get back to just a measure of the angle you need to divide this value, (-16), by the product of the magnitudes of the two vectors.

|A| * |B| = 4 * Sqrt(17) = 16.4924...

So the measure of sin(AB) = -16 / 16.4924 = -.97014...

This is a measure of whether the next segment after the vertex has bent to the left or right, and by how much. There is no need to take arc-sine. All we will care about is its magnitude, and of course its sign (positive or negative)!

Do this for each of the other 4 points around the closed path, and add up the values from this calculation at each vertex..

If final sum is positive, you went clockwise, negative, counterclockwise.

@Charles Bretana 2013-04-06 14:29:57

Actually, this solution is a different solution than the accepted solution. Whether they are equivalent or not is a question I am investigating, but I suspect they are not... The accepted answer calculates the area of the polygon, by taking the difference between the area under the top edge of the polygon and the area under the bottom edge of the polygon. One will be negative (the one where you are traversing from left to right), and the other will be negative. When traversing clockwise, The upper edge is traversed left to right and is larger, so the total is positive.

@Charles Bretana 2013-04-06 14:30:35

My solution measures the sum of sines of the changes in edge angles at each vertex. This will be positive when traversing clockwise and negative when traversing counter clockwise.

@agentp 2013-05-30 16:54:12

It seems with this approach you DO need to take the arcsin, unless you assume convexity (in which case you need only check one vertex)

@Charles Bretana 2015-05-05 23:53:52

No, you do not need to take the arcsine. You do need to divide by the magnitudes of the two vector, to eliminate scaling, (where one very large pair od segments that bends one way a small angle contributes more tot the sum than many smaller segments that bend the other way by a larger angle. But the arcsine just converts one measure of an angle to another.

@marsbear 2018-01-19 20:24:08

Why all the expensive math functions? Wouldn't it be faster to only calculate the z-part of a 3D cross product by doing a1 * b2 - a2 * b1 (as the x- and y-parts will be zero anyway) and checking the sign of the result? No trigonometric functions are needed and no sqrt.

@Charles Bretana 2018-01-19 21:51:00

@marsbear, That will only work if all vectors are the same magnitude (or, equivalently, have first been normalized). If they are of different lengths, then you will get the wrong answer., because a clockwise angle between two long segments of the polygon will contribute much more to the overall sum than a counterclockwise angle between two very small segments.

@marsbear 2018-01-24 10:25:24

But shouldn't it be enough to "add up" the signs of the cross products?

@Charles Bretana 2018-01-24 15:52:45

No, all that would do is "count" the number of turns to the right minus the number of turns to the left. If you start pointing East, and take ninety one degree turns to the right, (totaling 90 degrees of right), each an inch long, you will be pointing south, and then take five 90 degree turns to the left, you will be pointing east again, and with appropriate adjustment of the length of each segment, (make the lefts longer) reconnect to the start point. The overall closed curve is clearly counterclockwise, but just counting the signs would count 90 to the right and only 5 to the left.

@nbro 2018-02-24 23:46:35

I have reformulated your answer (to make it more precise mathematically) until your sentence "The magnitude of this value is a measure of the sine of the angle between the 2 original vectors". You should edit that sentence and the remaining ones to take into account my changes. The reason I didn't edit that sentence and the remaining ones is that I am not fully understanding what you mean those.

@Luke Hutchison 2018-06-08 09:07:54

You DO need to take the arcsin. Try it on a bunch of random non-convex polygons, and you will find the test will fail for some polygons if you don't take the arcsin.

@Charles Bretana 2018-06-12 14:48:53

No, You don't. Whatever test you are running, you are running it wrong. The value before taking the arcsine IS the sin of the original angle of rotation, converting it from a value between -1 and +1 to an actual angle between -Pi and +Pi does nothing to change the answer. And as long as you divide by the product of the magnitudes of the two vectors, the magnitudes will be normalized, and their sum will indicate the overall direction of curvature.

@ToolmakerSteve 2019-02-05 20:03:08

@CharlesBretana - while I have not run Luke's test, I believe he is correct. That is the nature of summing combined with a nonlinear scale [without arcsin vs. with arcsin]. Consider what marsbear suggested, that you correctly rejected. He suggested that you "just count", and you pointed out that a handful of large values could outweigh a large number of small values. Now consider arcsin of each value vs not. Isn't it still the case that failing to take arcsin gives incorrect weight to each value, therefore has same flaw (though much less so)?

@MSS 2014-11-18 10:56:17

This is the implemented function for OpenLayers 2. The condition for having a clockwise polygon is area < 0, it confirmed by this reference.

function IsClockwise(feature)
{
    if(feature.geometry == null)
        return -1;

    var vertices = feature.geometry.getVertices();
    var area = 0;

    for (var i = 0; i < (vertices.length); i++) {
        j = (i + 1) % vertices.length;

        area += vertices[i].x * vertices[j].y;
        area -= vertices[j].x * vertices[i].y;
        // console.log(area);
    }

    return (area < 0);
}

@MSS 2018-02-26 07:29:27

Openlayers is javascript based map management library like googlemaps and it's wrote and used in openlayers 2.

@nbro 2018-02-26 10:25:42

Can you explain a little bit what your code does, and why you're doing it?

@allez l'OM 2019-12-03 14:31:08

@nbro this code implements the lhf answer. It is easy to keep the non OpenLayer part in a pure javascript function by having vertices directly as parameter. It works well, and could be adapted to the case of multiPolygon.

@Ian 2014-05-01 17:33:38

As also explained in this Wikipedia article Curve orientation, given 3 points p, q and r on the plane (i.e. with x and y coordinates), you can calculate the sign of the following determinant

enter image description here

If the determinant is negative (i.e. Orient(p, q, r) < 0), then the polygon is oriented clockwise (CW). If the determinant is positive (i.e. Orient(p, q, r) > 0), the polygon is oriented counterclockwise (CCW). The determinant is zero (i.e. Orient(p, q, r) == 0) if points p, q and r are collinear.

In the formula above, we prepend the ones in front of the coordinates of p, q and r because we are using homogeneous coordinates.

@nbro 2018-02-25 11:34:39

@tibetty Can you explain why this method wouldn't work in many situations if the polygon is concave?

@tibetty 2018-02-28 01:20:19

Please look at the last table in the wiki item reference in your post. It's easy for me to give a false example but hard to prove it.

@tibetty 2018-02-28 01:22:49

Please look at the last table in the wiki item reference in your post. It's easy for me to give a false example but hard to prove it.

@ToolmakerSteve 2019-02-05 21:07:12

@tibetty is correct. You can't simply take any three points along the polygon; you might be in either a convex or concave region of that polygon. Reading wiki carefully, one must take three points along the convex hull that encloses the polygon. From "practical considerations": "One does not need to construct the convex hull of a polygon to find a suitable vertex. A common choice is the vertex of the polygon with the smallest X-coordinate. If there are several of them, the one with the smallest Y-coordinate is picked. It is guaranteed to be [a] vertex of the convex hull of the polygon."

@ToolmakerSteve 2019-02-05 21:11:13

Hence lhf's earlier answer, which is similar, and references the same wiki article, but specifies such a point. [Apparently it doesn't matter whether one takes smallest or largest, x or y, as long as one avoids being in the middle; effectively one is working from one edge of the bounding box around the polygon, to guarantee in a concave region.]

@Marjan Moderc 2017-09-08 08:50:46

After testing several unreliable implementations, the algorithm that provided satisfactory results regarding the CW/CCW orientation out of the box was the one, posted by OP in this thread (shoelace_formula_3).

As always, a positive number represents a CW orientation, whereas a negative number CCW.

@Venkata Goli 2018-01-21 22:04:05

A much computationally simpler method, if you already know a point inside the polygon:

  1. Choose any line segment from the original polygon, points and their coordinates in that order.

  2. Add a known "inside" point, and form a triangle.

  3. Calculate CW or CCW as suggested here with those three points.

@ToolmakerSteve 2019-02-05 21:22:59

Maybe this works if the polygon is entirely convex. It definitely is not reliable if there are any concave regions - its easy to pick a point that is on the "wrong" side of one of the edges of the cave, then connect it to that edge. Will get wrong answer.

@Venkata Goli 2019-02-06 22:28:14

It works even if the polygon is concave. The point needs to be inside that concave polygon. However I am not sure about complex polygon (did not test.)

@ToolmakerSteve 2019-02-07 18:29:28

"It works even if the polygon is concave." - Counterexample: poly (0,0), (1,1), (0,2), (2,2), (2,0). Line segment (1,1), (0, 2). If you pick an interior point within (1,1), (0,2), (1,2) to form triangle -> (1,1), (0,2), (0.5,1.5)), you get opposite winding than if you pick an interior point within (0,0), (1,1), (1,0) > (1,1), (0,2),(0.5,0.5). Those are both interior to the original polygon, yet have opposite windings. Therefore, one of them gives the wrong answer.

@ToolmakerSteve 2019-02-07 18:42:50

In general, if a polygon has any concave region, pick a segment in the concave region. Because it is concave, you can find two "interior" points that are on opposite sides of that line. Because they are on opposite sides of that line, the triangles formed have opposite windings. End of proof.

@Frederick 2013-10-17 14:10:42

If you use Matlab, the function ispolycw returns true if the polygon vertices are in clockwise order.

@Gianni Spear 2013-03-06 15:25:29

This is my solution using the explanations in the other answers:

def segments(poly):
    """A sequence of (x,y) numeric coordinates pairs """
    return zip(poly, poly[1:] + [poly[0]])

def check_clockwise(poly):
    clockwise = False
    if (sum(x0*y1 - x1*y0 for ((x0, y0), (x1, y1)) in segments(poly))) < 0:
        clockwise = not clockwise
    return clockwise

poly = [(2,2),(6,2),(6,6),(2,6)]
check_clockwise(poly)
False

poly = [(2, 6), (6, 6), (6, 2), (2, 2)]
check_clockwise(poly)
True

@nbro 2018-02-25 11:15:13

Can you specify which other answers exactly this answer is based on?

@mpen 2017-10-06 20:36:03

An implementation of Sean's answer in JavaScript:

function calcArea(poly) {
    if(!poly || poly.length < 3) return null;
    let end = poly.length - 1;
    let sum = poly[end][0]*poly[0][1] - poly[0][0]*poly[end][1];
    for(let i=0; i<end; ++i) {
        const n=i+1;
        sum += poly[i][0]*poly[n][1] - poly[n][0]*poly[i][1];
    }
    return sum;
}

function isClockwise(poly) {
    return calcArea(poly) > 0;
}

let poly = [[352,168],[305,208],[312,256],[366,287],[434,248],[416,186]];

console.log(isClockwise(poly));

let poly2 = [[618,186],[650,170],[701,179],[716,207],[708,247],[666,259],[637,246],[615,219]];

console.log(isClockwise(poly2));

Pretty sure this is right. It seems to be working :-)

Those polygons look like this, if you're wondering:

@bradgonesurfing 2016-06-24 07:21:13

My C# / LINQ solution is based on the cross product advice of @charlesbretana is below. You can specify a reference normal for the winding. It should work as long as the curve is mostly in the plane defined by the up vector.

using System.Collections.Generic;
using System.Linq;
using System.Numerics;

namespace SolidworksAddinFramework.Geometry
{
    public static class PlanePolygon
    {
        /// <summary>
        /// Assumes that polygon is closed, ie first and last points are the same
        /// </summary>
       public static bool Orientation
           (this IEnumerable<Vector3> polygon, Vector3 up)
        {
            var sum = polygon
                .Buffer(2, 1) // from Interactive Extensions Nuget Pkg
                .Where(b => b.Count == 2)
                .Aggregate
                  ( Vector3.Zero
                  , (p, b) => p + Vector3.Cross(b[0], b[1])
                                  /b[0].Length()/b[1].Length());

            return Vector3.Dot(up, sum) > 0;

        } 

    }
}

with a unit test

namespace SolidworksAddinFramework.Spec.Geometry
{
    public class PlanePolygonSpec
    {
        [Fact]
        public void OrientationShouldWork()
        {

            var points = Sequences.LinSpace(0, Math.PI*2, 100)
                .Select(t => new Vector3((float) Math.Cos(t), (float) Math.Sin(t), 0))
                .ToList();

            points.Orientation(Vector3.UnitZ).Should().BeTrue();
            points.Reverse();
            points.Orientation(Vector3.UnitZ).Should().BeFalse();



        } 
    }
}

@Beta 2009-07-22 15:02:22

Some of the suggested methods will fail in the case of a non-convex polygon, such as a crescent. Here's a simple one that will work with non-convex polygons (it'll even work with a self-intersecting polygon like a figure-eight, telling you whether it's mostly clockwise).

Sum over the edges, (x2 − x1)(y2 + y1). If the result is positive the curve is clockwise, if it's negative the curve is counter-clockwise. (The result is twice the enclosed area, with a +/- convention.)

point[0] = (5,0)   edge[0]: (6-5)(4+0) =   4
point[1] = (6,4)   edge[1]: (4-6)(5+4) = -18
point[2] = (4,5)   edge[2]: (1-4)(5+5) = -30
point[3] = (1,5)   edge[3]: (1-1)(0+5) =   0
point[4] = (1,0)   edge[4]: (5-1)(0+0) =   0
                                         ---
                                         -44  counter-clockwise

@Beta 2009-07-27 14:20:26

It's calculus applied to a simple case. (I don't have the skill to post graphics.) The area under a line segment equals its average height (y2+y1)/2 times its horizontal length (x2-x1). Notice the sign convention in x. Try this with some triangles and you'll soon see how it works.

@Beta 2009-08-31 15:44:30

I'd suggest proving it first for a triangle, then showing that any polygon can be split into triangles in a way that doesn't change its "clockwiseness".

@Mark Ingram 2009-11-09 17:19:42

Is there a particular name to this method?

@Beta 2009-11-09 19:17:25

Now that I think of it, the method of calculating area using a counterclockwise curve can be derived from Green's Theorem. I saw the method (without derivation) years ago-- I may have come up with idea of using it as a clockwiseness test myself, but it's a pretty obvious jump.

@the Ritz 2013-03-12 13:47:06

Is there a similar elegant solution for polygons in 3D space?

@the Ritz 2013-03-12 13:54:55

@Michael 2013-08-27 20:39:41

It also generalizes nicely to 3-D polytopes: for each face compute the signed area in x-y plane as described in this answer, multiply by the mean value of z, and sum it up for all faces. The absolute value of the result will be the polytope's volume.

@LarsH 2013-10-11 20:49:22

A minor caveat: this answer assumes a normal Cartesian coordinate system. The reason that's worth mentioning is that some common contexts, like HTML5 canvas, use an inverted Y-axis. Then the rule has to be flipped: if the area is negative, the curve is clockwise.

@Mr.Qbs 2014-12-13 16:49:57

This method is wrong! I implemented it ad is works in 98% cases. But not always. What worked for me (in every case) is this link.

@Beta 2014-12-14 00:21:10

@Mr.Qbs: Please give us the simplest counterexample you know.

@Mr.Qbs 2014-12-14 11:38:15

I checked my implementation to show you that it's properly done and to give you counterexample. Then I saw your coment "I'd suggest proving it for a triangle" so I tried it and it works! When i have 3 points (2 edges of my big polygon) I need to calculate 3 edges to get properly result (ab, bc, ca - imaginary edge). But when I try it only for 2 "real" edges (that realy belongs to my polygon) result sometimes is bad (its hard to give example - i compute 3d objects and I can see the bad arc quickly on modified model - sometimes arc is made in bad direction).

@Beta 2014-12-14 14:03:30

@Mr.Qbs: So my method works, but if you skip a vital part, then it doesn't work. This is not news.

@Olivier Jacot-Descombes 2015-07-20 14:45:09

@Mr.Qbs: You always have to link the last point to the first one. If you have N points numbered from 0 to N-1, then you must calculate: Sum( (x[(i+1) mod N] - x[i]) * (y[i] + y[(i+1) mod N]) ) for i = 0 to N-1. I.e., must must take the index Modulo N (N ≡ 0) The formula works only for closed polygons. Polygons have no imaginary edges.

@Dan M. 2016-08-02 17:58:55

Can the result be zero? And if it can, then what orientation is it? (It looks like some special cease).

@Beta 2016-08-03 02:23:24

@DanM. If the result is zero, that means that the positive and negative areas cancel out, as in a figure-eight.

@Dan M. 2016-08-13 22:30:03

Btw, It seems you only need to calculate sum of x<sub>i + 1</sub>*y<sub>i</sub> - x<sub>i</sub>*y<sub>i + 1</sub> to determine the answer.

@Beta 2016-08-14 19:54:49

@DanM.: True, but notice that my formula has only one multiplication, while yours has two.

@Matt Fiocca 2017-01-12 02:58:53

I'm trying this out with a clockwise diamond, and getting a negative (counter-clockwise) sum result, and a positive sum for counter-clockwise. CLOCKWISE DIAMOND: point[0] = 0,5 (5-0)(0+5) = 25, point[1] = 5,0 (10-5)(5+0) = 25, point[2] = 10,5 (5-10)(5+10) = -75, point[3] = 5,10 (0-5)(10+5) = -75 == -100. COUNTER DIAMOND: point[0] = 0,5 (5-0)(5+10) = 75, point[1] = 5,10 (10-5)(10+5) = 75, point[2] = 10,5 (5-10)(5+0) = -25, point[3] = 5,0 (0-5)(5+0) = -25 == +100 Is this a special case?

@Beta 2017-01-12 03:09:47

@MattFiocca: You are using the words "clockwise" and "counter-clockwise" in the opposite sense of everyone else.

@Matt Fiocca 2017-01-12 05:28:31

@Beta: just so i'm not going crazy, is everyone else using clockwise as counter-clockface? Here's a link to my diagram with calculations. s3-us-west-2.amazonaws.com/mattfiocca/polydirection.png.

@Beta 2017-01-12 23:41:12

@MattFiocca: The convention in mathematics is that the Y coordinate increases upward; your Y coordinates increase downward.

@Matt Fiocca 2017-01-13 03:33:49

@Beta: ah thanks! I knew I was missing something vital, having come from the backward world of web and apps.

@Sreeraj 2017-03-17 07:03:13

This works only with all positive co-ordinates. I tried with Polygons in X,-Z Plane (all Z values are negative in this case.). Only way I could get this logic to work in this case is by translating the polygons to X,+Z plane. i.e, by offsetting the Z values of all co-ordinates with a constant such that all values are positive. Please correct me if am wrong or if there is any simpler alternative.

@Beta 2017-03-18 22:08:36

@Sreeraj: You are mistaken, but a comment thread is not the place to discuss it. If I have time I'll make a chat room.

@David Zorychta 2017-05-17 01:16:03

This blog.element84.com/polygon-winding.html explains in simple english why this solution works.

@RunninglVlan 2017-11-26 16:14:04

Alternatively it can be calculated as sum of (point.x + nextPoint.x) * (point.y - nextPoint.y).

@Sandy Gifford 2018-04-16 11:36:36

Would this work for lat/lon coordinates assuming you're not near/crossing a pole, and all polygons that cross 180 are defined with numbers greater than 180?

@Ans 2018-11-27 14:51:18

@Macmee: the link given by you is broken, archive: web.archive.org/web/20171212175910/blog.element84.com/…

@Brian Carcich 2018-12-04 16:57:07

Ah, with a clever approach that cancels the extra XnYn terms over the sum, this is a cross product. Sweet!

@Chanwoo Park 2020-04-23 05:50:05

@LarsH unless it's Y coordinate is inverted resulting in a negative Y value, wouldn't positive sum of edges still means the polygon is clockwise?

@LarsH 2020-04-23 14:21:46

@ChanwooPark: Good question. I can see why that might appear to be the case, since y1 + y2 wouldn't change sign just from inverting the direction of the y axis. But I don't think so. Try it for a simple triangle. (I don't have time to try it now.)

@LarsH 2020-04-23 14:28:41

Ok I tried it with the example given in the answer. If you flip the y-axis to go from 5 to 0 instead of 0 to 5, you end up with 45 instead of -44. So no, a positive sum doesn't mean the polygon is clockwise, when the y-axis is flipped.

@LarsH 2020-04-23 15:28:52

@ChanwooPark Think of how the algorithm works: if edges that go leftward (x2 < x1) have higher y-values overall than edges that go rightward (x2 > x1), the negative terms will outweigh the positive, and in a Cartesian coordinate system this means the polygon goes counterclockwise.

@LarsH 2020-04-23 15:31:04

Whereas in an inverted y-coordinate system, edges that would have had higher y-values than the others now have lower y-values. So for a counterclockwise polygon, the leftward edges now have lower y-values, and the positive terms of the rightward edges will outweigh them, producing a positive result.

@Chanwoo Park 2020-04-24 07:49:53

@LarsH Thx for the answer it helped me alot

@LarsH 2020-04-24 15:24:48

@ChanwooPark Thx for the question. I hadn't really understood the algorithm until you provoked me to think it thru.

@Sean the Bean 2012-04-24 13:19:30

I guess this is a pretty old question, but I'm going to throw out another solution anyway, because it's straightforward and not mathematically intensive - it just uses basic algebra. Calculate the signed area of the polygon. If it's negative the points are in clockwise order, if it's positive they are counterclockwise. (This is very similar to Beta's solution.)

Calculate the signed area: A = 1/2 * (x1*y2 - x2*y1 + x2*y3 - x3*y2 + ... + xn*y1 - x1*yn)

Or in pseudo-code:

signedArea = 0
for each point in points:
    x1 = point[0]
    y1 = point[1]
    if point is last point
        x2 = firstPoint[0]
        y2 = firstPoint[1]
    else
        x2 = nextPoint[0]
        y2 = nextPoint[1]
    end if

    signedArea += (x1 * y2 - x2 * y1)
end for
return signedArea / 2

Note that if you are only checking the ordering, you don't need to bother dividing by 2.

Sources: http://mathworld.wolfram.com/PolygonArea.html

@Michael Eric Oberlin 2015-05-25 20:25:38

Was that a typo in your signed area formula above? It ends with "xn*y1 - x1*yn"; when I believe it should be "x_n y_{n+1} - y_n x_{n-1}" (in LaTeX, at least). On the other hand, it's been ten years since I took any linear algebra classes.

@Sean the Bean 2015-05-26 17:38:43

Nope. If you check the source, you'll see that the formula does in fact reference the first point again in the last term (y1 and x1). (Sorry, I'm not very familiar with LaTeX, but I formatted the subscripts to make them more readable.)

@Eric Fortier 2019-01-16 02:03:29

I used this solution and it worked perfectly for my use. Note that if you can plan ahead and spare and extra two vectors in your array, you can get rid of the comparison (or %) by adding the first vector at the tail of the array. That way you simply loop over all the elements, except the last one (length-2 instead of length-1).

@ToolmakerSteve 2019-02-05 20:16:19

@EricFortier - FWIW, rather than resize a possibly large array, an efficient alternative is for each iteration to save its point as previousPoint for next iteration. Before starting loop, set previousPoint to array's last point. Trade off is extra local variable copy but fewer array accesses. And most importantly, don't have to touch the input array.

@ToolmakerSteve 2019-02-05 20:20:07

@MichaelEricOberlin - its necessary to close the polygon, by including the line segment from last point to first point. (A correct calculation will be the same, no matter which point starts the closed polygon.)

@Eric Fortier 2019-02-10 16:45:21

@ToolmakerSteve - Good point. Unless you know the size beforehand and can create the array with the extra element, your solution is fast and prevents resizing.

@daniel 2015-02-24 12:41:31

I think in order for some points to be given clockwise all edges need to be positive not only the sum of edges. If one edge is negative than at least 3 points are given counter-clockwise.

@ToolmakerSteve 2019-02-05 21:40:49

True, but you misunderstand the concept of a polygon's winding order (clockwise or counter-clockwise). In an entirely convex polygon, the angle at all points will be clockwise or all will be counter-clockwise [as in your first sentence]. In a polygon with concave region(s), the "caves" will be in the opposite direction, but the polygon as a whole still has a well-defined interior, and is considered clockwise or counter-clockwise accordingly. See en.wikipedia.org/wiki/…

@Steve Jansen 2013-12-20 05:35:34

For what it is worth, I used this mixin to calculate the winding order for Google Maps API v3 apps.

The code leverages the side effect of polygon areas: a clockwise winding order of vertexes yields a positive area, while a counter-clockwise winding order of the same vertexes produces the same area as a negative value. The code also uses a sort of private API in the Google Maps geometry library. I felt comfortable using it - use at your own risk.

Sample usage:

var myPolygon = new google.maps.Polygon({/*options*/});
var isCW = myPolygon.isPathClockwise();

Full example with unit tests @ http://jsfiddle.net/stevejansen/bq2ec/

/** Mixin to extend the behavior of the Google Maps JS API Polygon type
 *  to determine if a polygon path has clockwise of counter-clockwise winding order.
 *  
 *  Tested against v3.14 of the GMaps API.
 *
 *  @author  [email protected]
 *
 *  @license http://opensource.org/licenses/MIT
 *
 *  @version 1.0
 *
 *  @mixin
 *  
 *  @param {(number|Array|google.maps.MVCArray)} [path] - an optional polygon path; defaults to the first path of the polygon
 *  @returns {boolean} true if the path is clockwise; false if the path is counter-clockwise
 */
(function() {
  var category = 'google.maps.Polygon.isPathClockwise';
     // check that the GMaps API was already loaded
  if (null == google || null == google.maps || null == google.maps.Polygon) {
    console.error(category, 'Google Maps API not found');
    return;
  }
  if (typeof(google.maps.geometry.spherical.computeArea) !== 'function') {
    console.error(category, 'Google Maps geometry library not found');
    return;
  }

  if (typeof(google.maps.geometry.spherical.computeSignedArea) !== 'function') {
    console.error(category, 'Google Maps geometry library private function computeSignedArea() is missing; this may break this mixin');
  }

  function isPathClockwise(path) {
    var self = this,
        isCounterClockwise;

    if (null === path)
      throw new Error('Path is optional, but cannot be null');

    // default to the first path
    if (arguments.length === 0)
        path = self.getPath();

    // support for passing an index number to a path
    if (typeof(path) === 'number')
        path = self.getPaths().getAt(path);

    if (!path instanceof Array && !path instanceof google.maps.MVCArray)
      throw new Error('Path must be an Array or MVCArray');

    // negative polygon areas have counter-clockwise paths
    isCounterClockwise = (google.maps.geometry.spherical.computeSignedArea(path) < 0);

    return (!isCounterClockwise);
  }

  if (typeof(google.maps.Polygon.prototype.isPathClockwise) !== 'function') {
    google.maps.Polygon.prototype.isPathClockwise = isPathClockwise;
  }
})();

@Cameron Roberts 2018-09-27 20:16:40

Trying this I get exactly the opposite result, a polygon drawn in clockwise order yields a negative area, while one drawn counter clockwise yields positive. In either case, this snippet is still super useful 5yrs on, thank you.

@allez l'OM 2019-12-03 13:46:21

@CameronRoberts The norm (see IETF in particular for geoJson) is to follow the 'right-hand rule'. I guess that Google is complaining with. In that case outer ring must be counterclockwise (yielding positive area), and inner rings (holes) are winding clockwise (negative area to be removed from main area).

@ufukgun 2009-07-22 14:38:14

find the center of mass of these points.

suppose there are lines from this point to your points.

find the angle between two lines for line0 line1

than do it for line1 and line2

...

...

if this angle is monotonically increasing than it is counterclockwise ,

else if monotonically decreasing it is clockwise

else (it is not monotonical)

you cant decide, so it is not wise

@Vicky Chijwani 2012-10-09 04:34:08

by "center of mass" I think you mean "centroid"?

@ToolmakerSteve 2019-02-05 21:24:21

Probably works if polygon is entirely convex. But better to instead use an answer that will work for non-convex polygons.

@Steve Gilham 2009-07-22 14:31:08

Start at one of the vertices, and compute the angle subtended by each side.

The first and the last will be zero (so skip those); for the rest, the sine of the angle will be given by the cross product of the normalizations to unit length of (point[n]-point[0]) and (point[n-1]-point[0]).

If the sum of the values is positive, then your polygon is drawn in the anti-clockwise sense.

@ReaperUnreal 2009-07-22 14:44:54

Seeing as how the cross product basically boils down to a positive scaling factor times the sine of the angle, it's probably better to just do a cross product. It'll be faster and less complicated.

Related Questions

Sponsored Content

48 Answered Questions

36 Answered Questions

12 Answered Questions

5 Answered Questions

[SOLVED] Sort points in clockwise order?

35 Answered Questions

3 Answered Questions

[SOLVED] sorting a list of 3D points in clockwise order

1 Answered Questions

2 Answered Questions

[SOLVED] How to determine ring points' orientation correctly?

1 Answered Questions

[SOLVED] How to order 3 points counter clockwise around normal

Sponsored Content