Подтвердить что ты не робот

Как определить, находится ли 2D-точка внутри многоугольника?

Я пытаюсь создать быструю двумерную точку внутри алгоритма многоугольника, для использования при хит-тестировании (например, Polygon.contains(p:Point)). Предложения по эффективным методам будут оценены.

4b9b3361

Ответ 1

Для графики я предпочел бы не использовать целые числа. Многие системы используют целые числа для рисования UI (пикселы - это ints), но macOS, например, использует float для всего. macOS знает только точки и точка может перевести на один пиксель, но в зависимости от разрешения монитора он может перевести что-то еще. На сетчатых экранах половина точки (0,5/0,5) - это пиксель. Тем не менее, я никогда не замечал, что пользовательские интерфейсы MacOS значительно медленнее, чем другие пользовательские интерфейсы. После того, как все 3D API (OpenGL или Direct3D) также работают с поплавками, а современные графические библиотеки очень часто используют ускорение GPU.

Теперь вы сказали, что скорость - это ваша главная забота, хорошо, отпустите скорость. Прежде чем запускать какой-либо сложный алгоритм, сначала выполните простой тест. Создайте ось с выравниванием по оси вокруг вашего многоугольника. Это очень просто, быстро и уже может вас сэкономить много расчетов. Как это работает? Перейдем по всем точкам многоугольника и найдем значения min/max для X и Y.

например. у вас есть точки (9/1), (4/3), (2/7), (8/2), (3/6). Это означает, что Xmin равно 2, Xmax равно 9, Ymin равно 1, а Ymax равно 7. Точка вне прямоугольника с двумя ребрами (2/1) и (9/7) не может находиться внутри многоугольника.

// 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!
}

Это первый тест для любой точки. Как вы можете видеть, этот тест очень быстрый, но он также очень грубый. Чтобы обрабатывать точки, находящиеся внутри ограничивающего прямоугольника, нам нужен более сложный алгоритм. Есть несколько способов, как это можно вычислить. Какой метод работает также зависит от того, может ли многоугольник иметь отверстия или всегда будет твердым. Вот примеры твердых (один выпуклый, один вогнутый):

Polygon without hole

И вот один с отверстием:

Polygon with hole

У зеленого есть отверстие посередине!

Самый простой алгоритм, который может обрабатывать все три случая выше и по-прежнему довольно быстро, называется ray casting. Идея алгоритма довольно проста: нарисуйте виртуальный луч из любой точки вне полигона до своей точки и подсчитайте, как часто он попадает на сторону многоугольника. Если количество обращений четное, это вне полигона, если оно нечетное, оно внутри.

Demonstrating how the ray cuts through a polygon

Алгоритм число обмоток был бы альтернативой, он более точным для точек, очень близких к линии многоугольников, но также намного медленнее. Рэй-кастинг может оказаться неудачным для точек, слишком близких к многоугольной стороне из-за ограниченной точности с плавающей запятой и вопросов округления, но на самом деле это вряд ли проблема, так как точка зрения находится близко к стороне, она часто визуально даже не возможна для зрителю, чтобы узнать, находится ли он уже внутри или снаружи.

У вас все еще есть ограничивающая рамка сверху, помните? Просто выберите точку за пределами рамки и используйте ее в качестве отправной точки для своего луча. Например. точка (Xmin - e/p.y) находится за пределами многоугольника.

Но что такое e? Ну, e (фактически epsilon) дает ограничивающее поле некоторую отступность. Как я уже сказал, трассировка лучей терпит неудачу, если мы начнем слишком близко к линии многоугольника. Так как ограничивающая рамка может равняться многоугольнику (если многоугольник представляет собой прямоугольник, выровненный по оси, ограничивающий прямоугольник равен самому многоугольнику!), Нам нужно некоторое дополнение, чтобы сделать это безопасным, все. Насколько большой вы должны выбрать e? Не слишком большой. Это зависит от масштаба системы координат, который вы используете для рисования. Если ширина шага пикселя равна 1,0, то просто выберите 1.0 (еще 0,1 работал бы также)

