Como posso determinar se um ponto 2D está dentro de um polígono?
https://stackoverflow.com/questions/217578
-
03-07-2019 - |
italiano
english
français
española
中国
日本の
العربية
Deutsch
한국어
Português
RussianPergunta
Eu estou tentando criar um rápido 2D algoritmo polígono ponto dentro, para uso em hit-teste (por exemplo Polygon.contains(p:Point)). Sugestões para técnicas eficazes seria apreciada.
Solução
Para os gráficos, eu prefiro não preferem inteiros. Muitos sistemas usam números inteiros para a pintura UI (pixels são ints depois de tudo), mas MacOS por exemplo, usa flutuar para tudo. MacOS só sabe pontos e um ponto pode se traduzir em um pixel, mas dependendo da resolução do monitor, pode traduzir para outra coisa. Em telas de retina metade de um ponto (0,5 / 0,5) é de pixel. Ainda assim, eu nunca percebi que o MacOS UIs são significativamente mais lento do que outras UIs. Afinal APIs 3D (OpenGL ou Direct3D) também trabalha com carros alegóricos e gráficos modernos bibliotecas, muitas vezes tirar proveito da aceleração GPU.
Agora você disse que a velocidade é a sua principal preocupação, ok, vamos para a velocidade. Antes de executar qualquer algoritmo sofisticado, primeiro fazer um teste simples. Criar uma caixa de eixo alinhado delimitadora em torno de seu polígono. Isto é muito fácil, rápida e já pode segura-lhe um monte de cálculos. Como isso funciona? Iterar sobre todos os pontos do polígono e encontrar os valores mínimo / máximo de X e Y.
por exemplo. você tem a pontos (9/1), (4/3), (2/7), (8/2), (3/6). Este meio é Xmin 2, Xmax é 9, Ymin é 1 e Ymax é 7. Um ponto fora de um rectângulo com as duas bordas (2/1) e (9/7) não pode estar dentro do polígono.
// p is your point, p.x is the x coord, p.y is the y coord
if (p.x < Xmin || p.x > Xmax || p.y < Ymin || p.y > Ymax) {
// Definitely not within the polygon!
}
Este é o primeiro teste a ser executado para qualquer ponto. Como você pode ver, este teste é ultra-rápido, mas também é muito grossa. Para lidar com pontos que estão dentro do retângulo delimitador, precisamos de um algoritmo mais sofisticado. Há um par de maneiras como isso pode ser calculado. Qual método funciona também depende do fato se o polígono pode ter furos ou será sempre sólida. Aqui estão alguns exemplos de uns sólidos (um convexas, um côncavo):
E um de aqui com um buraco:
O verde tem um buraco no meio!
O algoritmo mais fácil, que pode lidar com todos os três casos acima e ainda está muito rápido é chamado de fundição ray . A idéia do algoritmo é bastante simples: desenhar um raio virtual a partir de qualquer lugar fora do polígono para o seu ponto e contar quantas vezes ele bate um lado do polígono. Se o número de visitas é ainda, é fora do polígono, se é estranho, ele está dentro.
O número de enrolamento algoritmo seria uma alternativa, é mais preciso para pontos que estão sendo muito perto de uma linha de polígono mas também é muito mais lento. Ray fundição pode falhar por pontos muito perto de um lado polígono por causa do limitado flutuante precisão de ponto e arredondamento problemas, mas na realidade que é quase um problema, como se um ponto está tão perto de um lado, é muitas vezes visualmente nem mesmo possível para um espectador a reconhecer se ele já está dentro ou fora ainda.
Você ainda tem a caixa delimitadora de acima, lembra? Basta escolher um ponto fora da caixa delimitadora e usá-lo como ponto de partida para o seu raio. Por exemplo. o (Xmin - e/p.y) ponto está fora do polígono com certeza.
Mas o que é e? Bem, e (na verdade epsilon) dá a caixa delimitadora alguns estofamento . Como eu disse, o traçado de raios falhar se começarmos muito perto de uma linha de polígono. Desde a caixa delimitadora pode igualar o polígono (se o polígono é um retângulo eixo alinhado, a caixa delimitadora é igual ao próprio! Polígono), precisamos de algum estofo para fazer este seguro, isso é tudo. Como grande você deve escolher e? Não muito grande. Depende do sistema de coordenadas escala que você usa para desenhar. Se o seu passo de pixels de largura é de 1,0, então basta escolher o 1.0 (ainda 0,1 teria funcionado tão bem)
Agora que temos o raio com suas coordenadas iniciais e finais, as mudanças problema " é o ponto dentro do polígono " para " Com que freqüência o cruza ray um lado polígono ". Portanto, não podemos apenas trabalhocom os pontos de polígono como antes, agora precisamos dos lados reais. Um lado é sempre definida por dois pontos.
side 1: (X1/Y1)-(X2/Y2)
side 2: (X2/Y2)-(X3/Y3)
side 3: (X3/Y3)-(X4/Y4)
:
Você precisa testar o raio contra todos os lados. Considere o raio para ser um vetor e cada lado para ser um vetor. O raio tem que bater cada lado exatamente uma vez ou nunca. Não pode bater o mesmo lado duas vezes. Duas linhas em espaço 2D será sempre cruzam exatamente uma vez, a menos que eles são paralelos, caso em que eles nunca se cruzam. Contudo, desde vetores têm um comprimento limitado, dois vetores não pode ser paralelo e ainda nunca se cruzam, porque eles são muito curta para nunca conhecer uns aos outros.
// Test the ray against all sides
int intersections = 0;
for (side = 0; side < numberOfSides; side++) {
// Test if current side intersects with ray.
// If yes, intersections++;
}
if ((intersections & 1) == 1) {
// Inside of polygon
} else {
// Outside of polygon
}
Até aí tudo bem, mas como você testar se dois vetores se cruzam? Aqui está algum código C (não testada), que deve fazer o truque:
#define NO 0
#define YES 1
#define COLLINEAR 2
int areIntersecting(
float v1x1, float v1y1, float v1x2, float v1y2,
float v2x1, float v2y1, float v2x2, float v2y2
) {
float d1, d2;
float a1, a2, b1, b2, c1, c2;
// Convert vector 1 to a line (line 1) of infinite length.
// We want the line in linear equation standard form: A*x + B*y + C = 0
// See: http://en.wikipedia.org/wiki/Linear_equation
a1 = v1y2 - v1y1;
b1 = v1x1 - v1x2;
c1 = (v1x2 * v1y1) - (v1x1 * v1y2);
// Every point (x,y), that solves the equation above, is on the line,
// every point that does not solve it, is not. The equation will have a
// positive result if it is on one side of the line and a negative one
// if is on the other side of it. We insert (x1,y1) and (x2,y2) of vector
// 2 into the equation above.
d1 = (a1 * v2x1) + (b1 * v2y1) + c1;
d2 = (a1 * v2x2) + (b1 * v2y2) + c1;
// If d1 and d2 both have the same sign, they are both on the same side
// of our line 1 and in that case no intersection is possible. Careful,
// 0 is a special case, that's why we don't test ">=" and "<=",
// but "<" and ">".
if (d1 > 0 && d2 > 0) return NO;
if (d1 < 0 && d2 < 0) return NO;
// The fact that vector 2 intersected the infinite line 1 above doesn't
// mean it also intersects the vector 1. Vector 1 is only a subset of that
// infinite line 1, so it may have intersected that line before the vector
// started or after it ended. To know for sure, we have to repeat the
// the same test the other way round. We start by calculating the
// infinite line 2 in linear equation standard form.
a2 = v2y2 - v2y1;
b2 = v2x1 - v2x2;
c2 = (v2x2 * v2y1) - (v2x1 * v2y2);
// Calculate d1 and d2 again, this time using points of vector 1.
d1 = (a2 * v1x1) + (b2 * v1y1) + c2;
d2 = (a2 * v1x2) + (b2 * v1y2) + c2;
// Again, if both have the same sign (and neither one is 0),
// no intersection is possible.
if (d1 > 0 && d2 > 0) return NO;
if (d1 < 0 && d2 < 0) return NO;
// If we get here, only two possibilities are left. Either the two
// vectors intersect in exactly one point or they are collinear, which
// means they intersect in any number of points from zero to infinite.
if ((a1 * b2) - (a2 * b1) == 0.0f) return COLLINEAR;
// If they are not collinear, they must intersect in exactly one point.
return YES;
}
Os valores de entrada são os dois pontos finais de um vector (v1x1/v1y1 e v1x2/v1y2) e o vector 2 (v2x1/v2y1 e v2x2/v2y2). Então você tem 2 vetores, 4 pontos, 8 coordenadas. YES e NO são claras. YES aumenta cruzamentos, NO não faz nada.
E sobre colineares? Isso significa ambos os vectores encontram-se na mesma linha infinito, dependendo da posição e comprimento, eles não se cruzam em todos eles ou se cruzam em um infindável número de pontos. Não estou absolutamente certo de como lidar com este caso, eu não contá-lo como interseção de qualquer maneira. Bem, neste caso, é bastante raro na prática de qualquer maneira por causa de erros pontuais arredondamento flutuante; um código melhor, provavelmente, não teste para == 0.0f mas em vez de algo como < epsilon, onde epsilon é um número bastante pequeno.
Se você precisa testar um número maior de pontos, você pode certamente acelerar a coisa toda um pouco, mantendo a equação linear formulários dos lados do polígono na memória, para que você não tem que recalcular estes cada vez. Isto vai poupar duas multiplicações de vírgula flutuante e três flutuantes subtrações ponto em cada teste, em troca de armazenar três valores de ponto flutuante por lado polígono na memória. É uma lembrança típica vs comércio tempo de computação off.
Por último, mas não menos importante: Se você pode usar o hardware 3D para resolver o problema, não é uma alternativa interessante. Apenas deixe o GPU fazer todo o trabalho para você. Criar uma superfície de pintura que é fora da tela. Preenchê-lo completamente com a cor preta. Agora vamos OpenGL ou Direct3D pintar o seu polígono (ou mesmo todos os seus polígonos se você quiser apenas para testar se o ponto está dentro de qualquer um deles, mas você não se importa com qual) e preencher o polígono (s) com um diferente cor, por exemplo, branco. Para verificar se um ponto está dentro do polígono, obter a cor deste ponto da superfície de desenho. Esta é apenas uma O (1) Memória buscar.
É claro que este método só é útil se a sua superfície de desenho não tem que ser enorme. Se ele não pode caber na memória GPU, este método é mais lento do que fazê-lo na CPU. Se ele teria que ser enorme e sua GPU suporta shaders modernos, você ainda pode usar a GPU através da aplicação do elenco ray mostrado acima como um shader GPU, o que é absolutamente possível. Para um número maior de polígonos ou um grande número de pontos de teste, isso vai pagar, considere algumas GPUs será capaz de teste de 64 a 256 pontos em paralelo. Note, porém, que a transferência de dados de CPU para GPU e volta é sempre caro, por isso apenas para testar um par de pontos contra um par de polígonos simples, onde tanto os pontos ou os polígonos são dinâmicos e mudam com freqüência, uma abordagem GPU raramente vai pagar off.
Outras dicas
Eu acho o seguinte pedaço de código é a melhor solução (retirado aqui ):
int pnpoly(int nvert, float *vertx, float *verty, float testx, float testy)
{
int i, j, c = 0;
for (i = 0, j = nvert-1; i < nvert; j = i++) {
if ( ((verty[i]>testy) != (verty[j]>testy)) &&
(testx < (vertx[j]-vertx[i]) * (testy-verty[i]) / (verty[j]-verty[i]) + vertx[i]) )
c = !c;
}
return c;
}
Argumentos
- nvert : Número de vértices do polígono. Se repetir o primeiro vértice no final foi discutida no artigo mencionado acima.
- VertX, verty :. Arrays contendo os X e Y coordenadas dos vértices do polígono
- testx, irritado :. X e Y coordenadas do ponto de teste
É curto e eficiente e funciona tanto para convexo e polígonos côncavos. Como sugerido antes, você deve verificar o retângulo delimitador primeira e buracos polígono tratar separadamente.
A idéia por trás disso é muito simples. O autor descreve-o como segue:
Eu corro um raio semi-infinita na horizontal (aumento de x, y fixo) a partir do ponto de teste, e contar quantas bordas cruza. Em cada passagem, o raio muda entre o interior e exterior. Isso é chamado a curva teorema de Jordan.
A variável c é de comutação de 0 a 1 e 1-0 cada vez que o raio horizontal cruza qualquer borda. Então, basicamente, é manter o controle sobre se o número de arestas cruzados são par ou ímpar. 0 significa mesmo e 1 significa ímpares.
Aqui está uma versão C # da resposta dada pelo nirg , que vem de este RPI professora . Note-se que o uso do código a partir dessa fonte RPI requer atribuição.
Uma verificação caixa delimitadora foi adicionado no topo. No entanto, como James Brown aponta, o código principal é quase tão rápido como a caixa delimitadora verificar-se, por isso, a verificação de caixa delimitadora pode realmente retardar o funcionamento global, no caso em que a maioria dos pontos que está conferindo estão dentro da caixa delimitadora . Então você pode deixar a caixa delimitadora check-out, ou uma alternativa seria a de precompute as caixas delimitadoras dos seus polígonos se eles não mudam de forma muitas vezes.
public bool IsPointInPolygon( Point p, Point[] polygon )
{
double minX = polygon[ 0 ].X;
double maxX = polygon[ 0 ].X;
double minY = polygon[ 0 ].Y;
double maxY = polygon[ 0 ].Y;
for ( int i = 1 ; i < polygon.Length ; i++ )
{
Point q = polygon[ i ];
minX = Math.Min( q.X, minX );
maxX = Math.Max( q.X, maxX );
minY = Math.Min( q.Y, minY );
maxY = Math.Max( q.Y, maxY );
}
if ( p.X < minX || p.X > maxX || p.Y < minY || p.Y > maxY )
{
return false;
}
// http://www.ecse.rpi.edu/Homepages/wrf/Research/Short_Notes/pnpoly.html
bool inside = false;
for ( int i = 0, j = polygon.Length - 1 ; i < polygon.Length ; j = i++ )
{
if ( ( polygon[ i ].Y > p.Y ) != ( polygon[ j ].Y > p.Y ) &&
p.X < ( polygon[ j ].X - polygon[ i ].X ) * ( p.Y - polygon[ i ].Y ) / ( polygon[ j ].Y - polygon[ i ].Y ) + polygon[ i ].X )
{
inside = !inside;
}
}
return inside;
}
Aqui está uma variante JavaScript da resposta por M. Katz baseado em abordagem de Nirg:
function pointIsInPoly(p, polygon) {
var isInside = false;
var minX = polygon[0].x, maxX = polygon[0].x;
var minY = polygon[0].y, maxY = polygon[0].y;
for (var n = 1; n < polygon.length; n++) {
var q = polygon[n];
minX = Math.min(q.x, minX);
maxX = Math.max(q.x, maxX);
minY = Math.min(q.y, minY);
maxY = Math.max(q.y, maxY);
}
if (p.x < minX || p.x > maxX || p.y < minY || p.y > maxY) {
return false;
}
var i = 0, j = polygon.length - 1;
for (i, j; i < polygon.length; j = i++) {
if ( (polygon[i].y > p.y) != (polygon[j].y > p.y) &&
p.x < (polygon[j].x - polygon[i].x) * (p.y - polygon[i].y) / (polygon[j].y - polygon[i].y) + polygon[i].x ) {
isInside = !isInside;
}
}
return isInside;
}
calcular a soma orientada de ângulos entre o ponto P e cada um dos ápices polígono. Se o ângulo orientado total é de 360 ??graus, o ponto está dentro. Se o total for 0, o ponto está fora.
Eu gosto deste método melhor porque é mais robusto e menos dependente de precisão numérica.
Os métodos que uniformidade de computação de número de cruzamentos são limitados porque você pode 'hit' um ápice durante o cálculo do número de cruzamentos.
EDIT:. By The Way, este método funciona com polígonos côncavos e convexos
EDIT:. Eu encontrei recentemente um todo Wikipedia artigo sobre o tema
O Eric Haines artigo citado por bobobobo é realmente excelente. Particularmente interessantes são as tabelas que comparam o desempenho dos algoritmos; o método ângulo somatório é muito ruim em comparação com os outros. Também interessante é que as otimizações como o uso de uma grade de pesquisa para mais subdividir o polígono em "in" e "out" setores pode fazer o teste incrivelmente rápido mesmo em polígonos com> 1000 lados.
De qualquer forma, é cedo, mas meu voto vai para o método "cruzamentos", que é muito bonito o que Mecki descreve eu acho. No entanto eu achei mais sucintamente descrito e codificado por David Bourke . Eu amo que não há trigonometria verdadeiro necessário, e ele funciona para convexas e côncavas, e ele executa razoavelmente bem como o número de lados aumenta.
A propósito, aqui está uma das mesas de desempenho do artigo, o Eric Haines' para o interesse, testando em polígonos aleatórios.
number of edges per polygon
3 4 10 100 1000
MacMartin 2.9 3.2 5.9 50.6 485
Crossings 3.1 3.4 6.8 60.0 624
Triangle Fan+edge sort 1.1 1.8 6.5 77.6 787
Triangle Fan 1.2 2.1 7.3 85.4 865
Barycentric 2.1 3.8 13.8 160.7 1665
Angle Summation 56.2 70.4 153.6 1403.8 14693
Grid (100x100) 1.5 1.5 1.6 2.1 9.8
Grid (20x20) 1.7 1.7 1.9 5.7 42.2
Bins (100) 1.8 1.9 2.7 15.1 117
Bins (20) 2.1 2.2 3.7 26.3 278
Esta questão é tão interessante. Eu tenho uma outra idéia diferente viável de outras respostas deste post.
A idéia é usar a soma dos ângulos de decidir se o alvo está dentro ou fora. Se o alvo está dentro de uma área, a soma de forma ângulo pelo alvo e a cada dois pontos de fronteira será 360. Se o destino estiver fora, a soma não será 360. O ângulo tem sentido. Se o ângulo indo para trás, o ângulo é negativa. Isso é como calcular o número enrolamento .
Consulte esta imagem para obter uma compreensão básica da idéia:
O meu algoritmo assume o horário é o sentido positivo. Aqui está uma entrada potencial:
[[-122.402015, 48.225216], [-117.032049, 48.999931], [-116.919132, 45.995175], [-124.079107, 46.267259], [-124.717175, 48.377557], [-122.92315, 47.047963], [-122.402015, 48.225216]]
A seguir está o código python que implementa a idéia:
def isInside(self, border, target):
degree = 0
for i in range(len(border) - 1):
a = border[i]
b = border[i + 1]
# calculate distance of vector
A = getDistance(a[0], a[1], b[0], b[1]);
B = getDistance(target[0], target[1], a[0], a[1])
C = getDistance(target[0], target[1], b[0], b[1])
# calculate direction of vector
ta_x = a[0] - target[0]
ta_y = a[1] - target[1]
tb_x = b[0] - target[0]
tb_y = b[1] - target[1]
cross = tb_y * ta_x - tb_x * ta_y
clockwise = cross < 0
# calculate sum of angles
if(clockwise):
degree = degree + math.degrees(math.acos((B * B + C * C - A * A) / (2.0 * B * C)))
else:
degree = degree - math.degrees(math.acos((B * B + C * C - A * A) / (2.0 * B * C)))
if(abs(round(degree) - 360) <= 3):
return True
return False
Versão Swift da resposta href="https://stackoverflow.com/a/2922778/3310339"> :
extension CGPoint {
func isInsidePolygon(vertices: [CGPoint]) -> Bool {
guard !vertices.isEmpty else { return false }
var j = vertices.last!, c = false
for i in vertices {
let a = (i.y > y) != (j.y > y)
let b = (x < (j.x - i.x) * (y - i.y) / (j.y - i.y) + i.x)
if a && b { c = !c }
j = i
}
return c
}
}
Eu fiz alguns trabalhos sobre esta volta quando eu era um investigador no âmbito Michael Stonebraker - você sabe, o professor que veio com Ingres , PostgreSQL , etc
Nós percebemos que a maneira mais rápida foi primeiro fazer uma caixa delimitadora porque é SUPER rápido. Se é fora da caixa delimitadora, é lá fora. Caso contrário, você faz o trabalho mais difícil ...
Se você quiser uma grande algoritmo, olhar para o código-fonte PostgreSQL projeto de código aberto para o trabalho geo ...
Eu quero salientar, nunca tivemos qualquer conhecimento sobre direito vs lateralidade esquerda (também expresso como um "dentro" vs problema "fora" ...
Atualização
link do BKB desde um bom número de algoritmos razoáveis. Eu estava trabalhando em problemas de Ciências da Terra e, portanto, precisava de uma solução que funciona em latitude / longitude, e tem o problema peculiar de destreza manual - é a área dentro da área menor ou a maior área? A resposta é que a "direção" das matérias vértices - é tanto canhoto ou destro e, desta forma você pode indicar qualquer área como "dentro" qualquer polígono. Como tal, a minha solução trabalho usado três enumerados nessa página.
Além disso, meu trabalho usado funções separadas por "na linha" testes.
... Desde que alguém perguntou: descobrimos que limitam os testes caixa foram melhor quando o número de vértices foi além de um número - fazer um teste muito rápido antes de fazer o teste mais longo, se necessário ... uma caixa delimitadora é criado por simplesmente tomar o maior x, menores x, y maior e menor y e colocá-los juntos para fazer quatro pontos de uma caixa de ...
Outra dica para aqueles que seguem: nós fizemos todo o nosso mais sofisticado e "light-dimming" computação em um espaço de grade tudo em pontos positivos em um avião e, em seguida, re-projetado de volta em "real" de longitude / latitude, evitando assim possíveis erros de envolvimento ao redor de quando uma linha cruzada 180 de longitude e ao manusear regiões polares. Trabalhou ótimo!
Realmente como a solução postado por Nirg e editado por bobobobo. Eu só fiz isso JavaScript amigável e um pouco mais legível para meu uso:
function insidePoly(poly, pointx, pointy) {
var i, j;
var inside = false;
for (i = 0, j = poly.length - 1; i < poly.length; j = i++) {
if(((poly[i].y > pointy) != (poly[j].y > pointy)) && (pointx < (poly[j].x-poly[i].x) * (pointy-poly[i].y) / (poly[j].y-poly[i].y) + poly[i].x) ) inside = !inside;
}
return inside;
}
A resposta de David Segond é praticamente a resposta geral padrão, e Richard T de é a otimização mais comum, embora therre são alguns outros. Outras otimizações fortes são baseados em soluções menos gerais. Por exemplo, se você estiver indo para verificar o mesmo polígono com lotes de pontos, triangulação do polígono pode acelerar as coisas imensamente, pois há uma série de algoritmos TIN pesquisando muito rapidamente. Outra é se o polígono e os pontos estão em um plano limitado em baixa resolução, dizer uma tela, você pode pintar o polígono em um buffer de exibição memória mapeada em uma determinada cor, e verificar a cor de um determinado pixel para ver se ele mentiras nos polígonos.
Como muitas otimizações, estes são baseados em específico em vez de casos gerais, e beneifits rendimento com base no tempo amortizado em vez de uso único.
Trabalho neste campo, eu encontrei Joeseph O'Rourkes 'Computação Geometria em C' ISBN 0-521-44034-3 para ser uma grande ajuda.
A solução trivial seria dividir o polígono de triângulos e teste de ocorrência dos triângulos como explicado aqui
Se o seu polígono é CONVEX pode haver uma melhor abordagem embora. Olhe para o polígono como uma coleção de linhas infinitas. Cada linha que divide o espaço em dois. para cada ponto, é fácil dizer se a sua de um lado ou do outro lado da linha. Se um ponto é no mesmo lado de todas as linhas, em seguida, é no interior do polígono.
Sei que isso é velho, mas aqui é um algoritmo de lançamento de raios implementado em Cocoa, no caso de alguém está interessado. Não tenho certeza que é a maneira mais eficiente de fazer as coisas, mas pode ajudar alguém.
- (BOOL)shape:(NSBezierPath *)path containsPoint:(NSPoint)point
{
NSBezierPath *currentPath = [path bezierPathByFlatteningPath];
BOOL result;
float aggregateX = 0; //I use these to calculate the centroid of the shape
float aggregateY = 0;
NSPoint firstPoint[1];
[currentPath elementAtIndex:0 associatedPoints:firstPoint];
float olderX = firstPoint[0].x;
float olderY = firstPoint[0].y;
NSPoint interPoint;
int noOfIntersections = 0;
for (int n = 0; n < [currentPath elementCount]; n++) {
NSPoint points[1];
[currentPath elementAtIndex:n associatedPoints:points];
aggregateX += points[0].x;
aggregateY += points[0].y;
}
for (int n = 0; n < [currentPath elementCount]; n++) {
NSPoint points[1];
[currentPath elementAtIndex:n associatedPoints:points];
//line equations in Ax + By = C form
float _A_FOO = (aggregateY/[currentPath elementCount]) - point.y;
float _B_FOO = point.x - (aggregateX/[currentPath elementCount]);
float _C_FOO = (_A_FOO * point.x) + (_B_FOO * point.y);
float _A_BAR = olderY - points[0].y;
float _B_BAR = points[0].x - olderX;
float _C_BAR = (_A_BAR * olderX) + (_B_BAR * olderY);
float det = (_A_FOO * _B_BAR) - (_A_BAR * _B_FOO);
if (det != 0) {
//intersection points with the edges
float xIntersectionPoint = ((_B_BAR * _C_FOO) - (_B_FOO * _C_BAR)) / det;
float yIntersectionPoint = ((_A_FOO * _C_BAR) - (_A_BAR * _C_FOO)) / det;
interPoint = NSMakePoint(xIntersectionPoint, yIntersectionPoint);
if (olderX <= points[0].x) {
//doesn't matter in which direction the ray goes, so I send it right-ward.
if ((interPoint.x >= olderX && interPoint.x <= points[0].x) && (interPoint.x > point.x)) {
noOfIntersections++;
}
} else {
if ((interPoint.x >= points[0].x && interPoint.x <= olderX) && (interPoint.x > point.x)) {
noOfIntersections++;
}
}
}
olderX = points[0].x;
olderY = points[0].y;
}
if (noOfIntersections % 2 == 0) {
result = FALSE;
} else {
result = TRUE;
}
return result;
}
Obj-C da resposta de nirg com o método de amostra para testar pontos. A resposta de Nirg funcionou bem para mim.
- (BOOL)isPointInPolygon:(NSArray *)vertices point:(CGPoint)test {
NSUInteger nvert = [vertices count];
NSInteger i, j, c = 0;
CGPoint verti, vertj;
for (i = 0, j = nvert-1; i < nvert; j = i++) {
verti = [(NSValue *)[vertices objectAtIndex:i] CGPointValue];
vertj = [(NSValue *)[vertices objectAtIndex:j] CGPointValue];
if (( (verti.y > test.y) != (vertj.y > test.y) ) &&
( test.x < ( vertj.x - verti.x ) * ( test.y - verti.y ) / ( vertj.y - verti.y ) + verti.x) )
c = !c;
}
return (c ? YES : NO);
}
- (void)testPoint {
NSArray *polygonVertices = [NSArray arrayWithObjects:
[NSValue valueWithCGPoint:CGPointMake(13.5, 41.5)],
[NSValue valueWithCGPoint:CGPointMake(42.5, 56.5)],
[NSValue valueWithCGPoint:CGPointMake(39.5, 69.5)],
[NSValue valueWithCGPoint:CGPointMake(42.5, 84.5)],
[NSValue valueWithCGPoint:CGPointMake(13.5, 100.0)],
[NSValue valueWithCGPoint:CGPointMake(6.0, 70.5)],
nil
];
CGPoint tappedPoint = CGPointMake(23.0, 70.0);
if ([self isPointInPolygon:polygonVertices point:tappedPoint]) {
NSLog(@"YES");
} else {
NSLog(@"NO");
}
}
C # versão da resposta de nirg está aqui: Eu só vou compartilhar o código. Pode salvar alguém algum tempo.
public static bool IsPointInPolygon(IList<Point> polygon, Point testPoint) {
bool result = false;
int j = polygon.Count() - 1;
for (int i = 0; i < polygon.Count(); i++) {
if (polygon[i].Y < testPoint.Y && polygon[j].Y >= testPoint.Y || polygon[j].Y < testPoint.Y && polygon[i].Y >= testPoint.Y) {
if (polygon[i].X + (testPoint.Y - polygon[i].Y) / (polygon[j].Y - polygon[i].Y) * (polygon[j].X - polygon[i].X) < testPoint.X) {
result = !result;
}
}
j = i;
}
return result;
}
Não há nada mais linda do que uma definição indutiva de um problema. Por uma questão de exaustividade, aqui você tem uma versão em prólogo que também pode esclarecer os thoughs atrás fundição ray :
Com base na simulação do algoritmo simplicidade em http: // www.ecse.rpi.edu/Homepages/wrf/Research/Short_Notes/pnpoly.html
Alguns predicados auxiliares:
exor(A,B):- \+A,B;A,\+B.
in_range(Coordinate,CA,CB) :- exor((CA>Coordinate),(CB>Coordinate)).
inside(false).
inside(_,[_|[]]).
inside(X:Y, [X1:Y1,X2:Y2|R]) :- in_range(Y,Y1,Y2), X > ( ((X2-X1)*(Y-Y1))/(Y2-Y1) + X1),toggle_ray, inside(X:Y, [X2:Y2|R]); inside(X:Y, [X2:Y2|R]).
get_line(_,_,[]).
get_line([XA:YA,XB:YB],[X1:Y1,X2:Y2|R]):- [XA:YA,XB:YB]=[X1:Y1,X2:Y2]; get_line([XA:YA,XB:YB],[X2:Y2|R]).
A equação de uma linha dada 2 pontos A e B (linha (A, B)) é:
(YB-YA)
Y - YA = ------- * (X - XA)
(XB-YB)
É importante que o sentido de rotação para a linha é setado para ponteiros do relógio para fronteiras e-anti-relógio sábio para furos. Estamos indo para verificar se o ponto (X, Y), ou seja o ponto testado está à esquerda semi-plano da nossa linha (é uma questão de gosto, ele também poderia ser do lado direito, mas também a direção de limites linhas tem que ser mudado em nesse caso), este é projetar o raio do ponto para a direita (ou esquerda) e reconhecer a interseção com a linha. Optamos por projeto o raio na direção horizontal (novamente, é uma questão de gosto, ele também pode ser feito em vertical com restrições semelhantes), por isso temos:
(XB-XA)
X < ------- * (Y - YA) + XA
(YB-YA)
Agora precisamos saber se o ponto está à esquerda (ou direita) lado o segmento de linha só, não todo o plano, por isso precisamos restringir a pesquisa apenas para esse segmento, mas isso é fácil, já que para estar dentro do segmento de um e apenas um ponto na linha pode ser mais elevada que Y no eixo vertical. Como se trata de uma restrição que mais forte necessidades de ser o primeiro a verificar, assim que nós tomamos primeiros somente as linhas atender esta exigência e, em seguida, verificar a sua possition. Pela Jordan Curva teorema qualquer raio projetado para uma intersecção polígono must em um mesmo número de linhas. Portanto, estamos a fazer, vamos lançar o ray ao direita e, em seguida, cada vez que cruza uma linha, alternar seu estado. No entanto, em nossa implementação estamos goint para verificar o comprimento da o saco de soluções que satisfaçam as restrições dadas e decidir o innership sobre ela. para cada linha no polígono isso tem de ser feito.
is_left_half_plane(_,[],[],_).
is_left_half_plane(X:Y,[XA:YA,XB:YB], [[X1:Y1,X2:Y2]|R], Test) :- [XA:YA, XB:YB] = [X1:Y1, X2:Y2], call(Test, X , (((XB - XA) * (Y - YA)) / (YB - YA) + XA));
is_left_half_plane(X:Y, [XA:YA, XB:YB], R, Test).
in_y_range_at_poly(Y,[XA:YA,XB:YB],Polygon) :- get_line([XA:YA,XB:YB],Polygon), in_range(Y,YA,YB).
all_in_range(Coordinate,Polygon,Lines) :- aggregate(bag(Line), in_y_range_at_poly(Coordinate,Line,Polygon), Lines).
traverses_ray(X:Y, Lines, Count) :- aggregate(bag(Line), is_left_half_plane(X:Y, Line, Lines, <), IntersectingLines), length(IntersectingLines, Count).
% This is the entry point predicate
inside_poly(X:Y,Polygon,Answer) :- all_in_range(Y,Polygon,Lines), traverses_ray(X:Y, Lines, Count), (1 is mod(Count,2)->Answer=inside;Answer=outside).
Java Versão:
public class Geocode {
private float latitude;
private float longitude;
public Geocode() {
}
public Geocode(float latitude, float longitude) {
this.latitude = latitude;
this.longitude = longitude;
}
public float getLatitude() {
return latitude;
}
public void setLatitude(float latitude) {
this.latitude = latitude;
}
public float getLongitude() {
return longitude;
}
public void setLongitude(float longitude) {
this.longitude = longitude;
}
}
public class GeoPolygon {
private ArrayList<Geocode> points;
public GeoPolygon() {
this.points = new ArrayList<Geocode>();
}
public GeoPolygon(ArrayList<Geocode> points) {
this.points = points;
}
public GeoPolygon add(Geocode geo) {
points.add(geo);
return this;
}
public boolean inside(Geocode geo) {
int i, j;
boolean c = false;
for (i = 0, j = points.size() - 1; i < points.size(); j = i++) {
if (((points.get(i).getLongitude() > geo.getLongitude()) != (points.get(j).getLongitude() > geo.getLongitude())) &&
(geo.getLatitude() < (points.get(j).getLatitude() - points.get(i).getLatitude()) * (geo.getLongitude() - points.get(i).getLongitude()) / (points.get(j).getLongitude() - points.get(i).getLongitude()) + points.get(i).getLatitude()))
c = !c;
}
return c;
}
}
.Net port:
static void Main(string[] args)
{
Console.Write("Hola");
List<double> vertx = new List<double>();
List<double> verty = new List<double>();
int i, j, c = 0;
vertx.Add(1);
vertx.Add(2);
vertx.Add(1);
vertx.Add(4);
vertx.Add(4);
vertx.Add(1);
verty.Add(1);
verty.Add(2);
verty.Add(4);
verty.Add(4);
verty.Add(1);
verty.Add(1);
int nvert = 6; //Vértices del poligono
double testx = 2;
double testy = 5;
for (i = 0, j = nvert - 1; i < nvert; j = i++)
{
if (((verty[i] > testy) != (verty[j] > testy)) &&
(testx < (vertx[j] - vertx[i]) * (testy - verty[i]) / (verty[j] - verty[i]) + vertx[i]))
c = 1;
}
}
VBA VERSÃO:
Nota: Lembre-se que se o seu polígono é uma área dentro de um mapa que Latitude / Longitude são y / x valores em oposição a X / Y (Latitude = Y, Longitude = X), devido à do que eu entendo são implicações históricas de caminho de volta quando Longitude não era uma medição.
módulo de classe: CPoint
Private pXValue As Double
Private pYValue As Double
'''''X Value Property'''''
Public Property Get X() As Double
X = pXValue
End Property
Public Property Let X(Value As Double)
pXValue = Value
End Property
'''''Y Value Property'''''
Public Property Get Y() As Double
Y = pYValue
End Property
Public Property Let Y(Value As Double)
pYValue = Value
End Property
MÓDULO:
Public Function isPointInPolygon(p As CPoint, polygon() As CPoint) As Boolean
Dim i As Integer
Dim j As Integer
Dim q As Object
Dim minX As Double
Dim maxX As Double
Dim minY As Double
Dim maxY As Double
minX = polygon(0).X
maxX = polygon(0).X
minY = polygon(0).Y
maxY = polygon(0).Y
For i = 1 To UBound(polygon)
Set q = polygon(i)
minX = vbMin(q.X, minX)
maxX = vbMax(q.X, maxX)
minY = vbMin(q.Y, minY)
maxY = vbMax(q.Y, maxY)
Next i
If p.X < minX Or p.X > maxX Or p.Y < minY Or p.Y > maxY Then
isPointInPolygon = False
Exit Function
End If
' SOURCE: http://www.ecse.rpi.edu/Homepages/wrf/Research/Short_Notes/pnpoly.html
isPointInPolygon = False
i = 0
j = UBound(polygon)
Do While i < UBound(polygon) + 1
If (polygon(i).Y > p.Y) Then
If (polygon(j).Y < p.Y) Then
If p.X < (polygon(j).X - polygon(i).X) * (p.Y - polygon(i).Y) / (polygon(j).Y - polygon(i).Y) + polygon(i).X Then
isPointInPolygon = True
Exit Function
End If
End If
ElseIf (polygon(i).Y < p.Y) Then
If (polygon(j).Y > p.Y) Then
If p.X < (polygon(j).X - polygon(i).X) * (p.Y - polygon(i).Y) / (polygon(j).Y - polygon(i).Y) + polygon(i).X Then
isPointInPolygon = True
Exit Function
End If
End If
End If
j = i
i = i + 1
Loop
End Function
Function vbMax(n1, n2) As Double
vbMax = IIf(n1 > n2, n1, n2)
End Function
Function vbMin(n1, n2) As Double
vbMin = IIf(n1 > n2, n2, n1)
End Function
Sub TestPointInPolygon()
Dim i As Integer
Dim InPolygon As Boolean
' MARKER Object
Dim p As CPoint
Set p = New CPoint
p.X = <ENTER X VALUE HERE>
p.Y = <ENTER Y VALUE HERE>
' POLYGON OBJECT
Dim polygon() As CPoint
ReDim polygon(<ENTER VALUE HERE>) 'Amount of vertices in polygon - 1
For i = 0 To <ENTER VALUE HERE> 'Same value as above
Set polygon(i) = New CPoint
polygon(i).X = <ASSIGN X VALUE HERE> 'Source a list of values that can be looped through
polgyon(i).Y = <ASSIGN Y VALUE HERE> 'Source a list of values that can be looped through
Next i
InPolygon = isPointInPolygon(p, polygon)
MsgBox InPolygon
End Sub
Eu fiz uma implementação Python de de nirg c ++ código :
Entradas
- bounding_points:. nós que compõem o polígono
-
bounding_box_positions: pontos candidatos filtrar. (Em minha implementação criado a partir da caixa delimitadora.
(As entradas são listas de tuplas no formato:
[(xcord, ycord), ...])
Retorna
- Todos os pontos que estão dentro do polígono.
def polygon_ray_casting(self, bounding_points, bounding_box_positions):
# Arrays containing the x- and y-coordinates of the polygon's vertices.
vertx = [point[0] for point in bounding_points]
verty = [point[1] for point in bounding_points]
# Number of vertices in the polygon
nvert = len(bounding_points)
# Points that are inside
points_inside = []
# For every candidate position within the bounding box
for idx, pos in enumerate(bounding_box_positions):
testx, testy = (pos[0], pos[1])
c = 0
for i in range(0, nvert):
j = i - 1 if i != 0 else nvert - 1
if( ((verty[i] > testy ) != (verty[j] > testy)) and
(testx < (vertx[j] - vertx[i]) * (testy - verty[i]) / (verty[j] - verty[i]) + vertx[i]) ):
c += 1
# If odd, that means that we are inside the polygon
if c % 2 == 1:
points_inside.append(pos)
return points_inside
Novamente, a idéia é retirado aqui
Aqui está um ponto no teste de polígono em C que não está usando raios-casting. E ele pode trabalhar para áreas (autointersecções) sobreposição, consulte o argumento use_holes.
/* math lib (defined below) */
static float dot_v2v2(const float a[2], const float b[2]);
static float angle_signed_v2v2(const float v1[2], const float v2[2]);
static void copy_v2_v2(float r[2], const float a[2]);
/* intersection function */
bool isect_point_poly_v2(const float pt[2], const float verts[][2], const unsigned int nr,
const bool use_holes)
{
/* we do the angle rule, define that all added angles should be about zero or (2 * PI) */
float angletot = 0.0;
float fp1[2], fp2[2];
unsigned int i;
const float *p1, *p2;
p1 = verts[nr - 1];
/* first vector */
fp1[0] = p1[0] - pt[0];
fp1[1] = p1[1] - pt[1];
for (i = 0; i < nr; i++) {
p2 = verts[i];
/* second vector */
fp2[0] = p2[0] - pt[0];
fp2[1] = p2[1] - pt[1];
/* dot and angle and cross */
angletot += angle_signed_v2v2(fp1, fp2);
/* circulate */
copy_v2_v2(fp1, fp2);
p1 = p2;
}
angletot = fabsf(angletot);
if (use_holes) {
const float nested = floorf((angletot / (float)(M_PI * 2.0)) + 0.00001f);
angletot -= nested * (float)(M_PI * 2.0);
return (angletot > 4.0f) != ((int)nested % 2);
}
else {
return (angletot > 4.0f);
}
}
/* math lib */
static float dot_v2v2(const float a[2], const float b[2])
{
return a[0] * b[0] + a[1] * b[1];
}
static float angle_signed_v2v2(const float v1[2], const float v2[2])
{
const float perp_dot = (v1[1] * v2[0]) - (v1[0] * v2[1]);
return atan2f(perp_dot, dot_v2v2(v1, v2));
}
static void copy_v2_v2(float r[2], const float a[2])
{
r[0] = a[0];
r[1] = a[1];
}
Nota: este é um dos métodos menos ideal, uma vez que inclui uma série de chamadas para atan2f, mas pode ser de interesse para desenvolvedores que lêem este thread (nos meus testes o seu ~ 23x mais lento, em seguida, usando o método de intersecção linha).
Para detectar bater em Polygon precisamos testar duas coisas:
- Se Point é área de polígono dentro. (Pode ser conseguida por Raio-Fundição Algorithm)
- Se Ponto está na fronteira polígono (pode ser conseguida por mesmo algoritmo que é usado para a detecção do ponto em poligonal (linha)).
Para lidar com os seguintes casos especiais em Ray :
- O raio sobrepõe um dos lados do polígono.
- O ponto está dentro do polígono e o raio passa através de um vértice do polígono.
- O ponto está fora do polígono eo ray apenas toca uma das ângulo do polígono.
Verifique determinar se um ponto está dentro de um complexo polígono . O artigo fornece uma maneira fácil de resolvê-los assim não haverá nenhum tratamento especial necessário para os casos acima.
Você pode fazer isso verificando se a área formada por liga o ponto desejado para os vértices de seu polígono corresponde à área do próprio polígono.
Ou você pode verificar se a soma dos ângulos internos de seu ponto para cada par de dois vértices do polígono consecutivos ao seu ponto de verificação somas a 360, mas tenho a sensação de que a primeira opção é mais rápido porque ele não envolve divisões nem cálculos do inverso de funções trigonométricas.
Eu não sei o que acontece se o seu polígono tem um buraco dentro dela, mas parece-me que a idéia principal pode ser adaptado a esta situação
Você pode também postar a pergunta em uma comunidade matemática. Aposto que eles têm um milhão de maneiras de fazer isso
Se você está procurando uma biblioteca javascript há um google javascript mapeia extensão v3 para a classe Polygon para detectar ou não a reside ponto dentro dele.
var polygon = new google.maps.Polygon([], "#000000", 1, 1, "#336699", 0.3);
var isWithinPolygon = polygon.containsLatLng(40, -90);
Ao usar qt (Qt 4.3 +) , pode-se usar a função de QPolygon containsPoint
A resposta depende se você tem o simples ou polígonos complexos. polígono simples não deve ter quaisquer cruzamentos segmento de linha. Assim, eles podem ter os buracos, mas as linhas não podem cruzar entre si. regiões complexas podem ter os cruzamentos de linha -. para que possam ter as regiões sobrepostas, ou regiões que se tocam apenas por um único ponto
Para polígono simples o melhor algoritmo é fundição Ray (número Travessia) algoritmo. Para polígonos complexos, este algoritmo não detectar pontos que estão dentro das regiões sobrepostas. Assim, para polígonos complexos você tem que usar Winding algoritmo número.
Aqui está um excelente artigo com a implementação C de ambos os algoritmos. Eu tentei-los e eles funcionam bem.
versão Scala de solução por nirg (assume delimitadora rectângulo pré-verificação é feita separadamente):
def inside(p: Point, polygon: Array[Point], bounds: Bounds): Boolean = {
val length = polygon.length
@tailrec
def oddIntersections(i: Int, j: Int, tracker: Boolean): Boolean = {
if (i == length)
tracker
else {
val intersects = (polygon(i).y > p.y) != (polygon(j).y > p.y) && p.x < (polygon(j).x - polygon(i).x) * (p.y - polygon(i).y) / (polygon(j).y - polygon(i).y) + polygon(i).x
oddIntersections(i + 1, i, if (intersects) !tracker else tracker)
}
}
oddIntersections(0, length - 1, tracker = false)
}
Aqui é a versão golang de resposta @nirg (inspirado por código C # por @@ m-Katz)
func isPointInPolygon(polygon []point, testp point) bool {
minX := polygon[0].X
maxX := polygon[0].X
minY := polygon[0].Y
maxY := polygon[0].Y
for _, p := range polygon {
minX = min(p.X, minX)
maxX = max(p.X, maxX)
minY = min(p.Y, minY)
maxY = max(p.Y, maxY)
}
if testp.X < minX || testp.X > maxX || testp.Y < minY || testp.Y > maxY {
return false
}
inside := false
j := len(polygon) - 1
for i := 0; i < len(polygon); i++ {
if (polygon[i].Y > testp.Y) != (polygon[j].Y > testp.Y) && testp.X < (polygon[j].X-polygon[i].X)*(testp.Y-polygon[i].Y)/(polygon[j].Y-polygon[i].Y)+polygon[i].X {
inside = !inside
}
j = i
}
return inside
}
surpreso que ninguém trouxe isso antes, mas para os pragmáticos que requerem um banco de dados:. MongoDB tem excelente suporte para Geo consulta incluindo este
O que você está procurando é:
db.neighborhoods.findOne ({geometria: {$ geoIntersects: {$ geometria: { digite: "Point", coordenadas: [ "longitude", "latitude"]}}} })
Neighborhoods é a coleção que armazena um ou mais polígonos em formato GeoJSON padrão. Se a consulta retorna nulo não é cortada caso contrário ele é.
Muito bem documentado aqui: https://docs.mongodb.com/manual/tutorial/geospatial-tutorial/
O desempenho por mais de 6.000 pontos classificados em uma grade polígono irregular 330 foi inferior a um minuto sem otimização em tudo e incluindo o tempo para atualizar documentos com suas respectivas polígono.
