Como determinar se uma lista de pontos de polígono está no sentido horário?
https://stackoverflow.com/questions/1165647
-
19-09-2019 - |
italiano
english
français
española
中国
日本の
العربية
Deutsch
한국어
Português
RussianPergunta
Tendo uma lista de pontos, como encontro se eles estão no sentido horário?
Por exemplo:
point[0] = (5,0)
point[1] = (6,4)
point[2] = (4,5)
point[3] = (1,5)
point[4] = (1,0)
diria que é anti-horário (ou no sentido anti-horário, para algumas pessoas).
Solução
Alguns dos métodos sugeridos falharão no caso de um polígono não convexo, como um crescente. Aqui está um simples que funcionará com polígonos não convexos (ele até funcionará com um polígono auto-intercalado como uma figura oito, dizendo se está majoritariamente sentido horário).
Soma sobre as bordas, (x2 - x1) (y2 + y1). Se o resultado for positivo, a curva é no sentido horário, se for negativo, a curva será no sentido anti-horário. (O resultado é o dobro da área fechada, com uma convenção +/-.)
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
Outras dicas
o produto cruzado mede o grau de perpendicular-estação de dois vetores. Imagine que cada borda do seu polígono é um vetor no plano XY de um espaço XYZ tridimensional (3-D). Em seguida, o produto cruzado de duas arestas sucessivas é um vetor na direção z (direção z positiva se o segundo segmento for no sentido horário, menos a direção z se for no sentido anti-horário). A magnitude desse vetor é proporcional ao seno do ângulo entre as duas bordas originais, por isso atinge um máximo quando é perpendicular e diminui para desaparecer quando as bordas são colineares (paralelas).
Portanto, para cada vértice (ponto) do polígono, calcule a magnitude do produto cruzado das duas bordas adjacentes:
Using your data:
point[0] = (5, 0)
point[1] = (6, 4)
point[2] = (4, 5)
point[3] = (1, 5)
point[4] = (1, 0)
Então rotule as bordas consecutivamente como
edgeA é o segmento de point0 para point1 e
edgeB entre point1 para point2
...
edgeE está entre point4 e point0.
Então vértice a (point0) está entre
edgeE A partir de point4 para point0]
edgeA A partir de point0 para `ponto1 '
Essas duas arestas são próprios vetores, cujas coordenadas X e Y podem ser determinadas subtraindo as coordenadas de seus pontos de início e final:
edgeE = point0 - point4 = (1, 0) - (5, 0) = (-4, 0) e
edgeA = point1 - point0 = (6, 4) - (1, 0) = (5, 4) e
E o produto cruzado dessas duas arestas adjacentes é calculado usando o determinante da matriz a seguir, que é construída colocando as coordenadas dos dois vetores abaixo dos símbolos que representam o eixo de três coordenadas (i, j, & k). A terceira coordenada com valor (zero) está lá porque o conceito de produto cruzado é um construto 3D e, portanto, estendemos esses vetores em 3D para 3D, a fim de aplicar o produto cruzado:
i j k
-4 0 0
1 4 0
Dado que todos os produtos cruzados produzem um vetor perpendicular ao plano de dois vetores sendo multiplicados, o determinante da matriz acima só tem um k, componente (ou zis z).
A fórmula para calcular a magnitude do k ou componente do eixo z é
a1*b2 - a2*b1 = -4* 4 - 0* 1 = -16
A magnitude desse valor (-16), é uma medida do seno do ângulo entre os 2 vetores originais, multiplicados pelo produto das magnitudes dos 2 vetores.
Na verdade, outra fórmula para seu valor é
A X B (Cross Product) = |A| * |B| * sin(AB).
Então, para voltar a apenas uma medida do ângulo, você precisa dividir esse valor, (-16), pelo produto das magnitudes dos dois vetores.
|A| * |B| = 4 * Sqrt(17) = 16.4924...
Então a medida do pecado (AB) = -16 / 16.4924 = -.97014...
Essa é uma medida de se o próximo segmento depois que o vértice está dobrado para a esquerda ou direita e quanto. Não há necessidade de tomar arco-sine. Tudo o que nos importamos é com sua magnitude e, claro, seu sinal (positivo ou negativo)!
Faça isso para cada um dos outros 4 pontos ao redor do caminho fechado e adicione os valores deste cálculo em cada vértice.
Se a soma final for positiva, você ficou no sentido horário, negativo, no sentido anti -horário.
Eu acho que essa é uma pergunta bastante antiga, mas vou jogar fora outra solução de qualquer maneira, porque é direta e não matematicamente intensiva - apenas usa álgebra básica. Calcule a área assinada do polígono. Se for negativo, os pontos estão no sentido horário, se for positivo, eles são no sentido anti -horário. (Isso é muito semelhante à solução da versão beta.)
Calcule a área assinada: a = 1/2 * (x1*y2 - x2*y1 + x2*y3 - x3*y2 + ... + xn*y1 - x1*yn)
Ou em pseudo-código:
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
Observe que, se você estiver verificando apenas o pedido, não precisa se preocupar em dividir por 2.
Encontre o vértice com o menor y (e o maior x se houver laços). Deixe o vértice ser A e o vértice anterior da lista é B e o próximo vértice da lista é C. Agora calcule o sinal do produto cruzado de AB e AC.
Referências:
Como encontro a orientação de um polígono simples? dentroPerguntas frequentes: comp.graphics.algorithms.
Orientação da curva na Wikipedia.
Aqui está uma simples implementação C# do algoritmo baseado em esta resposta.
Vamos supor que tenhamos um Vector Tipo tendo X e Y propriedades do tipo 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;
}
% é o módulo ou operador restante executando a operação do módulo que (De acordo com a Wikipedia) encontra o restante após a divisão de um número por outro.
Comece em um dos vértices e calcule o ângulo subtendido por cada lado.
O primeiro e o último será zero (então pule isso); Para restos, o seno do ângulo será dado pelo produto cruzado das normalizações ao comprimento da unidade de (ponto [n] -po [0]) e (ponto [n-1] -po [0]).
Se a soma dos valores for positiva, seu polígono será desenhado no sentido anti-horário.
Pelo que vale a pena, usei esse mixin para calcular a ordem de enrolamento para aplicativos API V3 do Google Maps.
O código aproveita o efeito colateral das áreas de polígono: uma ordem de enrolamento no sentido horário de vértices produz uma área positiva, enquanto uma ordem de enrolamento no sentido anti-horário dos mesmos vértexes produz a mesma área que um valor negativo. O código também usa uma espécie de API privada na biblioteca de geometria do Google Maps. Eu me senti confortável em usá -lo - use por sua conta e risco.
Uso de amostra:
var myPolygon = new google.maps.Polygon({/*options*/});
var isCW = myPolygon.isPathClockwise();
Exemplo completo com testes de unidade @ 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;
}
})();
Uma implementação de Resposta de Sean em 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));Tenho certeza de que está certo. Parece que está funcionando :-)
Esses polígonos são assim, se você está se perguntando:
Se você usa o MATLAB, a função ispolycw Retorna true se os vértices de polígono estiverem no sentido horário.
Esta é a função implementada para OpenLayers 2. A condição para ter um polígono no sentido horário é area < 0, confirmou por esta referência.
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);
}
Como também explicado neste artigo da Wikipedia Orientação da curva, dado 3 pontos p, q e r No avião (ou seja, com coordenadas X e Y), você pode calcular o sinal do seguinte determinante
Se o determinante for negativo (ou seja, Orient(p, q, r) < 0), então o polígono é orientado no sentido horário (CW). Se o determinante for positivo (ou seja, Orient(p, q, r) > 0), o polígono é orientado no sentido anti -horário (CCW). O determinante é zero (ou seja, Orient(p, q, r) == 0) se pontos p, q e r são Collinear.
Na fórmula acima, precendemos aqueles em frente às coordenadas de p, q e r Porque estamos usando coordenadas homogêneas.
Eu acho que para que alguns pontos sejam dados no sentido horário, todas as arestas precisam ser positivas não apenas a soma das arestas. Se uma borda for negativa, pelo menos 3 pontos forem recebidos no sentido anti-horário.
Minha solução C# / LINQ é baseada nos conselhos cruzados do @Charlesbretana está abaixo. Você pode especificar uma referência normal para o enrolamento. Deve funcionar enquanto a curva estiver principalmente no plano definido pelo vetor da UP.
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;
}
}
}
com um teste de unidade
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();
}
}
}
Esta é a minha solução usando as explicações nas outras respostas:
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
Um método muito mais simples, computacionalmente mais simples, Se você já conhece um ponto dentro do polígono:
Escolha qualquer segmento de linha no polígono original, pontos e suas coordenadas nessa ordem.
Adicione um ponto "interno" conhecido e forme um triângulo.
Calcule CW ou CCW conforme sugerido aqui com esses três pontos.
Depois de testar várias implementações não confiáveis, o algoritmo que forneceu resultados satisfatórios em relação à orientação CW/CCW fora da caixa foi o, publicado por OP em isto fio (shoelace_formula_3).
Como sempre, um número positivo representa uma orientação CW, enquanto um número negativo CCW.
Aqui está a solução Swift 3.0 com base nas respostas acima:
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
Outra solução para isso;
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;
}
Pegue todos os vértices como uma matriz como essa;
const vertices = [{x: 5, y: 0}, {x: 6, y: 4}, {x: 4, y: 5}, {x: 1, y: 5}, {x: 1, y: 0}];
isClockwise(vertices);
Solução para r para determinar a direção e reverter se no sentido horário (achado necessário para objetos 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
Embora essas respostas estejam corretas, elas são mais matematicamente intensas do que o necessário. Suponha que as coordenadas do mapa, onde o ponto mais norte é o ponto mais alto do mapa. Encontre o mais North Point e, se 2 pontos empatarem, é o mais norte do que o mais leste (este é o ponto que LHF usa em sua resposta). Em seus pontos,
ponto [0] = (5,0)
ponto [1] = (6,4)
ponto [2] = (4,5)
ponto [3] = (1,5)
ponto [4] = (1,0)
Se assumirmos que P2 é o mais norte do que o East Point, ou o próximo ponto anterior ou o próximo, determine no sentido horário, CW ou CCW. Como o ponto mais norte está na face norte, se o P1 (anterior) a P2 se mover para o leste, a direção for cw. Nesse caso, ele se move para o oeste, então a direção é CCW como diz a resposta aceita. Se o ponto anterior não tiver movimento horizontal, o mesmo sistema se aplica ao próximo ponto, p3. Se p3 estiver a oeste de P2, é, então o movimento é CCW. Se o movimento P2 a P3 é leste, é oeste neste caso, o movimento é CW. Suponha que NTE, P2 em seus dados, seja o mais norte do que o leste e o PRV é o ponto anterior, P1 em seus dados, e o NXT é o próximo ponto, p3 em seus dados e [0] é horizontal ou leste/ Oeste, onde oeste é menor que o leste e [1] é 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)
Código C# para implementar A resposta do LHF:
// 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;
}
Encontre o centro de massa desses pontos.
Suponha que haja linhas deste ponto para seus pontos.
Encontre o ângulo entre duas linhas para a linha0 line1
do que fazer isso para linha1 e linha2
...
...
Se esse ângulo estiver monotonicamente aumentando do que é no sentido anti -horário,
caso contrário, se diminuir monotonicamente, está no sentido horário
caso contrário (não é monotônico)
você não pode decidir, então não é sábio