Теперь, когда у нас есть луч с его начальными и конечными координатами, проблема сдвигается от "является точкой внутри полигона" до "как часто луч пересекает сторону многоугольника". Поэтому мы не можем просто работать с точками полигона по-прежнему, теперь нам нужны фактические стороны. Сторона всегда определяется двумя точками.

side 1: (X1/Y1)-(X2/Y2)
side 2: (X2/Y2)-(X3/Y3)
side 3: (X3/Y3)-(X4/Y4)
:

Вам нужно протестировать луч со всех сторон. Рассмотрим луч как вектор, а каждая сторона - вектор. Луч должен ударить каждую сторону ровно один раз или никогда вообще. Он не может ударить по одной и той же стороне дважды. Две строки в 2D-пространстве всегда будут пересекаться ровно один раз, если они не параллельны, и в этом случае они никогда не пересекаются. Однако, поскольку векторы имеют ограниченную длину, два вектора могут быть не параллельными и никогда не пересекаться, потому что они слишком коротки, чтобы встречаться друг с другом.

// 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
}

До сих пор так хорошо, но как вы проверяете, пересекаются ли два вектора? Здесь некоторый C-код (не тестировался), который должен делать трюк:

#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 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;
}

Входными значениями являются две конечные точки вектора 1 (v1x1/v1y1 и v1x2/v1y2) и вектор 2 (v2x1/v2y1 и v2x2/v2y2). Таким образом, у вас есть 2 вектора, 4 точки, 8 координат. YES и NO понятны. YES увеличивает пересечения, NO ничего не делает.

Как насчет COLLINEAR? Это означает, что оба вектора лежат на одной и той же бесконечной линии, в зависимости от положения и длины, они не пересекаются вообще или пересекаются в бесконечном числе точек. Я не совсем уверен, как справиться с этим делом, я бы не стал считать его пересечением в любом случае. Ну, этот случай на практике довольно редок в любом случае из-за ошибок округления с плавающей запятой; лучший код, вероятно, не будет тестировать для == 0.0f, а вместо чего для чего-то вроде < epsilon, где epsilon - довольно небольшое число.

Если вам нужно проверить большее количество очков, вы можете, конечно, немного ускорить все это, сохранив линейные уравнения стандартных форм сторон многоугольника в памяти, поэтому вам не нужно их пересчитывать каждый раз. Это сэкономит вам два умножения с плавающей запятой и три вычитания с плавающей запятой для каждого теста в обмен на сохранение трех значений с плавающей запятой на каждую сторону полигона в памяти. Это типичная ошибка памяти и время вычисления.

И последнее, но не менее важное: если вы можете использовать 3D-оборудование для решения проблемы, есть интересная альтернатива. Просто позвольте GPU делать всю работу за вас. Создайте поверхность для рисования, которая выключена. Заполните его полностью черным цветом. Теперь пусть OpenGL или Direct3D нарисуют ваш полигон (или даже все ваши полигоны, если вы просто хотите проверить, находится ли точка в любом из них, но вам все равно, какой) и заполнить многоугольник другим цвет, например белый. Чтобы проверить, находится ли точка внутри полигона, получите цвет этой точки с поверхности чертежа. Это всего лишь выборка памяти O (1).

Конечно, этот метод можно использовать только в том случае, если ваша поверхность рисования не должна быть огромной. Если он не может вписаться в память графического процессора, этот метод работает медленнее, чем при работе с ЦП. Если это должно быть огромным, и ваш GPU поддерживает современные шейдеры, вы все равно можете использовать графический процессор, реализуя лучевое кастинг, показанный выше, как графический шейдер, что абсолютно возможно. Для большего количества полигонов или большого количества точек для тестирования это окупится, подумайте, что некоторые графические процессоры смогут протестировать от 64 до 256 точек параллельно. Обратите внимание, однако, что передача данных с CPU на GPU и обратно всегда стоит дорого, поэтому для простого тестирования пары точек против нескольких простых полигонов, где точки или полигоны являются динамическими и будут часто меняться, подход GPU редко будет платить выкл.

