быстрый предикат геометрической близости

Question

У меня есть 3 точки (A, B и X) и расстояние (d). Мне нужно сделать функцию, которая проверяет, находится ли точка X ближе, чем расстояние d к любой точке на отрезке AB.

Вопрос, во-первых, правильное ли мое решение, а затем найти лучшее (более быстрое) решение.

Мой первый проход выглядит следующим образом

AX = X-A
BX = X-B
AB = A-B

    // closer than d to A (done squared to avoid needing to compute the sqrt in mag)
If d^2 > AX.mag^2  return true

    // closer than d to B
If d^2 > BX.mag^2 return true

    // "beyond"  B
If (dot(BX,AB) < 0) return false

    // "beyond"  A
If (dot(AX,AB) > 0) return false

    // find component of BX perpendicular to AB
Return (BX.mag)^2 - (dot(AB,BX)/AB.mag)^2 < d^2

Этот код будет выполняться для большого набора P и большого набора триплетов A / B / d с целью найти все P, которые проходят по крайней мере для одного A / B / d, поэтому я подозреваю, что это способ снизить общую стоимость, основанную на этом, но я еще не рассматривал это.

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

---

РЕДАКТИРОВАТЬ: это двумерная проблема (но решение, которое обобщается на 3D, было бы круто

Редактировать: при дальнейшем осмыслении я ожидаю, что коэффициент попадания будет около 50% Точка X может быть сформирована во вложенной иерархии, поэтому я ожидаю, что смогу обрезать большие поддеревья как все проходные и все неудачные. A / B / d, ограничивающий триплеты, будет больше хитростью.

Редактировать: d в том же порядке, что и AB.

---

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

• Первый бит Артелиуса, предварительно рассчитать матрицу, в которой AB будет располагаться по центру в начале координат и выровнено по оси X. (2 добавления, 4 муль и 2 добавления)
• сложить все в 1-й квадрант (2 абс)
• масштаб в X & Y, чтобы сделать центральную часть зоны единым квадратом (2 мул)
• проверить, находится ли точка в этом квадрате (2 теста), так что выйти
• проверить заглушку (вернуться к немасштабированным значениям
  • перевести в х, чтобы поместить конец в начало координат (1 добавление)
  • квадрат и сложение (2 муль, 1 сложение)
  • сравните с d ^ 2 (1 cmp)

Я вполне уверен, что это лучше моего решения.

(если ничего лучше не приходит, Артелиус получает «приз»)

Answer 1

Если ваш набор (A, B, d) фиксирован, вы можете рассчитать пару матриц для каждой, чтобы перевести систему координат, чтобы линия AB стала осью X, а средняя точка AB была происхождение.

Я думаю это простой способ построения матриц:

trans = - ((A + B) / 2)        // translate midpoint of AB to origin

rot.col1 = AB / AB.mag         // unit vector in AB direction

                        0 -1    
rot.col2 = rot.col1 * (      ) // unit vector perp to AB
                        1  0 

rot = rot.inverse()            // but it needs to be done in reverse

Тогда вы просто берете каждый X и делаете rot * (X + trans).

Рассматриваемая область теперь фактически симметрична как по осям x, так и по оси y, поэтому вы можете взять абсолютное значение координаты x и координаты y.

Тогда вам просто нужно проверить:

y < d && x < AB.mag/2            //"along" the line segment
|| (x - AB.mag/2)^2 + y^2 < d^2  // the "end cap".

Вы можете сделать еще один трюк; матрица может уменьшаться с коэффициентом AB.mag/2 (помните, что матрицы рассчитываются только один раз за (A, B), что означает, что лучше их найти медленнее, чем сама проверка) Это означает, что ваш чек становится:

y < 2*d/AB.mag && x < 1
|| (x - 1)^2 + y^2 < (2*d/AB.mag)^2

Заменив два экземпляра AB.mag / 2 на константу 1, это может быть касание быстрее. И, конечно же, вы можете предварительно рассчитать 2*d/AB.mag и (2*d/AB.mag)^2.

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

Answer 2

Если я правильно прочитал, то это почти то же самое, что и классический тест пересечения луча / сферы, используемый в 3D-трассировке лучей.

В этом случае у вас есть сфера в точке X радиуса d, и вы пытаетесь выяснить, пересекает ли линия AB сферу. Единственное отличие от трассировки лучей состоит в том, что в этом случае у вас есть конкретная линия AB, тогда как при трассировке лучей луч обычно обобщается как origin + distance * direction, и вам все равно, как далеко вдоль бесконечной линии AB+ это так.

В псевдокоде от моего собственного трассировщика лучей (на основе алгоритма, приведенного в «Введение в трассировку лучей» (ред. Гласснер):

Vector v = X - A
Vector d = normalise(B - A)  // unit direction vector of AB
double b = dot(v, B - A)
double discrim = b^2 - dot(v, v) + d^2
if (discrim < 0)
     return false            // definitely no intersection

Если у вас так далеко, то есть некоторый шанс, что ваше условие выполнено. Вы просто должны выяснить, находится ли пересечение (я) на линии AB:

discrim = sqrt(discrim)
double t2 = b + discrim
if (t2 <= 0)
    return false             // intersection is before A

double t1 = b  - discrim

result = (t1 < length(AB) || (t2 < length(AB))

Answer 3

Хммммммм ... Какова частота попаданий? Как часто точка «Х» соответствует требованиям близости?

Я думаю, что ваш существующий алгоритм хорош, и любые дополнительные оптимизации придут от настройки на реальные данные. Например, если точка «X» соответствует тесту на близость 99% времени, то я думаю, что ваша стратегия оптимизации должна быть значительно иной, чем если бы она проходила тест только 1% времени.

---

Кстати, когда вы дойдете до точки, где вы запускаете этот алгоритм с тысячами точек, вы должны организовать все точки в K-мерное дерево (или KDTree ). Это значительно упрощает расчет "ближайшего соседа".

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

Answer 4

Если количество наборов A / B / d велико и вы определенно находитесь в 2D, рассмотрите возможность использования R-деревьев (или их восьмиугольных эквивалентов), где каждая запись в R-дереве является минимальной ограничительной рамкой тройки A / B / d. Это позволит вам избавиться от необходимости тестировать множество тройных A / B / d и сфокусировать мощность вашего процессора только на тех немногих, где каждая точка X находится в ограничивающих прямоугольниках тройки A / B / d. Затем проведите более подробный тест, подобный тому, который вы упомянули.

Answer 5

Этот код будет выполняться для большого набора P и большого набора триплетов A / B / d с целью найти все P, которые проходят по крайней мере для одного A / B / d, поэтому я подозреваю, что это способ снизить общую стоимость, основанную на этом, но я еще не рассматривал это.

В случае, когда d ~ AB, для данной точки X вы можете сначала проверить, принадлежит ли X одной из множества сфер радиуса d и центру Ai или Bi. Посмотрите на картинку:

     ......        .....
  ...........   ...........
 ...........................
.......A-------------B.......
 ...........................
  ...........   ...........
     .....         .....

Первые два теста

If d^2 > AX.mag^2 return true
If d^2 > BX.mag^2 return true

являются самыми быстрыми, и, если d ~ AB, они также имеют наибольшую вероятность успеха (с учетом того, что тест пройден вообще). Учитывая X, вы можете сначала выполнить все «тесты сферы», а затем все последние:

Return (BX.mag)^2 - (dot(AB,BX)/AB.mag)^2