Какой самый эффективный способ найти барицентрические координаты?

45

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

Это следует из метода Ширли , где вы вычисляете площадь треугольников, образованных путем вложения точки P внутри треугольника.

Bary

Код:

Vector Triangle::getBarycentricCoordinatesAt( const Vector & P ) const
{
  Vector bary ;

  // The area of a triangle is 
  real areaABC = DOT( normal, CROSS( (b - a), (c - a) )  ) ;
  real areaPBC = DOT( normal, CROSS( (b - P), (c - P) )  ) ;
  real areaPCA = DOT( normal, CROSS( (c - P), (a - P) )  ) ;

  bary.x = areaPBC / areaABC ; // alpha
  bary.y = areaPCA / areaABC ; // beta
  bary.z = 1.0f - bary.x - bary.y ; // gamma

  return bary ;
}

Этот метод работает, но я ищу более эффективный!

bobobobo
источник
2
Помните, что наиболее эффективные решения могут быть наименее точными.
Питер Тейлор
Я предлагаю вам сделать модульный тест, чтобы вызвать этот метод ~ 100k раз (или что-то подобное) и измерить производительность. Вы можете написать тест, который гарантирует, что оно меньше некоторого значения (например, 10 с), или вы можете использовать его просто для сравнения старой и новой реализации.
ashes999

Ответы:

54

Переписано из « Обнаружения столкновений в реальном времени» Кристера Эриксона (что, кстати, является отличной книгой):

// Compute barycentric coordinates (u, v, w) for
// point p with respect to triangle (a, b, c)
void Barycentric(Point p, Point a, Point b, Point c, float &u, float &v, float &w)
{
    Vector v0 = b - a, v1 = c - a, v2 = p - a;
    float d00 = Dot(v0, v0);
    float d01 = Dot(v0, v1);
    float d11 = Dot(v1, v1);
    float d20 = Dot(v2, v0);
    float d21 = Dot(v2, v1);
    float denom = d00 * d11 - d01 * d01;
    v = (d11 * d20 - d01 * d21) / denom;
    w = (d00 * d21 - d01 * d20) / denom;
    u = 1.0f - v - w;
}

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

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

Джон Калсбек
источник
7
# операций может быть красная сельдь. От того, насколько они зависимы и насколько важны графики для современных процессоров. всегда проверяйте предположения и производительность "улучшения".
Шон Мидлдитч
1
Две рассматриваемые версии имеют почти одинаковую задержку на критическом пути, если вы смотрите только на скалярные математические операции. Что мне нравится в этом, так это то, что, заплатив место за два поплавка, вы можете сбрить одно вычитание и одно деление с критического пути. Это того стоит? Только тест производительности знает наверняка ...
Джон Калсбек
1
Он описывает, как он получил это на странице 137-138 с разделом «Ближайшая точка треугольника к точке»
бобобо
1
Небольшое примечание: нет аргумента pдля этой функции.
Барт
2
Незначительное замечание по реализации: если все 3 точки расположены друг над другом, вы получите ошибку «деление на 0», поэтому обязательно проверьте этот случай в фактическом коде.
frodo2975
9

Правило Крамера должно быть лучшим способом его решения. Я не графический парень, но мне было интересно, почему в книге «Обнаружение столкновений в реальном времени» они не делают следующую более простую вещь:

// Compute barycentric coordinates (u, v, w) for
// point p with respect to triangle (a, b, c)
void Barycentric(Point p, Point a, Point b, Point c, float &u, float &v, float &w)
{
    Vector v0 = b - a, v1 = c - a, v2 = p - a;
    den = v0.x * v1.y - v1.x * v0.y;
    v = (v2.x * v1.y - v1.x * v2.y) / den;
    w = (v0.x * v2.y - v2.x * v0.y) / den;
    u = 1.0f - v - w;
}

Это напрямую решает линейную систему 2x2

v v0 + w v1 = v2

в то время как метод из книги решает систему

(v v0 + w v1) dot v0 = v2 dot v0
(v v0 + w v1) dot v1 = v2 dot v1
user5302
источник
Разве ваше предлагаемое решение не делает предположений о третьем ( .z) измерении (в частности, о том, что оно не существует)?
Кукурузные початки
1
Это лучший метод, если вы работаете в 2D. Просто небольшое улучшение: нужно вычислить обратную величину знаменателя, чтобы использовать два умножения и одно деление вместо двух делений.
Рубик
8

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

// Compute barycentric coordinates (u, v, w) for
// point p with respect to triangle (a, b, c)
void Barycentric(Point a, Point b, Point c, float &u, float &v, float &w)
{
    Vector v0 = b - a, v1 = c - a, v2 = p - a;
    float d00 = Dot(v0, v0);
    float d01 = Dot(v0, v1);
    float d11 = Dot(v1, v1);
    float d20 = Dot(v2, v0);
    float d21 = Dot(v2, v1);
    float invDenom = 1.0 / (d00 * d11 - d01 * d01);
    v = (d11 * d20 - d01 * d21) * invDenom;
    w = (d00 * d21 - d01 * d20) * invDenom;
    u = 1.0f - v - w;
}

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

Vector v0;
Vector v1;
float d00;
float d01;
float d11;
float invDenom;

Итак, окончательный код выглядит так:

// Compute barycentric coordinates (u, v, w) for
// point p with respect to triangle (a, b, c)
void Barycentric(Point a, Point b, Point c, float &u, float &v, float &w)
{
    Vector v2 = p - a;
    float d20 = Dot(v2, v0);
    float d21 = Dot(v2, v1);
    v = (d11 * d20 - d01 * d21) * invDenom;
    w = (d00 * d21 - d01 * d20) * invDenom;
    u = 1.0f - v - w;
}
NielW
источник
2

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

В качестве альтернативы вы можете вычислить несколько барицентрических координат одновременно, используя sse или avx для ускорения в 4 или 8 раз.

Crowley9
источник
1

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

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

Герт
источник
1

Я пытался скопировать код @ NielW в C ++, но не получил правильных результатов.

Проще было прочитать https://en.wikipedia.org/wiki/Barycentric_coordinate_system#Barycentric_coordinates_on_triangles и вычислить lambda1 / 2/3, как указано там (векторные функции не нужны).

Если p (0..2) - это Точки треугольника с x / y / z:

Precalc для треугольника:

double invDET = 1./((p(1).y-p(2).y) * (p(0).x-p(2).x) + 
                   (p(2).x-p(1).x) * (p(0).y-p(2).y));

тогда лямбды для точки «точка»

double l1 = ((p(1).y-p(2).y) * (point.x-p(2).x) + (p(2).x-p(1).x) * (point.y-p(2).y)) * invDET; 
double l2 = ((p(2).y-p(0).y) * (point.x-p(2).x) + (p(0).x-p(2).x) * (point.y-p(2).y)) * invDET; 
double l3 = 1. - l1 - l2;
user1712200
источник
0

Для заданной точки N внутри треугольника ABC вы можете получить барицентрический вес точки C, разделив площадь подтреугольника ABN на общую площадь треугольника AB C.

ловкач
источник