Ответ 2

Я думаю, что следующий фрагмент кода является наилучшим решением (взят из здесь):

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;
}

Аргументы

  • nvert: количество вершин в многоangularьнике. Повторять ли первую вершину в конце, обсуждали в статье, упомянутой выше.
  • vertx, verty: массивы, содержащие x- и y-координаты вершин многоangularьника.
  • testx, testy: x- и y-координата контрольной точки.

Он короткий и эффективный и работает как для выпуклых, так и для вогнутых многоangularьников. Как предлагалось ранее, вы должны сначала проверить ограничивающий прямоangularьник и обработать отверстия многоangularьника отдельно.

Идея этого довольно проста. Автор описывает это следующим образом:

Я запускаю полубесконечный луч по горизонтали (с увеличением x, фиксированный y) из контрольной точки и подсчитываю, сколько ребер он пересекает. На каждом пересечении луч переключается между внутренним и внешним. Это называется теорема кривой Джордана.

Переменная с переключается с 0 на 1 и с 1 на 0 каждый раз, когда горизонтальный луч пересекает любое ребро. Так что в основном он отслеживает, является ли количество пересеченных ребер четным или нечетным. 0 означает четное, а 1 означает нечетное.

Ответ 3

Вот версия С# ответа предоставленного nirg, которая исходит от этого профессора RPI. Обратите внимание, что использование кода из этого источника RPI требует атрибуции.

В верхней части была добавлена ​​проверка рамки. Тем не менее, как отмечает Джеймс Браун, основной код почти так же быстро, как и ограничивающий блок, сам проверяет, поэтому проверка ограничивающего блока может фактически замедлить общую операцию в том случае, если большинство точек, которые вы проверяете, находятся в ограничивающей рамке, Таким образом, вы можете оставить контрольную коробку, или альтернативой было бы прекомпотировать ограничивающие поля ваших полигонов, если они не меняют форму слишком часто.

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;
}

Ответ 4

Вот вариант ответа М. Каца на JavaScript, основанный на подходе Нирга:

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;
}

Ответ 5

Вычислить ориентированную сумму углов между точкой p и каждой из вершин многоугольника. Если общий ориентированный угол равен 360 градусам, точка внутри. Если сумма равна 0, точка находится вне.

Мне нравится этот метод лучше, потому что он более устойчив и менее зависит от числовой точности.

Методы, которые вычисляют четность числа пересечений, ограничены, потому что вы можете "ударить" вершину при вычислении числа пересечений.

EDIT: By The Way, этот метод работает с вогнутыми и выпуклыми многоугольниками.

EDIT: Недавно я нашел целую статью Википедии по этой теме.

Ответ 6

статья Эрика Хейнса, приведенная bobobobo, действительно отличная. Особенно интересны таблицы, сравнивающие производительность алгоритмов; метод суммирования углов действительно плох по сравнению с другими. Также интересно, что оптимизация, например, использование сетки поиска для дальнейшего разбиения полигона на "in" и "out" сектора, может сделать тест невероятно быстрым даже на многоугольниках s > 1000 сторонами.

В любом случае, это в первые дни, но мой голос идет на метод "пересечений", что в значительной степени отражает то, что описывает Мекки. Однако я нашел его наиболее сугубо описанный и кодифицированный Дэвидом Бурком. Мне нравится, что нет реальной тригонометрии, и она работает на выпуклую и вогнутую, и она достаточно хорошо работает, когда увеличивается количество сторон.

Кстати, вот одна из таблиц производительности из статьи Эрика Хейнса для интереса, тестирование на случайных полигонах.

                       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

Ответ 7

Этот вопрос так интересен. У меня есть другая работающая идея, отличная от других ответов этого сообщения.

Идея состоит в том, чтобы использовать сумму углов, чтобы решить, находится ли цель внутри или снаружи. Если цель находится внутри области, сумма формы угла по цели и каждые две граничные точки будет равна 360. Если цель находится снаружи, сумма не будет равна 360. Угол имеет направление. Если угол идет назад, угол является отрицательным. Это точно так же, как расчет числа введите описание изображения здесь

