ポイントのリストがある場合、それらが時計回りの順序であるかどうかを確認するにはどうすればよいですか?
例えば:
point[0] = (5,0)
point[1] = (6,4)
point[2] = (4,5)
point[3] = (1,5)
point[4] = (1,0)
それは反時計回り(または一部の人にとっては反時計回り)であると言うでしょう。
ポイントのリストがある場合、それらが時計回りの順序であるかどうかを確認するにはどうすればよいですか?
例えば:
point[0] = (5,0)
point[1] = (6,4)
point[2] = (4,5)
point[3] = (1,5)
point[4] = (1,0)
それは反時計回り(または一部の人にとっては反時計回り)であると言うでしょう。
提案された方法のいくつかは、三日月などの非凸多角形の場合に失敗します。これは、非凸多角形で機能する単純なものです(8の字のような自己交差する多角形でも機能し、ほとんど時計回りであるかどうかを示します)。
エッジを合計します(x 2 − x 1)(y 2 + y 1)。結果が正の場合、曲線は時計回りになり、負の場合、曲線は反時計回りになります。(結果は、+/-規則で、囲まれた領域の2倍になります。)
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
yが最小(同点の場合はxが最大)の頂点を見つけます。頂点A
をとし、リスト内の前の頂点をとし、リストB
内の次の頂点をとしますC
。次に、との交差積の符号を計算します。AB
AC
参照:
単純なポリゴンの方向を見つけるにはどうすればよいですか?よくある 質問:comp.graphics.algorithms。
ウィキペディアでの曲線の向き。
単純で数学的に集中的ではないため、別のソリューションを破棄します。これは、基本的な代数を使用するだけです。ポリゴンの符号付き領域を計算します。負の場合、ポイントは時計回りの順序であり、正の場合、ポイントは反時計回りです。(これはベータ版のソリューションと非常によく似ています。)
符号付き領域を計算します:A = 1/2 *(x 1 * y 2 -x 2 * y 1 + x 2 * y 3 -x 3 * y 2 + ... + x n * y 1 -x 1 * y n)
または擬似コードで:
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
順序を確認するだけの場合は、わざわざ2で割る必要はありません。
出典: http: //mathworld.wolfram.com/PolygonArea.html
外積は、2つのベクトルの垂直度を測定します。ポリゴンの各エッジが、3次元(3-D)xyz空間のxy平面内のベクトルであると想像してください。次に、2つの連続するエッジの外積は、z方向のベクトルです(2番目のセグメントが時計回りの場合は正のz方向、反時計回りの場合はz方向を引いたもの)。このベクトルの大きさは、元の2つのエッジ間の角度の正弦に比例するため、2つの元のエッジが垂直になると最大になり、エッジが同一線上(平行)になると先細りになって消えます。
したがって、ポリゴンの各頂点(ポイント)について、隣接する2つのエッジの外積の大きさを計算します。
Using your data:
point[0] = (5, 0)
point[1] = (6, 4)
point[2] = (4, 5)
point[3] = (1, 5)
point[4] = (1, 0)
したがって、 toとbetween to
...の間
edgeA
のセグメントがとの間であるように、エッジに連続してラベルを付けます。 point0
point1
edgeB
point1
point2
edgeE
point4
point0
次に、頂点A(point0
)は
edgeE
[From point4
to point0
]
edgeA
[From point0
to`point1'の間にあります
これらの2つのエッジはそれ自体がベクトルであり、そのx座標とy座標は、開始点と終了点の座標を引くことで決定できます。
edgeE
= point0
- point4
== および(1, 0) - (5, 0)
= - == および_ _(-4, 0)
edgeA
point1
point0
(6, 4) - (1, 0)
(5, 4)
そして、これら2つの隣接するエッジの外積は、次の行列式を使用して計算されます。この行列式は、2つのベクトルの座標を3つの座標軸(、、、および)を表す記号の下に配置することによって作成さi
れj
ますk
。クロス積の概念は3D構造であるため、3番目の(ゼロ)値の座標があります。クロス積を適用するために、これらの2Dベクトルを3Dに拡張します。
i j k
-4 0 0
1 4 0
すべての外積が、乗算される2つのベクトルの平面に垂直なベクトルを生成するとすると、上記の行列の行列式にはk
、(またはz軸)成分しかありません。 またはz 軸成分
の大きさを計算するための式は次のとおりです。k
a1*b2 - a2*b1 = -4* 4 - 0* 1
-16
この値の大きさ(-16
)は、2つの元のベクトル間の角度の正弦に、2つのベクトルの大きさの積を掛けたものです。
実際、その値の別の式はです
A X B (Cross Product) = |A| * |B| * sin(AB)
。
したがって、角度の測定値に戻るには、この値(-16
)を2つのベクトルの大きさの積で割る必要があります。
|A| * |B|
= 4 * Sqrt(17)
=16.4924...
したがって、sin(AB)の測度= -16 / 16.4924
=-.97014...
これは、頂点の次のセグメントが左または右に曲がっているかどうか、およびその程度の尺度です。アークサインを取る必要はありません。私たちが気にするのはその大きさ、そしてもちろんその符号(正または負)だけです!
閉じたパスの周りの他の4つのポイントのそれぞれに対してこれを行い、各頂点でこの計算からの値を合計します。
最終的な合計が正の場合、時計回り、負、反時計回りに移動します。
これは、 @ Betaの回答に基づくアルゴリズムの単純なC#実装です。
Vector
タイプを持つタイプX
とY
タイプのプロパティがあると仮定しましょう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;
}
%
は、モジュロ演算を実行するモジュロ演算子または剰余演算子であり、(ウィキペディアによると)ある数値を別の数値で除算した後に剰余を求めます。
@MichelRouzicのコメントによると最適化されたバージョン:
double sum = 0.0;
Vector v1 = vertices[vertices.Count - 1]; // or vertices[^1] with
// C# 8.0+ and .NET Core
for (int i = 0; i < vertices.Count; i++) {
Vector v2 = vertices[i];
sum += (v2.X - v1.X) * (v2.Y + v1.Y);
v1 = v2;
}
return sum > 0.0;
これにより、モジュロ演算だけでなく%
、配列のインデックス付けも節約できます。
頂点の1つから開始し、各辺のなす角を計算します。
最初と最後はゼロになります(したがって、それらをスキップします)。残りの部分では、角度の正弦は、(point [n] -point [0])と(point [n-1] -point [0])の単位長さへの正規化の外積によって与えられます。
値の合計が正の場合、ポリゴンは反時計回りに描画されます。
JavaScriptでのSeanの回答の実装:
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));
これが正しいことはかなり確かです。動作しているようです:-)
疑問に思っている場合、これらのポリゴンは次のようになります。
その価値については、このミックスインを使用して、Google MapsAPIv3アプリの巻き取り順序を計算しました。
このコードは、ポリゴン領域の副作用を利用しています。頂点の時計回りの巻き順序は正の領域を生成し、同じ頂点の反時計回りの巻き順序は負の値と同じ領域を生成します。このコードは、Googleマップジオメトリライブラリの一種のプライベートAPIも使用しています。私はそれを快適に使用できました-自己責任で使用してください。
使用例:
var myPolygon = new google.maps.Polygon({/*options*/});
var isCW = myPolygon.isPathClockwise();
ユニットテストを使用した完全な例@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 stevejansen_github@icloud.com
*
* @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;
}
})();
これは、 OpenLayers2に実装されている関数です。時計回りの多角形を持つための条件は、この参照area < 0
によって確認されます。
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);
}
lhfの答えを実装するためのC#コード:
// 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;
}
Matlabを使用するispolycw
場合、ポリゴンの頂点が時計回りの順序である場合、関数はtrueを返します。
このウィキペディアの記事で説明されているように、3つの点が与えられ、平面上(つまり、x座標とy座標)の曲線の向きを使用すると、次の行列式の符号を計算できます。p
q
r
行列式が負の場合(つまりOrient(p, q, r) < 0
)、ポリゴンは時計回り(CW)に向けられます。行列式が正の場合(つまりOrient(p, q, r) > 0
)、ポリゴンは反時計回り(CCW)に方向付けられます。行列式は、点の場合はゼロ(つまりOrient(p, q, r) == 0
)p
でq
ありr
、同一線上にあります。
p
上記の式では、の座標の前にあるものを追加します。これは、同次座標を使用しているためq
です。r
これは、他の回答の説明を使用した私の解決策です。
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
一部のポイントを時計回りに指定するには、エッジの合計だけでなく、すべてのエッジが正である必要があると思います。1つのエッジが負の場合、少なくとも3つのポイントが反時計回りに与えられます。
私のC#/ LINQソリューションは、以下の@charlesbretanaのクロス積アドバイスに基づいています。巻線の基準法線を指定できます。曲線が主に上向きベクトルによって定義された平面内にある限り、機能するはずです。
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;
}
}
}
ユニットテスト付き
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();
}
}
}
いくつかの信頼性の低い実装をテストした後、箱から出してCW / CCWの向きに関して満足のいく結果を提供したアルゴリズムは、このスレッドでOPによって投稿されたものでした(shoelace_formula_3
)。
いつものように、正の数はCW方向を表し、負の数はCCWを表します。
ポリゴン内のポイントをすでに知っている場合は、はるかに計算が簡単な方法です。
元のポリゴン、ポイント、およびそれらの座標から任意の線分をこの順序で選択します。
既知の「内側」の点を追加し、三角形を形成します。
これらの3つのポイントを使用して、ここで提案されているようにCWまたはCCWを計算します。
上記の回答に基づく迅速な3.0ソリューションは次のとおりです。
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
これに対する別の解決策。
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;
}
すべての頂点をこのような配列として取ります。
const vertices = [{x: 5, y: 0}, {x: 6, y: 4}, {x: 4, y: 5}, {x: 1, y: 5}, {x: 1, y: 0}];
isClockwise(vertices);
Rが方向を決定し、時計回りの場合は反転するためのソリューション(owinオブジェクトに必要であることがわかりました):
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
これらの答えは正しいですが、必要以上に数学的に強烈です。マップの座標を想定します。ここで、最も北のポイントがマップの最も高いポイントです。最も北のポイントを見つけます。2つのポイントが同点の場合、それは最も北で、次に最も東です(これは、lhfが彼の回答で使用するポイントです)。あなたのポイントでは、
point [0] =(5,0)
point [1] =(6,4)
point [2] =(4,5)
point [3] =(1,5)
ポイント[4]=(1,0)
P2が最も北にあると仮定すると、前のポイントまたは次のポイントのいずれかの東のポイントが時計回り、CW、またはCCWを決定します。最も北のポイントは北面にあるため、P1(前)からP2が東に移動すると、方向はCWになります。この場合、それは西に移動するので、受け入れられた答えが言うように、方向はCCWです。前のポイントに水平方向の動きがない場合、同じシステムが次のポイントP3に適用されます。P3がP2の西にある場合、その場合、移動はCCWです。P2からP3への移動が東、この場合は西、移動はCWです。nte、データのP2が最も北、次に東のポイントであり、prvが前のポイント、データのP1、nxtが次のポイント、データのP3、[0]が水平または東/であると仮定します。西は東よりも小さく、[1]は垂直です。
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)
これは、この回答に基づく単純なPython 3の実装です(これは、受け入れられた回答で提案されたソリューションに基づいています) 。
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
これらの点の重心を見つけます。
このポイントからあなたのポイントまでの線があると仮定します。
line0line1の2本の線の間の角度を求めます
line1とline2よりも
..。
..。
この角度が反時計回りよりも単調に増加している場合、
それ以外の場合、単調に減少する場合は時計回りです
それ以外(単調ではありません)
あなたは決めることができないので、それは賢明ではありません