Ответ 1

Результатом вершинного шейдера является четырехкомпонентный вектор 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 шейдеры

Вот шейдеры 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):

Ответ 2

Формула, которую вы найдете в спецификации 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. Обратите внимание на несколько вещей для эффективной реализации этого. Инверсия 1/w_i может быть предварительно рассчитана для каждой вершины (позвольте в дальнейшем назвать их q_i), ее не нужно переоценивать для каждого фрагмента. И это абсолютно бесплатно, так как мы w любом случае делим на w при переходе в пространство NDC, поэтому мы можем сохранить это значение. Спецификация GL никогда не описывает, как определенная функция должна быть реализована внутри, но тот факт, что координаты экранного пространства будут доступны в glFragCoord.xyz, а gl_FragCoord.w гарантированно предоставит (линеаризованное интерполированное) 1/w пространство клипа координата довольно показательна здесь. Это значение 1_w каждого фрагмента фактически является знаменателем приведенной выше формулы.

2. Коэффициенты 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'
   

Таким образом, перспективная интерполяция сводится к

• 3 дополнительных умножения,
• 2 дополнительных дополнения и
• 1 дополнительное подразделение

за фрагмент.