Мой алгоритм предполагает, что по часовой стрелке есть положительное направление. Вот потенциальный ввод:

[[-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]]

Ниже приведен код python, который реализует идею:

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

Ответ 8

Быстрая версия ответа by nirg:

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
    }
}

Ответ 9

Очень нравится решение, опубликованное Nirg и отредактированное bobobobo. Я только сделал его дружественным к JavaScript и немного более разборчивым для моего использования:

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;
}

Ответ 10

Я немного поработал над этим, когда я был исследователем в Майкл Стоунбрейкер - вы знаете, профессор, который придумал Ingres, PostgreSQL и т.д.

Мы поняли, что самый быстрый способ состоял в том, чтобы сначала сделать ограничивающий прямоугольник, потому что он СУПЕР быстр. Если он находится за пределами рамки, он находится снаружи. В противном случае вы делаете более трудную работу...

Если вам нужен отличный алгоритм, посмотрите на исходный код PostgreSQL с открытым исходным кодом для гео-работы...

Я хочу указать, что мы никогда не понимаем права на левую руку (также выражаемая как проблема "внутри" и "снаружи"...


ОБНОВЛЕНИЕ

Ссылка BKB обеспечила большое количество разумных алгоритмов. Я работал над проблемами науки о Земле и поэтому нуждался в решении, которое работает в широте/долготе, и у него есть своеобразная проблема ручности - это область внутри меньшей площади или большей площади? Ответ заключается в том, что "направление" вершин имеет значение - оно либо левое, либо правое, и таким образом вы можете указать любую область как "внутри" любого заданного полигона. Таким образом, моя работа использовала решение, перечисленное на этой странице.

Кроме того, моя работа использовала отдельные функции для тестов "на линии".

... Так как кто-то спросил: мы выяснили, что тесты с ограничивающей коробкой лучше всего, когда количество вершин выходит за рамки некоторого числа - сделайте очень быстрый тест, прежде чем выполнять более длительный тест, если это необходимо... Ограничительная коробка создается просто беря наибольшее x, наименьшее x, наибольшее y и наименьшее y и складывая их вместе, чтобы сделать четыре точки окна...

Еще один совет для следующих: мы сделали все наши более сложные и "световые тускнеющие" вычисления в пространстве сетки все в положительных точках на плоскости, а затем снова проецировались обратно в "настоящую" долготу/широту, избегая при этом возможные ошибки обертывания вокруг одной перекрестной линии 180 долготы и при обработке полярных областей. Отлично работает!

Ответ 11

Ответ Дэвид Сегонд в значительной степени является стандартным общим ответом, а Ричард Т - наиболее распространенная оптимизация, хотя некоторые другие. Другие сильные оптимизации основаны на менее общих решениях. Например, если вы собираетесь проверять один и тот же полигон с большим количеством точек, триангуляция полигона может ускорить работу, так как существует множество очень быстрых алгоритмов поиска TIN. Другим является то, что многоугольник и точки находятся на ограниченной плоскости при низком разрешении, например, на экране, вы можете рисовать многоугольник в буфер отображения отображаемой памяти в заданном цвете и проверять цвет данного пикселя, чтобы увидеть, в многоугольниках.

Как и многие оптимизации, они основаны на конкретных, а не общих случаях, и дают преимущества, основанные на амортизированном времени, а не на одном использовании.

Работая в этом поле, я нашел полезную геометрию вычислений Джозефа О'Руркса в C 'ISBN 0-521-44034-3.

Ответ 12

Тривиальное решение состоит в том, чтобы разделить многоугольник на треугольники и нажать тест треугольников, как описано здесь

Если ваш многоугольник CONVEX, возможно, существует лучший подход. Посмотрите на многоугольник как совокупность бесконечных линий. Каждая строка разделяет пространство на две. для каждой точки легко сказать, находится ли она с одной или с другой стороны линии. Если точка находится на одной стороне всех линий, она находится внутри многоугольника.

Ответ 13

Я тоже считал, что суммирование суммируется только для выпуклых многоугольников, но это неверно.

На этом сайте есть хорошая диаграмма, показывающая именно это, и хорошее объяснение при тестировании ударов: Gamasutra - Сбой в Новый год: обнаружение столкновений

Ответ 14

Я понимаю, что это старо, но вот алгоритм кастинга лучей, реализованный в Cocoa, в случае, если кто-то заинтересован. Не уверен, что это самый эффективный способ сделать что-то, но это может помочь кому-то.

- (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;
}

Ответ 15

Obj-C версия ответа nirg с образцом метода для тестирования точек. Ответ Нирга работал хорошо для меня.

- (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");
    }
}

