В большинстве систем 3D-рендеринга, включая OpenGL, базовыми примитивами являются выпуклые многоугольники. Напомню, что многоугольник называется выпуклым, если любые две точки, принадлежащие данному многоугольнику, могут быть соединены отрезком, все точки которого тоже принадлежат данному многоугольнику. Примерами выпуклых многоугольников являются любой треугольник, а также любой правильный многоугольник. Своим широким применением в трехмерной графике выпуклые многоугольники (ВМ) обязаны характерным свойствам, некоторые из которых перечислены ниже:
Для низкоуровневого закрашивания и наложения текстур важно еще одно свойство: Если сторона ВМ не лежит на данной прямой, эта прямая
Многие свойства ВМ справедливы и для выпуклых многогранников. В частности, при сечении выпуклого многогранника плоскостью, образуются два выпуклых многогранника. Также если точка находится внутри выпуклого многогранника, полупрямая, исходящая из этой точки, пересечет только одну грань многогранника.
Большинство алгоритмов, рассмотренных далее, работает корректно только с ВМ, по этому, прежде всего, следует проверить,
Для того, чтобы решить эту задачу, будем сравнивать последовательно ориентацию треугольников, образованных вершинами смежных сторон. В ВМ ориентация всех таких треугольников одинакова и совпадает с ориентацией самого ВМ. У невыпуклых многоугольников ориентация треугольников различается. Данное замечание вообще говоря справедливо только для простых многоугольников, т. е. таких, стороны которых не пересекаются. Если в вашем приложении возможно появление пересекающихся многоугольников, следует сначала убедиться, что многоугольник - простой.
Рисунок 1 Иллюстрирует ситуацию с выпуклым и невыпуклым многоугольниками при направлении обхода вершин против часовой стрелки.
Прежде, чем представить алгоритм проверки на языке C++, введем вспомогательные типы данных:
typedef struct {
float x, y;
} Vect2D;
typedef struct {
float x, y, z;
} Vect3D;
и вспомогательную функцию, классифицирующую положение точки относительно прямой, заданной двумя точками:
inline int classify_point2D(Vect2D p0, Vect2D p1, Vect2D p2) {
return (p1.x-p0.x)*(p2.y-p0.y) - (p2.x-p0.x)*(p1.y-p0.y);
}
Функция classify_point2D возвращает положительное значение[2], если точка p2 лежит слева от прямой, проходящей через p0 и p1, (если смотреть в направлении, заданном вектором p1-p0), отрицательное значение, если точка p2 лежит справа от прямой и 0, если p2 принадлежит прямой. Эту же функцию можно использовать для проверки ориентации треугольника p0,p1, p2, как показано на рисунке 2.
Рисунок 2. Определение ориентации треугольника P0,P1,P2 с помощью classify_point2D. Если P2 лежит слева от отрезка P0,P1, треугольник ориентирован против часовой стрелки, если справа - то по часовой.
Как видим, classify_point2D - очень полезная функция, и мы воспользуемся ей неоднократно.
Ну а теперь сам алгоритм. Функция is_poly_convex возвращает true, если многоугольник выпуклый и false в противном случае.
bool is_poly_convex2D(Vect2D * poly, int np) {
// треугольник - всегда выпуклый.
if (np == 3) return true;
int res = SIGN(classify_point2D(poly[0], poly[1], poly[2]));
for (int i = 1; i < np; i++)
if (res != SIGN(classify_point2D(poly[i], poly[(i+1)%np], poly[(i+2) %np]))
return false;
return true;
}
Следующий алгоритм, определяет ориентацию простого многоугольника (выпуклого или невыпуклого) "в целом". Часто эту задачу решают, вычисляя площадь многоугольника со знаком. Мы же рассмотрим другой алгоритм, являющийся, на мой взгляд, самым быстрым.
int poly_orientation2D(Vect2D * poly, int np) {
int vmin = 0;
float xmin = poly[0].x;
float ymin = poly[0].y;
// Находим самую правую нижнюю вершину
for (int i = 1; i < np; i++) {
if (poly[i].y > ymin) continue;
if ((poly[i].y == ymin) && (poly[i].x < xmin)) continue;
vmin = i;
xmin = poly[i].x;
ymin = poly[i].y;
}
// Проверяем ориентацию треугольника
if (vmin == 0)
return classify_point2D(poly[np-1], poly[0], poly[1]);
else
return classify_point2D(poly[vmin-1], poly[vmin], poly[(vmin+1)%np]);
}
Функция возвращает положительное значение, если многоугольник ориентирован против часовой стрелки и отрицательное значение в противном случае. Объяснить этот алгоритм можно так: если многоугольник выпуклый, его ориентация совпадает с ориентацией любого из треугольников, образованных соседними вершинами. Если же треугольник невыпуклый, самая правая нижняя точка не может быть точкой "впадины" и ориентация треугольника (vmin-1), vmin, (vmin+1) совпадает с ориентацией всего многоугольника.
А как же быть в случае плоского многоугольника в трехмерном пространстве? Возникает соблазн использовать векторные операции - вычислять векторное произведение для вершин каждого треугольника и сравнивать затем направление векторов. Однако и в этом случае можно воспользоваться приведенной выше функцией, которая выполняет меньше вычислений, чем векторный метод. Обычно, после выполнения всех предварительных преобразований плоскость наблюдателя совпадает с одной из координатных плоскостей (обычно x, y). Определение ориентации многоугольника как правило имеет смысл только если плоскость многоугольника не перпендикулярна плоскости наблюдателя. Заметим, что в этом случае параллельная проекция многоугольника на плоскость наблюдателя будет ориентирована относительно наблюдателя также, как и сам многоугольник. Таким образом, нам следует спроецировать многоугольник на эту плоскость, просто отбросив "лишнюю" координату. Получившийся двумерный многоугольник будет иметь туже ориентацию, что и исходный трехмерный. Строго говоря, ориентации многоугольников совпадают, когда угол между плоскостью x,y и плоскостью наблюдателя лежит в интервале (-90°..90°).
Отсечение широко применяется в компьютерной графике. В этом разделе мы рассмотрим двумерный случай отсечения многоугольника прямой и трехмерный случай отсечения многоугольника плоскостью. Прежде всего нам потребуется еще одна вспомогательная функция:
bool intersection_point2D(Vect2D * dst, Vect2D * seg1, Vect2D * seg2) {
if ( (seg1[0].x == seg1[1].x) && (seg2[0].x == seg2[1].x) ) return false;
if ( (seg1[0].y == seg1[1].y) && (seg2[0].y == seg2[1].y) ) return false;
float d = (seg1[1].x-seg1[0].x)*(seg2[1].y-seg2[0].y) - (seg1[1].y-seg1[0].y)*(seg2[1].x-seg2[0].x);
if (d == 0) return false;
float r = ((seg1[0].y-seg2[0].y)*(seg2[1].x-seg2[0].x) - (seg1[0].x-seg2[0].x)*(seg2[1].y-seg2[0].y))/d;
dst->x = (seg1[1].x-seg1[0].x)*r + seg1[0].x;
dst->y = (seg1[1].y-seg1[0].y)*r + seg1[0].y;
return true;
}
Функция intersection_point2D находит точку пересечения прямых, заданных наборами точек seg1 и seg2. Функция возвращает false если прямые не пересекаются. Теперь объявим структуру, в которой будут передаваться аргументы для функции выполняющей двумерное отсечение ВМ.
typedef struct {
// Вершины результирующего многоугольника
Vect2D * dst;
// их число
int nd;
// вершины исходного многоугольника
Vect2D * poly;
// их число
int np;
// первая точка отсекающей прямой
Vect2D p1;
// вторая точка отсекающей прямой
Vect2D p2;
// Определяет, в какой полуплоскости лежит результирующий многоугольник
int cl;
} Clip2D;
Последний член структуры требует пояснения. При отсечении многоугольника прямой получаются два новых многоугольника. Параметр cl определяет, из какой полуплоскости (в терминах функции classify_point2D) следует брать результирующий многоугольник. Далее следует сама функция.
void clip2D(Clip2D * cp) {
Vect2D clippoint[2];
Vect2D seg1[2];
Vect2D seg2[2];
int pp[2];
int pcount = 0;
seg1[0] = cp->p1;
seg1[1] = cp->p2;
for(int i = 0; i < cp->np; i++)
if (classify_point2D(cp->p1, cp->p2, cp->poly[i]) !=
classify_point2D(cp->p1, cp->p2, cp->poly[(i+1)%cp->np])) {
pp[pcount] = (i + 1)%cp->np;
seg2[0] = cp->poly[i];
seg2[1] = cp->poly[(i+1)%cp->np];
intersection_point2D(&clippoint[pcount], seg1, seg2);
if ( ++pcount > 1 ) break;
}
if (pcount < 2) {
if ( SIGN(classify_point2D(cp->p1, cp->p2, poly[0])) == SIGN(cp->cl)) {
cp->nd = cp->np;
cp->dst = (Vect2D*) malloc(cp->nd * sizeof(Vect2D));
memcpy(cp->dst, cp->poly, cp->nd * sizeof(Vect2D));
} else {
cp->nd = 0;
cp->dst = NULL;
}
} else {
if (SIGN(classify_point2D(cp->p1, cp->p2, poly[0])) == SIGN(cp->cl)) {
cp->nd = (pp[1] == 0 ? cp->np - pp[0] : pp[1] - pp[0]) + 2;
cp->dst = (Vect2D*) malloc(cp->nd * sizeof(Vect2D));
for( int i = 0; i < (cp->nd - 2); i++) cp->dst[i] =
cp->poly[(i+pp[0])%cp->np];
cp->dst[cp->nd-2] = clippoint[1];
cp->dst[cp->nd-1] = clippoint[0];
} else {
cp->nd = (pp[1] == 0 ? pp[0] : cp->np - pp[1] + pp[0]) + 2;
cp->dst = (Vect2D*) malloc(cp->nd * sizeof(Vect2D));
for( int i = 0; i < (cp->nd - 2); i++) cp->dst[i] =
cp->poly[(i+pp[1])%cp->np];
cp->dst[cp->nd-2] = clippoint[0];
cp->dst[cp->nd-1] = clippoint[1];
}
}
}
Результат функции возвращается в параметрах cp->dst и cp->nd. Если весь исходный многоугольник лежит в полуплоскости, заданной cp->cl, возвращается исходный многоугольник. Если весь исходный многоугольник лежит в другой полуплоскости, cp->dst возвращает NULL, а cp->nd - 0.
В остальных случаях возвращается многоугольник, полученный в результате отсечения. Блок памяти, на который указывает cp->dst, должен освобождаться функцией free. С помощью функции отсечения можно реализовать простой алгоритм нахождения ВМ, получающегося в результате пересечения двух ВМ.
Для начала определим структуру для передачи параметров:
typedef struct {
// Вершины результирующего многоугольника
Vect2D * dst;
// их число
int nd;
// вершины первого многоугольника
Vect2D * poly1;
// их число
int np1;
// вершины второго многоугольника
Vect2D * poly2;
// их число
int np2;
} IntersectPoly;
Ниже следует функция, реализующая алгоритм:
void intersect_poly2D(IntersectPoly * ip) {
Clip2D cp;
cp.poly = (Vect2D*) malloc(ip->np1 * sizeof(Vect2D));
cp.np = ip->np1;
for (int i = 0; i < ip->np1; i++) cp.poly[i] = ip->poly1[i];
for (int i = 0; i < ip->np2; i++) {
cp.p1 = ip->poly2[i];
cp.p2 = ip->poly2[(i+1)%ip->np2];
cp.cl = classify_point2D(ip->poly2[i], ip->poly2[(i+1)%ip->np2],
ip->poly2[(i+2)%ip->np2]);
clip2D(&cp);
free(cp.poly);
if (cp.dst == NULL) break;
cp.poly = cp.dst;
cp.np = cp.nd;
}
ip->dst = cp.dst;
ip->nd = cp.nd;
}
Если многоугольники не пересекаются, cp->dst возвращает NULL, а cp->nd - 0. Блок памяти, на который указывает cp->dst, должен освобождаться функцией free.
Результирующий многоугольник получается в результате последовательного отсечения первого многоугольника прямыми, проходящими через стороны второго. Рисунок 3 иллюстрирует принцип работы алгоритма.
Данный алгоритм, безусловно, не является самым быстрым (время его работы пропорционально n2, где n - среднее число сторон многоугольников), но он достаточно компактен и легко справляется со всеми частными случаями (один многоугольник лежит внутри другого и т. п.).
Для решения трехмерной задачи об отсечении многоугольника, принадлежащего некоторой плоскости, другой плоскостью, нам потребуется вспомогательная процедура, вычисляющая точку пересечения плоскости и прямой, проходящей через две заданные точки. При решении этой задачи мы будем оперировать векторами, так что полезно ввести следующие операторы:
// Скалярное произведение векторов
float operator * (Vect3D a, Vect3D b) {
return a.x*b.x + a.y*b.y + a.z*b.z;
}
// Умножение вектора на число
Vect3D operator * (float a, Vect3D b) {
Vect3D res = {a*b.x, a*b.y, a*b.z};
return res;
}
Аналогично определим операторы для сложения и вычитания векторов. Точка пересечения прямой и плоскости вычисляется в следующей функции:
bool plane_line_intersection3D (Vect3D * dst, Vect3D v, Vect3D n, Vect3D p0, Vect3D p1) {
float s = n*(p1-p0);
if (s == 0) return false;
*dst = p0 + ((n*(v-p0))/s)*(p1-p0);
return true;
}
Функция возвращает false, если прямая параллельна плоскости. Если функция возвращает true, параметр dst содержит координаты точки пересечения. Параметр v - какая-либо точка на плоскости, n - вектор, перпендикулярный плоскости, p0, p1 - точки, через которые проходит прямая.
Ну и, наконец, функция, выполняющая отсечение в трехмерном пространстве.
typedef struct {
// Вершины результирующего многоугольника
Vect3D * dst;
// их число
int nd;
// нормаль к плоскости
Vect3D normal;
// точка на плоскости
Vect3D V;
// Вершины отсекаемого многоугольника
Vect3D * poly;
// их число
int np;
float cl;
} Clip3D;
Многоугольник poly отсекается плоскостью, заданной параметрами normal и v. Смысл параметра cl - такой же, как и в случае двумерного отсечения.
void clip3D(Clip3D * cp) {
int pp[2];
int pcount = 0;
Vect3D clippoint[2];
float d = cp->normal*cp->v;
for (int i = 0; i < cp->np; i++) {
if ((SIGN(cp->normal*cp->poly[i]-d)) != (SIGN(cp->normal*cp->poly[(i+1)%cp->np]-d))) {
pp[pcount] = (i+1)%cp->np;
plane_line_intersection3D(&clippoint[pcount++], cp->v, cp->n, cp->poly[i], cp->poly[(i+1)%cp->np]);
}
if (pcount > 1) break;
}
if ( pcount < 2 ) {
if (SIGN(cp->normal*cp->poly[0]-d) == SIGN(cp->cl)) {
cp->nd = cp->np;
cp->dst = (Vect3D*) malloc(cp->nd * sizeof(Vect3D));
memcpy(cp->dst, cp->poly, cp->nd * sizeof(Vect3D));
} else {
cp->nd = 0;
cp->dst = NULL;
}
} else {
if ((SIGN(cp->normal*clippoint[0]-d)) == (SIGN(cp->cl))) {
cp->nd = (pp[1] == 0 ? cp->np - pp[0] : pp[1] - pp[0]) + 2;
cp->dst = (Vect3D*) malloc(cp->nd * sizeof(Vect3D));
for (int i = 0; i < (cp->nd - 2); i++)
cp->dst[i] = cp->poly[(i+pp[0])%cp->np];
cp->dst[cp->nd-2] = clippoint[1];
cp->dst[cp->nd-1] = clippoint[0];
} else {
cp->nd = (pp[1] == 0 ? pp[0] : cp->np - pp[1] + pp[0]) + 2;
cp->dst = (Vect3D*) malloc(cp->nd * sizeof(Vect3D));
for (int i = 0; i < (cp->nd - 2); i++)
cp->dst[i] = cp->poly[(i+pp[1])%cp->np];
cp->dst[cp->nd-2] = clippoint[0];
cp->dst[cp->nd-1] = clippoint[1];
}
}
}
Блок памяти, на который указывает cp->dst, должен освобождаться функцией free.
Многие двумерные алгоритмы легко применить к многоугольникам в пространстве, просто спроецировав многоугольник на одну из координатных плоскостей путем отбрасывания третьей координаты. Если в процессе выполнения алгоритма появляются новые точки, мы можем "вернуть" их в 3D-базис, просто вычислив третью координату на исходной плоскости.
Однако если плоскость многоугольника не параллельна координатной плоскости, при таком проецировании возникают линейные искажения. Если такие искажения недопустимы, необходимо выполнить переход к 2D-базису. Одна из ситуаций, где такой переход может понадобиться - генерация текстурных координат. Самым очевидным способом перехода к 2D базису является поворот плоскости многоугольника до совпадения ее с координатной плоскостью. Однако такое решение требует довольно большого объема вычислений (определение угла между плоскостями, тригонометрические вычисления, перемножение матриц). Существует способ перехода к 2D-базису без вращения плоскостей. Этот способ заключается в задании на плоскости двух перпендикулярных векторов (ортов).
Рассмотрим, как это делается. Сначала нужно перенести плоскость таким образом, чтобы она проходила через начало координат. В нашем случае достаточно просто вычесть из координат всех вершин многоугольника координаты одной из вершин. Затем мы строим орты. Один из ортов можно вычислить из векторного уравнения
, где
- искомый базисный вектор (орт), а
- нормаль к плоскости.
Второй орт получается в результате векторного произведения
.
Для того чтобы перейти в новый базис, нужно скалярно перемножить каждую вершину многоугольника с каждым ортом. Полученные пары чисел и будут соответствовать координатам многоугольника в новом базисе. Если необходимо соблюдать масштаб, орты следует предварительно нормализовать. Рисунок 4 иллюстрирует данный подход.
Рисунок 4. Плоскость, проходящая через начало координат, нормаль к плоскости, новые орты на плоскости и точка плоскости в новом базисе.
В качестве примера приведу фрагмент функции, генерирующей текстурные координаты для многоугольника.
...
int np; // число вершин многоугольника
Vect3D * poly; // массив вершин многоугольника
Vect3D normal; // нормаль многоугольника
...
Vect3D Vmin;
Vect3D orts[2];
Vect2D * 2dcoords;
// Вычисляем орты
if (normal.x > 0) {
orts[0].x = -normal.z;
orts[0].y = 0;
orts[0].z = normal.x;
} else if (normal.y > 0) {
orts[0].x = 0;
orts[0].y = -normal.z;
orts[0].z = normal.y;
} else {
orts[0].x = 1;
orts[0].y = 0;
orts[0].z = 0;
}
orts[1] = orts[0]^normal;
normalize(&orts[0]);
normalize(&orts[1]);
// Находим наименьшие x и y и переносим плоскость в начало координат
Vmin.x = poly[0].x;
Vmin.y = poly[0].y;
for (int i = 1; i < np; i++) {
if (Vmin.x > poly[i].x) {
Vmin.x = poly[i].x;
if(normal.x == 0) Vmin.z = poly[i].z;
}
if (Vmin.y > poly[i].y) {
Vmin.y = poly[i].y;
if(normal.y == 0) Vmin.z = poly[i].z;
}
if (normal.z != 0)
Vmin.z = (normal*poly[0] - Vmin.x*normal.x - Vmin.y*normal.y)/normal.z;
for (int i = 0; i < np; i++) poly[i] = poly[i]-Vmin;
// Генерируем двумерные координаты
2dcoords = (Vect2D*) malloc(np*sizeof(Vect2D));
for (int i = 0; i < np; i++) {
2dcoords[i].x = poly[i]*orts[0];
2dcoords[i].y = poly[i]*orts[1];
}
...
Здесь символ ^ выполняет роль оператора векторного произведения, который определяется следующим образом:
Vect3D operator ^ (Vect3D a, Vect3D b) {
Vect3D res = {a.y * b.z - a.z * b.y, a.z * b.x - a.x * b.z, a.x * b.y - a.y * b.x;};
return res;
}
Функция normalize нормализует векторы (т. е. делит координаты вектора на значение его длины).
Отсечение по сектору обзора (view frustum clipping) является важной частью процесса нахождения видимых поверхностей в сложных трехмерных сценах. Обычно сектор обзора (view frustum) представляет собой усеченную прямую пирамиду с четырехугольным основанием.
Рассмотрим рис. 5. На нем сектор обзора представлен в параллельной проекции "вид сверху".
Тонкие синие линии - продолжения плоскостей граней сектора. Толстые синие отрезки - проекции многоугольников (для упрощения картинки все многоугольники перпендикулярны по отношению к нашей проекции). Многоугольник D полностью лежит в секторе обзора, Многоугольник B попадает в сектор обзора частично, а многоугольник A не попадает в сектор обзора.
В своем графическом движке вы, скорее всего, не будете выполнять отсечения многоугольников, частично попадающих в сектор обзора, так как это потребует большого объема вычислений на центральном процессоре, а предоставите эту операцию нижележащему уровню (например, OpenGL с аппаратным ускорением). Благодаря средствам анализа глубины (depth buffer) и отсечения, OpenGL может корректно сортировать и отображать все многоугольники, составляющие сцену. Тем не менее, возлагать на уровень OpenGL обработку всех многоугольников нецелесообразно, так как это опять же приведет к резкому снижению производительности приложения.
В связи с вышеизложенным возникает довольно тонкая проблема, касающаяся "разделения обязанностей" между уровнем приложения и уровнем OpenGL. С одной стороны, желательно уменьшить число многоугольников, передаваемых на нижний уровень. С другой стороны, алгоритм определения видимых многоугольников не должен быть слишком ресурсоемким, а значит, его нельзя сделать абсолютно точным. В современных приложениях 3D-графики применяются консервативные алгоритмы определения видимых поверхностей. Такие алгоритмы достаточно быстры, но их недостаток заключается в том, что они не все многоугольники, определенные консервативными алгоритмами как видимые, являются таковыми на самом деле.
Рассмотрим классический консервативный алгоритм отсечения по сектору обзора. Каждая грань сектора лежит на плоскости, нормаль которой обращена вовнутрь сектора. Классический алгоритм считает невидимым всякий многоугольник, полностью лежащий по другую сторону хотя бы одной из плоскостей относительно ее нормали. Алгоритм не ошибается в том смысле, что всякий такой многоугольник действительно будет невидим. Примером такого многоугольника на нашем рисунке является многоугольник A.
Однако классический алгоритм сочтет видимыми не только все полностью или частично видимые, но некоторые невидимые многоугольники (на рисунке это многоугольники C и E). В этом и заключается "консервативность" данного алгоритма. Конечно, для OpenGL эти лишние многоугольники не будут проблемой, однако нам хочется уменьшить количество лишних многоугольников, тем более что в некоторых сценах их число может быть весьма велико.
Для того чтобы усовершенствовать классический алгоритм, взглянем внимательно на рис. 5. Плоскости, в которых лежат многоугольники C и E не пересекают граней сектора. Очевидно, что если плоскость не пересекает граней сектора, лежащий на ней многоугольник не может быть видимым. Таким образом, алгоритм можно усовершенствовать, добавив к проверке расположения многоугольника относительно плоскостей граней сектора, проверку пересечения граней сектора и плоскостей многоугольника.
Насколько это усложнит алгоритм? У сектора 6 граней. Любая плоскость, пересекающая сектор, пересекает по меньшей мере 3 его грани. Значит, достаточно выполнять проверку только для четырех граней сектора . Далее, проверка выполняется только для тех многоугольников, которые признаны видимыми. Для большинства из этих многоугольников не нужно будет перебирать все 4 грани. Таким образом, усложнение алгоритма отсечения может быть оправдано. Разумеется, даже при этой дополнительной проверке алгоритм останется консервативным, так как некоторые невидимые многоугольники смогут пройти и ее.
Далее следует функция, демонстрирующая алгоритм:
// Сторона сектора, 5-й вектор - нормаль
typedef Vect3D FSide[5];
// Сектор
typedef FSide Frustum[6];
bool poly_in_frustum3D(Frustum frustum, Vect3D * poly, int np) {
// Сначала классический алгоритм
for (int i = 0; i < 6; i++)
for (int j = 0; j < np; j++) {
if ((poly[j]-frustum[i][0])*frustum[i][4] > 0) break;
if (j == np-1) return false;
}
// Проверяем пересечение плоскости и граней
float normal = (poly[1]-poly[0])^(poly[2]-poly[1]);
for(int i = 0; i < 4; i++) {
int res = SIGN((frustum[i][0]-poly[0])*normal);
for (int j = 1; j < 4; j++)
if (res != SIGN((frustum[i][j]-poly[0])*normal)
return true;
}
return false;
}
(c) 2002 Андрей Боровский, контакты: anb@symmetrica.net