Если на этапе растеризации в конвейере OpenGL происходит линейная интерполяция, и вершины уже были преобразованы в экранное пространство, где информация о глубине, используемая для корректно правильной интерполяции, исходит из?
Можно ли дать подробное описание того, как OpenGL переходит из примитивов экранного пространства в фрагменты с правильно интерполированными значениями?
Результатом вершинного шейдера является четырехкомпонентный вектор vec4 gl_Position. Из раздела 13.6 Преобразования координат ядра GL 4.4 spec:
Координаты
gl_Positionдля вершины получаются в результате выполнения шейдера, что дает координату вершиныgl_Position.Разделение перспективы по координатам клипа дает нормализованные координаты устройства, после чего следует преобразование области просмотра (см. Раздел 13.6.1) для преобразования этих координат в координаты окна.
OpenGL делит перспективу как
device.xyz = gl_Position.xyz / gl_Position.w
Но затем сохраняет 1/gl_Position.w как последний компонент gl_FragCoord:
gl_FragCoord.xyz = device.xyz scaled to viewport
gl_FragCoord.w = 1 / gl_Position.w
Это преобразование является биективным, поэтому информация о глубине не теряется. На самом деле, как мы видим ниже, 1/gl_Position.w имеет решающее значение для правильной интерполяции перспективы.
Для данного треугольника (P0, P1, P2) одним из способов параметризации точек внутри треугольника является выбор одной вершины (здесь P0) и выражение каждой другой точки в виде:
P(u,v) = P0 + (P1 - P0)u + (P2 - P0)v
где u> = 0, v> = 0 и u + v <= 1. Учитывая атрибут (f0, f1, f2) в вершинах треугольника, мы можем использовать u, v, чтобы интерполировать его по треугольнику
f(u,v) = f0 + (f1 - f0)u + (f2 - f0)v
Вся математика может быть выполнена с использованием описанной выше параметризации, и на самом деле иногда предпочтительнее из-за более быстрых вычислений. Однако это менее удобно и имеет числовые проблемы (например, P (1,0) может не совпадать с P1).
Вместо этого обычно используются барицентрические координаты. Каждая точка внутри треугольника является взвешенной суммой вершин:
P(b0,b1,b2) = P0*b0 + P1*b1 + P2*b2
f(b0,b1,b2) = f0*b0 + f1*b1 + f2*b2
где b0 + b1 + b2 = 1, b0> = 0, b1> = 0, b2> = 0 - барицентрические координаты точки в треугольнике. Каждое би можно представить как "сколько пи нужно смешать". Таким образом, b = (1,0,0), (0,1,0) и (0,0,1) - вершины треугольника, (1/3, 1/3, 1/3) - барицентр, и так далее.
Допустим, мы заполняем спроектированный 2D-треугольник на экране. Для каждого фрагмента у нас есть свои оконные координаты. Сначала мы вычисляем его барицентрические координаты, инвертируя функцию P(b0,b1,b2), которая является линейной функцией в оконных координатах. Это дает нам барицентрические координаты фрагмента на проекции 2D треугольника.
Правильная корректная интерполяция атрибута будет линейно изменяться в координатах клипа (и, соответственно, в мировых координатах). Для этого нам нужно получить барицентрические координаты фрагмента в пространстве клипа.
Как это происходит (см. [1] и [2]), глубина фрагмента не является линейной в координатах окна, но обратная глубина (1/gl_Position.w) равна. Соответственно, атрибуты и барицентрические координаты пространства клипа при взвешивании по обратной глубине изменяются линейно в координатах окна.
Поэтому мы вычисляем барицентрическую форму с поправкой на перспективу:
( b0 / gl_Position[0].w, b1 / gl_Position[1].w, b2 / gl_Position[2].w )
B = -------------------------------------------------------------------------
b0 / gl_Position[0].w + b1 / gl_Position[1].w + b2 / gl_Position[2].w
а затем использовать его для интерполяции атрибутов из вершин.
Примечание: GL_NV_fragment_shader_barycentric выставляет линейные барицентрические координаты устройства через gl_BaryCoordNoPerspNV а перспектива корректируется через gl_BaryCoordNV.
Вот код C++, который растеризует и затеняет треугольник на процессоре, подобно OpenGL. Я рекомендую вам сравнить его с шейдерами, перечисленными ниже:
struct Renderbuffer {
int w, h, ys;
void *data;
};
struct Vert {
vec4f position;
vec4f texcoord;
vec4f color;
};
struct Varying {
vec4f texcoord;
vec4f color;
};
void vertex_shader(const Vert &in, vec4f &gl_Position, Varying &out)
{
out.texcoord = in.texcoord;
out.color = in.color;
gl_Position = { in.position[0], in.position[1], -2*in.position[2] - 2*in.position[3], -in.position[2] };
}
void fragment_shader(vec4f &gl_FragCoord, const Varying &in, vec4f &out)
{
out = in.color;
vec2f wrapped = vec2f(in.texcoord - floor(in.texcoord));
bool brighter = (wrapped[0] < 0.5) != (wrapped[1] < 0.5);
if(!brighter)
(vec3f&)out = 0.5f*(vec3f&)out;
}
void store_color(Renderbuffer &buf, int x, int y, const vec4f &c)
{
// can do alpha composition here
uint8_t *p = (uint8_t*)buf.data + buf.ys*(buf.h - y - 1) + 4*x;
p[0] = linear_to_srgb8(c[0]);
p[1] = linear_to_srgb8(c[1]);
p[2] = linear_to_srgb8(c[2]);
p[3] = lrint(c[3]*255);
}
void draw_triangle(Renderbuffer &color_attachment, const box2f &viewport, const Vert *verts)
{
Varying perVertex[3];
vec4f gl_Position[3];
box2f aabbf = { viewport.hi, viewport.lo };
for(int i = 0; i < 3; ++i)
{
// invoke the vertex shader
vertex_shader(verts[i], gl_Position[i], perVertex[i]);
// convert to device coordinates by perspective division
gl_Position[i][3] = 1/gl_Position[i][3];
gl_Position[i][0] *= gl_Position[i][3];
gl_Position[i][1] *= gl_Position[i][3];
gl_Position[i][2] *= gl_Position[i][3];
// convert to window coordinates
auto &pos2 = (vec2f&)gl_Position[i];
pos2 = mix(viewport.lo, viewport.hi, 0.5f*(pos2 + vec2f(1)));
aabbf = join(aabbf, (const vec2f&)gl_Position[i]);
}
// precompute the affine transform from fragment coordinates to barycentric coordinates
const float denom = 1/((gl_Position[0][0] - gl_Position[2][0])*(gl_Position[1][1] - gl_Position[0][1]) - (gl_Position[0][0] - gl_Position[1][0])*(gl_Position[2][1] - gl_Position[0][1]));
const vec3f barycentric_d0 = denom*vec3f( gl_Position[1][1] - gl_Position[2][1], gl_Position[2][1] - gl_Position[0][1], gl_Position[0][1] - gl_Position[1][1] );
const vec3f barycentric_d1 = denom*vec3f( gl_Position[2][0] - gl_Position[1][0], gl_Position[0][0] - gl_Position[2][0], gl_Position[1][0] - gl_Position[0][0] );
const vec3f barycentric_0 = denom*vec3f(
gl_Position[1][0]*gl_Position[2][1] - gl_Position[2][0]*gl_Position[1][1],
gl_Position[2][0]*gl_Position[0][1] - gl_Position[0][0]*gl_Position[2][1],
gl_Position[0][0]*gl_Position[1][1] - gl_Position[1][0]*gl_Position[0][1]
);
// loop over all pixels in the rectangle bounding the triangle
const box2i aabb = lrint(aabbf);
for(int y = aabb.lo[1]; y < aabb.hi[1]; ++y)
for(int x = aabb.lo[0]; x < aabb.hi[0]; ++x)
{
vec4f gl_FragCoord;
gl_FragCoord[0] = x + 0.5;
gl_FragCoord[1] = y + 0.5;
// fragment barycentric coordinates in window coordinates
const vec3f barycentric = gl_FragCoord[0]*barycentric_d0 + gl_FragCoord[1]*barycentric_d1 + barycentric_0;
// discard fragment outside the triangle. this doesn't handle edges correctly.
if(barycentric[0] < 0 || barycentric[1] < 0 || barycentric[2] < 0)
continue;
// interpolate inverse depth linearly
gl_FragCoord[2] = dot(barycentric, vec3f(gl_Position[0][2], gl_Position[1][2], gl_Position[2][2]));
gl_FragCoord[3] = dot(barycentric, vec3f(gl_Position[0][3], gl_Position[1][3], gl_Position[2][3]));
// clip fragments to the near/far planes (as if by GL_ZERO_TO_ONE)
if(gl_FragCoord[2] < 0 || gl_FragCoord[2] > 1)
continue;
// convert to perspective correct (clip-space) barycentric
const vec3f perspective = 1/gl_FragCoord[3]*barycentric*vec3f(gl_Position[0][3], gl_Position[1][3], gl_Position[2][3]);
// interpolate the attributes using the perspective correct barycentric
Varying varying;
for(int i = 0; i < sizeof(Varying)/sizeof(float); ++i)
((float*)&varying)[i] = dot(perspective, vec3f(
((const float*)&perVertex[0])[i],
((const float*)&perVertex[1])[i],
((const float*)&perVertex[2])[i]
));
// invoke the fragment shader and store the result
vec4f color;
fragment_shader(gl_FragCoord, varying, color);
store_color(color_attachment, x, y, color);
}
}
int main()
{
Renderbuffer buffer = { 512, 512, 512*4 };
buffer.data = malloc(buffer.ys * buffer.h);
memset(buffer.data, 0, buffer.ys * buffer.h);
// interleaved attributes buffer
Vert verts[] = {
{ { -1, -1, -2, 1 }, { 0, 0, 0, 1 }, { 0, 0, 1, 1 } },
{ { 1, -1, -1, 1 }, { 10, 0, 0, 1 }, { 1, 0, 0, 1 } },
{ { 0, 1, -1, 1 }, { 0, 10, 0, 1 }, { 0, 1, 0, 1 } },
};
box2f viewport = { 0, 0, buffer.w, buffer.h };
draw_triangle(buffer, viewport, verts);
lodepng_encode32_file("out.png", (unsigned char*)buffer.data, buffer.w, buffer.h);
}
Вот шейдеры OpenGL, используемые для генерации эталонного изображения.
Вершинный шейдер:
#version 450 core
layout(location = 0) in vec4 position;
layout(location = 1) in vec4 texcoord;
layout(location = 2) in vec4 color;
out gl_PerVertex {
vec4 gl_Position;
};
layout(location = 0) out PerVertex {
vec4 texcoord;
vec4 color;
} OUT;
void main() {
OUT.texcoord = texcoord;
OUT.color = color;
gl_Position = vec4(position[0], position[1], -2*position[2] - 2*position[3], -position[2]);
}
Фрагмент шейдера:
#version 450 core
layout(location = 0) in PerVertex {
vec4 texcoord;
vec4 color;
} IN;
layout(location = 0) out vec4 OUT;
void main() {
OUT = IN.color;
vec2 wrapped = fract(IN.texcoord.xy);
bool brighter = (wrapped[0] < 0.5) != (wrapped[1] < 0.5);
if(!brighter)
OUT.rgb *= 0.5;
}
Вот почти идентичные изображения, сгенерированные кодом C++ (слева) и OpenGL (справа):
Различия вызваны различной точностью и режимами округления.
Для сравнения приведем пример неправильной перспективы (для интерполяции в вышеприведенном коде используется barycentric вместо perspective):
Формула, которую вы найдете в спецификации GL (см. Стр. 427; ссылка - текущая спецификация 4.4, но так было всегда) для интерполяции с коррекцией перспективы значения атрибута в треугольнике:
a * f_a / w_a + b * f_b / w_b + c * f_c / w_c
f=-----------------------------------------------------
a / w_a + b / w_b + c / w_c
где a,b,c обозначают барицентрические координаты точки в треугольнике, для которого мы интерполируем (a,b,c >=0, a+b+c = 1), f_i значение атрибута в вершине i и w_i пространство клипа w координатой вершины i. Обратите внимание, что барицентрические координаты рассчитываются только для 2D-проекции координат оконного пространства треугольника (поэтому z игнорируется).
Это то, что формулы, которые Ибунгалоубилл дал в своем прекрасном ответе, сводятся, в общем случае, к произвольной оси проекции. Фактически, последняя строка матрицы проекции определяет только ось проекции, к которой плоскость изображения будет ортогональна, а компонент пространства клипа w является просто точечным произведением между координатами вершины и этой осью.
В типичном случае матрица проекции имеет (0,0, -1, 0) в качестве последней строки, поэтому она w_clip = -z_eye так, что w_clip = -z_eye, и это то, что использовала ybungalowbill. Однако, поскольку w - это то, чем мы на самом деле будем заниматься делением (это единственный нелинейный шаг во всей цепочке преобразования), это будет работать для любой оси проекции. Это также будет работать в тривиальном случае ортогональных проекций, где w всегда равно 1 (или, по крайней мере, постоянно).
Обратите внимание на несколько вещей для эффективной реализации этого. Инверсия 1/w_i может быть предварительно рассчитана для каждой вершины (позвольте в дальнейшем назвать их q_i), ее не нужно переоценивать для каждого фрагмента. И это абсолютно бесплатно, так как мы w любом случае делим на w при переходе в пространство NDC, поэтому мы можем сохранить это значение. Спецификация GL никогда не описывает, как определенная функция должна быть реализована внутри, но тот факт, что координаты экранного пространства будут доступны в glFragCoord.xyz, а gl_FragCoord.w гарантированно предоставит (линеаризованное интерполированное) 1/w пространство клипа координата довольно показательна здесь. Это значение 1_w каждого фрагмента фактически является знаменателем приведенной выше формулы.
Коэффициенты a/w_a, b/w_b и c/w_c каждый используются в формуле два раза. И они также постоянны для любого значения атрибута, теперь независимо от того, сколько атрибутов нужно интерполировать. Таким образом, для каждого фрагмента вы можете вычислить a'=q_a * a, b'=q_b * b и c'=q_c и получить
a' * f_a + b' * f_b + c' * f_c
f=------------------------------
a' + b' + c'
Таким образом, перспективная интерполяция сводится к
за фрагмент.