sample polygon

Ответ 16

версия ответа nirg здесь: я просто передам код. Это может сэкономить некоторое время.

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;
        }

Ответ 17

Нет ничего более прекрасного, чем индуктивное определение проблемы. Для полноты здесь у вас есть версия в прологе, которая также может прояснить хотя бы за литье лучей:

На основе моделирования алгоритма простоты в http://www.ecse.rpi.edu/Homepages/wrf/Research/Short_Notes/pnpoly.html

Некоторые вспомогательные предикаты:

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]).

Уравнение прямой, заданной 2 точками A и B (Линия (A, B)):

                    (YB-YA)
           Y - YA = ------- * (X - XA) 
                    (XB-YB) 

Важно, что направление вращения для линии установленный по часовой стрелке для границ и против часовой стрелки для отверстий. Мы будем проверять, находится ли точка (X, Y), т.е. тестируемая точка слева полуплоскость нашей линии (это вопрос вкуса, это также может быть правая сторона, но также направление линий границ должно быть изменено в в этом случае), это проецирование луча из точки вправо (или влево) и подтвердите пересечение с линией. Мы выбрали проект луч в горизонтальном направлении (опять же это вопрос вкуса, это также можно было бы сделать вертикально с аналогичными ограничениями), поэтому мы имеем:

               (XB-XA)
           X < ------- * (Y - YA) + XA
               (YB-YA) 

Теперь нам нужно знать, находится ли точка слева (или справа) только сегмент линии, а не всю плоскость, поэтому нам нужно ограничьте поиск только этим сегментом, но это легко, так как быть внутри сегмента только одна точка в линии может быть выше чем Y по вертикальной оси. Поскольку это более сильное ограничение, должен быть первым, чтобы проверить, поэтому сначала берем только те линии выполнив это требование, а затем проверит его принадлежность. Иорданией Теорема о кривых любой луч, проецируемый на многоугольник, должен пересекаться в точке даже количество строк. Итак, мы закончили, мы выбросим луч в справа, а затем каждый раз, когда он пересекает линию, переключает ее состояние. Однако в нашей реализации мы собираемся проверить длину пакет решений, отвечающих заданным ограничениям, и взаимозависимость. для каждой линии в многоугольнике это должно быть сделано.

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).

Ответ 18

Версия Java:

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;
    }

}

Ответ 19

.Net-порт:

    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;
        }
    }

Ответ 20

VBA VERSION:

Примечание. Помните, что если ваш многоугольник является областью в карте, ширина/долгота - это значения Y/X, а не X/Y (Latitude = Y, Longitude = X) из-за того, что я понимаю, являются историческими последствиями из назад, когда долгота не была измерением.

МОДУЛЬ КЛАССА: 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

МОДУЛЬ:

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

Ответ 21

Я сделал реализацию Python из nirg c++ код:

входные

  • bounding_points: узлы, составляющие многоугольник.
  • bounding_box_positions: кандидат указывает на фильтрацию. (В моей реализации создано из ограничивающей рамки.

    (Входы представляют собой списки кортежей в формате: [(xcord, ycord),...])

Возвращает

  • Все точки, находящиеся внутри многоугольника.
def polygon_ray_casting(self, bounding_points, bounding_box_positions):
    # Arrays containing the x- and y-coordinates of the polygon 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

Опять же, идея взята здесь

Ответ 22

Вот точка в полигоне в C, которая не использует лучевое литье. И он может работать для перекрывающихся областей (самопересечений), см. Аргумент 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];
}

Примечание. Это один из менее оптимальных методов, поскольку он включает в себя множество вызовов atan2f, но может быть интересным для разработчиков, читающих этот поток (в моих тестах его ~ 23x медленнее, а затем с использованием метода пересечения строк).

Ответ 23

Для обнаружения удара по полигону нам нужно проверить две вещи:

  • Если точка находится внутри области полигона. (может выполняться алгоритмом Ray-Casting)
  • Если точка находится на границе многоугольника (может выполняться по тому же алгоритму, который используется для определения точки на полилинии (линия)).

Ответ 24

Чтобы справиться со следующими особыми случаями в алгоритм преобразования лучей:

  • Лук перекрывает одну из сторон многоугольника.
  • Точка находится внутри многоугольника, а луч проходит через вершину многоугольника.
  • Точка находится за пределами многоугольника, и луч просто касается одного из углов многоугольника.

Проверить Определить, находится ли точка внутри сложного многоугольника. В статье представлен простой способ их устранения, поэтому для вышеуказанных случаев не требуется специального лечения.

Ответ 25

Вы можете сделать это, проверив, соответствует ли область, образованная путем подключения нужной точки к вершинам вашего многоугольника, к области самого полигона.

Или вы могли бы проверить, равна ли сумма внутренних углов от вашей точки к каждой паре двух последовательных вершин многоугольников к вашей контрольной точке до 360, но у меня есть ощущение, что первый вариант быстрее, потому что он не включает делениям или вычислениям инверсных тригонометрических функций.

Я не знаю, что произойдет, если ваш полигон имеет отверстие внутри него, но мне кажется, что основная идея может быть адаптирована к этой ситуации.

Вы также можете задать вопрос в математическом сообществе. Уверен, у них есть миллион способов сделать это.

Ответ 26

Если вы ищете библиотеку java- script, там есть расширение javascript google maps v3 для класса Polygon для определения того, находится ли в нем точка.

var polygon = new google.maps.Polygon([], "#000000", 1, 1, "#336699", 0.3);
var isWithinPolygon = polygon.containsLatLng(40, -90);

Google Extention Github

Ответ 27

При использовании (Qt 4.3 +), можно использовать функцию QPolygon containsPoint

Ответ 28

Ответ зависит от того, есть ли у вас простые или сложные полигоны. У простых многоугольников не должно быть пересечений линий. Таким образом, они могут иметь отверстия, но линии не могут пересекать друг друга. Сложные области могут иметь пересечения линий, поэтому они могут иметь перекрывающиеся области или области, которые касаются друг друга только одной точкой.

Для простых многоугольников лучшим алгоритмом является алгоритм преобразования лучей (Crossing number). Для сложных полигонов этот алгоритм не обнаруживает точек, находящихся внутри перекрывающихся областей. Итак, для сложных полигонов вам нужно использовать алгоритм числа Winding.

Вот отличная статья с реализацией C обоих алгоритмов. Я попробовал их, и они хорошо работают.

http://geomalgorithms.com/a03-_inclusion.html

Ответ 29

Scala версия решения по nirg (предполагается, что предварительная проверка ограничивающего прямоугольника выполняется отдельно):

  def внутри (p: Point, polygon: Array [Point], bounds: Bounds): Boolean = {
 val length = polygon.length
 @tailrec def oddIntersections (i: Int, j: Int, tracker: Boolean): Boolean = {   если (i = длина)     трекер   else {     val пересекает = (многоугольник (i).y > p.y)!= (многоугольник (j).y > p.y) && & p.x < (многоугольник (j).x - многоугольник (i).x) * (p.y - polygon (i).y)/(многоугольник (j).y - многоугольник (i).y) + многоугольник (i).x     oddIntersections (i + 1, i, if (пересекает)! tracker else tracker)   } }
 oddIntersections (0, length - 1, tracker = false)
}
Код>

Ответ 30

Вот голанская версия ответа @nirg (вдохновленная кодом С# на @@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
}