Рецепт Trapezoidal Shadow Maps (TSM)

Статья описывающая метод TSM была опубликована на Симпозиуме по Рендерингу Eurographics 2004 (Eurographics Symposium on Rendering 2004) и может быть найдена здесь. По той же ссылке вы также можете найти видеоролики и презентацию EGSR2004, небольшое сравнение нашей техники с Perspective Shadow Maps, Bounding Box Approximation и Standard Shadow Maps. Эта страница предоставит конкретные детали, чтобы упростить реализацию нашего подхода. Она рассматривает вычисление N_T и сохранение polygon offset с помощью fragment programs. В этом документе мы предполагаем, что вы уже прочитали статью по TSM и можете вычислить четыре вершины трапеции.

Идея TSM

Мы лучше используем разрешение shadow map аппроксимируя усеченную пирамиду камеры (eye's frustum) видимую источником света до трапеции, а затем искривляем ее и накладываем на shadow map. Это увеличивает кол-во сэмплов для областей более близких к камере и таким образом приводит к более высокому качеству теней. Трапеция вычисляется таким образом, что достигается плавное изменение разрешения тени. Трапеция вычисляется основываясь только на восьми вершинах усеченной пирамиды камеры, а не на всей сцене (как например в PSM), что устраняет проблемы непрерывности проявляющиеся в динамических сценах. Более того, аппрокисмация трапеции это постоянная процедура (constant operation) и алгоритм хорошо масштабируется как и standard shadow map. Искривление содержит трансформацию перспективы, где появляется проблема polygon offset. Мы справляемся с ней используя fragment programs современных графических карт.

Аппроксимация трапеции

Раздел 6 в статье по TSM детально рассматривает аппроксимация трапеции. Для большей ясности мы приводим эту ссылку на небольшой видеоролик, который визуализирует алгоритм.

Трансформация трапеции

Мы хотим вычислить трансформацию N_T (матрица 4x4), которая наложит четыре угла трапеции t_0, t_1, t_2 и t_3 на переднюю сторону единичного куба, т.е. мы хотим вычислить N_T со следующими условиями:

(-1, -1, 1, 1)^T = N_T * t_0
(+1, -1, 1, 1)^T = N_T * t_1
(+1, +1, 1, 1)^T = N_T * t_2
(-1, +1, 1, 1)^T = N_T * t_3



Есть несколько способов добится этого. Прямолинейный метод это применить к трапеции вращение, сдвиг, искажена (shearing*), масштабирование и нормализацию, чтобы наложить ее на переднюю часть единичного куба. Это достижимо (в теории) вычислением восьми матриц 4x4 (T_1, R, T_2, H, S_1, N, T_3 и S_2):

1. Первым шагом, T_1 перенесем центр верхней стороны в начало координат (вектора u = (x_u, y_u, z_u, w_u) и v = (x_v, y_v, z_v, w_v) hold intermediate results):

Код:

u = (t_2 + t_3) / 2,       |1 0 0 -x_u| T_1 = |0 1 0 -y_u|       |0 0 1   0 |       |0 0 0   1 |

2. Потом поворачиваем трапецию T вокруг начала координат применяя R таким образом, что верхняя сторона становится коллинеарной оси X:

Код:

u = (t_2 - t_3) / |t_2 - t_3|,     |x_u  y_u 0 0| R = |y_u -x_u 0 0|     | 0    0  1 0|     | 0    0  0 1|

3. После шага 2, пересечение прямых содержащих боковые стороны (t_0, t_3) и (t_1, t_2), обозначенное как i, переносится применяя T_2 к началу координат:

Код:

u = R * T_1 * i,       |1 0 0 -x_u| T_2 = |0 1 0 -y_u|       |0 0 1   0 |       |0 0 0   1 |

4. Следующим шагом, трапеция должна быть искажена (sheared*) с помощью H так, что она станет симметрична оси Y, т.е. прямая проходящая через середины верхней и нижней сторон буде коллинеарна оси Y:

Код:

u = (T_2 * R * T_1 * (t_2 + t_3)) / 2,     |1 -x_u/y_u 0 0| H = |0     1    0 0|     |0     0    1 0|     |0     0    0 1|

5. Теперь применив S_1 масштабируем трапецию по осям так, что угол между прямыми содержащими боковые стороны (t_0, t_3) и (t_1, t_2) станет 90 градусов, и так, что расстояние между верхней гранью и осью X станет 1:

Код:

u = H * T_2 * R * T_1 * t_2,       |1/x_u   0   0 0| S_1 = |  0   1/y_u 0 0|       |  0     0   1 0|       |  0     0   0 1|

6. Следующим шагом применив N мы трансформируем трапецию в прямоугольник:

Код:

|1 0 0 0| N = |0 1 0 1|     |0 0 1 0|     |0 1 0 0|

7. Потом прямоугольник смещаем вдоль оси Y пока его центр не совпадет с началом координат. Это достигается применением T_3. После этой трансформации прямоугольник также симметричен и относительно оси X:

Код:

u = N * S_1 * H * T_2 * R * T_1 * t_0, v = N * S_1 * H * T_2 * R * T_1 * t_2,       |1 0 0          0          | T_3 = |0 1 0 -(y_u/w_u+y_v/w_v)/2|       |0 0 1          0          |       |0 0 0          1          |

8. В качестве последнего шага прямоугольник должен быть масштабирован S_2 по оси Y так, чтобы он покрыл переднюю часть единичного куба:

Код:

u = T_3 * N * S_1 * H * T_2 * R * T_1 * t_0,       |1     0     0 0| S_2 = |0 -w_u/y_u  0 0|       |0     0     1 0|       |0     0     0 1|

Итак, трансформация трапеции N_T может быть вычислена таким образом:

Код:

N_T = S_2 * T_3 * N * S_1 * H * T_2 * R * T_1 .

Далее приведен оптимизированный C-код для вычисления проецрования трапеции (trapezoidal mapping) описаного выше:

Код:

#define ASSIGN_MAT(M, u0, u3, u6, u1, u4, u7, u2, u5, u8) { \    M[0][0] = u0; M[0][1] = u3; M[0][2] = u6; \    M[1][0] = u1; M[1][1] = u4; M[1][2] = u7; \    M[2][0] = u2; M[2][1] = u5; M[2][2] = u8; \    } #define DET2(a, b, c, d) ((a) * (d) - (b) * (c)) #define DOT2(u, v) (u[0] * v[0] + u[1] * v[1]) void intersect(float i[2], float g0[3], float g1[3], float h0[3], float h1[3]) {    float a, b;    i[0] = i[1] = 1.0f / DET2(g0[0] - g1[0], g0[1] - g1[1], h0[0] - h1[0], h0[1] - h1[1]);    a = DET2(g0[0], g0[1], g1[0], g1[1]);    b = DET2(h0[0], h0[1], h1[0], h1[1]);    i[0] *=   DET2(a, g0[0] - g1[0], b, h0[0] - h1[0]);    i[1] *=   DET2(a, g0[1] - g1[1], b, h0[1] - h1[1]); } void map_Trapezoid_To_Square(float TR[3][3], float t0[3], float t1[3], float t2[3], float t3[3]) {    float i[2], a, b, c, d;    //M1 = R * T1    a = 0.5f * (t2[0] - t3[0]);    b = 0.5f * (t2[1] - t3[1]);    ASSIGN_MAT(TR, a  ,  b  , a * a + b * b,                   b  , -a  , a * b - b * a,                  0.0f, 0.0f, 1.0f);    //M2 = T2 * M1 = T2 * R * T1    intersect(i, t0, t3, t1, t2);    TR[0][2] = -DOT2(TR[0], i);    TR[1][2] = -DOT2(TR[1], i);    //M1 = H * M2 = H * T2 * R * T1    a = DOT2(TR[0], t2) + TR[0][2];    b = DOT2(TR[1], t2) + TR[1][2];    c = DOT2(TR[0], t3) + TR[0][2];    d = DOT2(TR[1], t3) + TR[1][2];    a = -(a + c) / (b + d);    TR[0][0] += TR[1][0] * a;    TR[0][1] += TR[1][1] * a;    TR[0][2] += TR[1][2] * a;    //M2 = S1 * M1 = S1 * H * T2 * R * T1    a = 1.0f / (DOT2(TR[0], t2) + TR[0][2]);    b = 1.0f / (DOT2(TR[1], t2) + TR[1][2]);    TR[0][0] *= a; TR[0][1] *= a; TR[0][2] *= a;    TR[1][0] *= b; TR[1][1] *= b; TR[1][2] *= b;    //M1 = N * M2 = N * S1 * H * T2 * R * T1    TR[2][0] = TR[1][0]; TR[2][1] = TR[1][1]; TR[2][2] = TR[1][2];    TR[1][2] += 1.0f;    //M2 = T3 * M1 = T3 * N * S1 * H * T2 * R * T1    a = DOT2(TR[1], t0) + TR[1][2];    b = DOT2(TR[2], t0) + TR[2][2];    c = DOT2(TR[1], t2) + TR[1][2];    d = DOT2(TR[2], t2) + TR[2][2];    a = -0.5f * (a / b + c / d);    TR[1][0] += TR[2][0] * a;    TR[1][1] += TR[2][1] * a;    TR[1][2] += TR[2][2] * a;    //M1 = S2 * M2 = S2 * T3 * N * S1 * H * T2 * R * T1    a = DOT2(TR[1], t0) + TR[1][2];    b = DOT2(TR[2], t0) + TR[2][2];    c = -b / a;    TR[1][0] *= c; TR[1][1] *= c; TR[1][2] *= c; }

TR это матрица 3x3. Чтобы использовать ее (например) в Opengl, ее надо сконвертировать в матрицу 4x4:

...
float TR[3][3];

map_Trapezoid_To_Square(TR, t0, t1, t2, t3);

GLdouble N_T[16];

N_T[0] = TR[0][0]; N_T[4] = TR[0][1]; N_T[ 8] = 0.0f; N_T[12] = TR[0][2];
N_T[1] = TR[1][0]; N_T[5] = TR[1][1]; N_T[ 9] = 0.0f; N_T[13] = TR[1][2];
N_T[2] = 0.0f; N_T[6] = 0.0f; N_T[10] = 1.0f; N_T[14] = 0.0f;
N_T[3] = TR[2][0]; N_T[7] = TR[2][1]; N_T[11] = 0.0f; N_T[15] = TR[2][2];
...

Другой более общий метод заключается в том, чтобы произвести вычисления quadrilateral to quad mapping, которые могут быть найдены в "Master Thesis: Fundamentals of Texture Mapping and Image Warping" Paul Heckbert'а на стр. 17-21, и исходный код в Appendix A.2.2, стр 70-72.

Аппроксимация трапеции и вычисление N_T - это единственная часть алгоритма, где используется CPU.

Работа над решением проблемы Polygon Offset

Трансформация трапеции включает в себя двумерное проецирование. Важным свойством этой трансформации является то, что z_T вершины в трапециоидном пространстве зависит от w_T. По факту, распределение значений z меняется по всей трапецевидной карте тени (shadow map) таким образом, что выбор неизменного polygon offset, такого как в методе standard shadow map, более не является приемлимым. Проблема в том, что конкретный polygon offset может быть слишком велик для пикселей расположенных ближе к камере или может быть слишком мал для пикселей, которые расположен дальше. Если polygon offset слишком велик может случиться так, что будут исчезать тени; с другой стороны, если он слишком мал могут проявится черные точки на поверхностях (surface acne - угри на поверхности - прим. пер.). В статье мы рекомендуем сохранять значения z в пространстве источника света после вычисления перспективы (post-perspecive space of the light**), чтобы мог быть найден polygon offset похожий на тот, что используется в методе standard shadow map.


Иллюстарция: N_T влияет на распределение значений z. Значение конкретного polygon offset слишком велико для областей ближе к камере, что приводит к исчезновению теней. Уменьшая polygon offset, будет привнесено неправильное самозатенение в областях дальше от камеры.

Иллюстрация: Сохраняя значения глубины в пространстве источника света после вычисления перспективы (post-perspecive space of the light**), можно указать постоянный polygon offset как в методе standard shadow map. Распределение останется равномерным.

Далее мы покажем как может быть реализована трансформация трапеции (используя компоненты NV_vertex_program и NV_fragment_program) таким образом, что можно избежать серьезных проблем с polygon offset. Мы покажем только ту часть, где генерируется TSM, т.е. первый проход алгоритма, потому что второй проход работает таким же образом (Трансформируем первую координату текстуры в трапециоидное пространство, а вторую - в пространство источника света после вычисления перспективы (post-perspecive space of the light**). Потом, на стадии работы с fragment заменяем z_T на z_L)***.
(Transform the first texture coordinate to trapezoidal space and the second to the post-perspective space of the light. Then, on the fragment stage replace z_T with z_L).

light's vertex program: ****

!!VP1.0

# c[0-3] : N_T
# c[8-11] : матрица проецирования света и modelview (light's projection and modelview matrix)
# c[12-15]: матрица мира (world matrix)
# c[16-19]: инвертированная матрица мира
# v[OPOS] : позиция объекта
# o[HPOS] : результирующая вершина (result vertex)

# трансфоримруем v[OPOS] в пространство мира и сохраняем результат в R4:
# R4 = W * v[OPOS]
DP4 R4.x, c[12], v[OPOS];
DP4 R4.y, c[13], v[OPOS];
DP4 R4.z, c[14], v[OPOS];
DP4 R4.w, c[15], v[OPOS];

# трансформируем R4 в пространство источника света после вычисления перспективы (light's post-perspective space**) и сохраняем результат в R1:
# R1 = P_L * C_L * (R4) = P_L * C_L * W * v[OPOS]
DP4 R1.x, c[8], R4;
DP4 R1.y, c[9], R4;
DP4 R1.z, c[10], R4;
DP4 R1.w, c[11], R4;

# сохраняем R1 в first texture coordinate***:
# o[TEX0] = P_L * C_L * W * v[OPOS]
MOV o[TEX0], R1;

# трансформируем R4 в трапециоидное пространство:
# o[HPOS] = N_T * P_L * C_L * W * v[OPOS]
DP4 o[HPOS].x, c[0], R1;
DP4 o[HPOS].y, c[1], R1;
DP4 o[HPOS].z, c[2], R1;
DP4 o[HPOS].w, c[3], R1;

MOV o[COL0], v[COL0];

END

light's fragment program:

!!FP1.0

# f[WPOS] : fragment в трапециоидном пространстве (положение окна (window position))
# f[TEX0] : положение fragment в пространстве источника света после вычисления перспективы (post-perspective space of the light**)

RCP R0.x, f[TEX0].w; # R0.x = 1 / w_L
MUL R1, f[TEX0], R0.x; # R1 = (x_L/w_L, y_L/w_L, z_L/w_L, 1)

MAD R2.z, R1.z, 0.5, 0.5; # R2.z = R1.z * 0.5 + 0.5; depth is now in the range [0;1]

MOV o[DEPR], R2.z; # replace "z_T" with "z_L"

MOV o[COLR], f[COL0];

END


Обратите внимание, что для ясности мы опустили два шейдера перел вычислением постоянного polygon offset, который добавляется к финальному значению глубины.

Может сложится неправильное впечатление, что TSM работает только на том железе, которое поддерживает "xxx"_fragment_program. Используя эти простые шейдеры мы работаем над решением проблемы polygon offset, которая ухудшилась из-за проективной природы N_T. (Заметим, что другие методы также используют 2D/3D трансформацию перспективы чтобы улучшить качество теней (как PSM). Поэтому, Для решения проблемы polygon offset в других методах тоже могут быть необходимы fragment/vertex programs. Однако в нашем случае проблема может быть решена очень просто.)

Если графическая карта не поддерживает "xxx"_fragment_program, может быть применен другой метод для "решения" проблемы polygon offset. Вместо fragment program можно использовать vertex program для сохранения значений z вершины в пространстве источника света после вычисления перспективы (post-perspective space of the light**). Это достигается в результате v_T = (x_T, y_T, (z_L*w_T)/w_L, w_T), т.е. v_T = (x_T/w_T, y_T/w_T, z_L/w_L, 1) после homogeneous division. Заметим, что это не общее решение, т.к. это ведет к неправильной интерполяции во внутренних частях треугольнико: Это корректно на уровне вершин, но не на уровне пикселей. However, если сцена достаточно тесселирована, то видимых артефактов не будет.

Используя эти шейдеры и трансформацию трапеции N_T как приведено выше, программа для отображения может выглядеть примерно так:

Код:

void display() {    // 1-ый проход: генерация TSM    // как описано в разделе Section 6 в статье TSM,    // в качестве первого шага мы должны вычеслить трапецию основываясь    // на камере E.    // Четыре вершины хранятся в t_0, t_1, t_2, t_3    calculateTrapezoid(&t_0, &t_1, &t_2, &t_3, P_L, C_L, E);    // ...после этого N_T вычисляется как описано выше    calculateTrapezoidalTransformation(&N_T, t_0, t_1, t_2, t_3);    glViewport(0, 0, shadowmapWidth, shadowmapHeight);    // Bind the above vertex program...    glBindProgramNV(GL_VERTEX_PROGRAM_NV, vpLight);    // ...and bind the above fragment program    glBindProgramNV(GL_FRAGMENT_PROGRAM_NV, fpLight);    glMatrixMode(GL_PROJECTION); // store N_T in GL_PROJECTION    glLoadMatrixf(N_T);    glTrackMatrixNV(GL_VERTEX_PROGRAM_NV, 0, GL_PROJECTION, GL_IDENTITY_NV);    glMatrixMode(GL_MATRIX1_NV); // store in P_L * C_L and track in GL_MATRIX1_NV    glLoadMatrixf(P_L);          // матрица световой проекции (light's projection matrix), полученная, например, от gluPerspective()    glMultMatrixf(C_L);          // матрица камеры света (light's camera matrix), полученная, например, от gluLookAt()    glTrackMatrixNV(GL_VERTEX_PROGRAM_NV, 8, GL_MATRIX1_NV, GL_IDENTITY_NV);    glMatrixMode(GL_MODELVIEW);  // сохраняем матрицу modelview в GL_MODELVIEW    glLoadIdentity();    glTrackMatrixNV(GL_VERTEX_PROGRAM_NV, 12, GL_MODELVIEW, GL_IDENTITY_NV);    renderScene();    // Как и в методе standard shadow map мы копируем depth buffer в текстуру:    CopyDepthBufferToTexture();    ...    // 2-ой проход: рендер сцены с позиции камеры и проецируем текстуру TSM на сцену    ...    SwapBuffers(); }

Tips and Tricks

• Обычный алгоритм проецирует focus region на линию 80%. Статья также упоминает итератиный метод для поиска оптимальной процентной линии в целом. Такой поиск должен быть затратным для современных CPU. С другой стороны, как указано в статье, можно заранее вычислить нужный процент для каждой конфигурации камера/свет и сохранить их в таблицу (это работает потому что это не зависит от сцены). Также можно сберечь часть вычислительной мощности ограничивая поиск (сильно уменьшив от 80%) после того, как (локальная) максимальная область была обнаружена. В наших тестовых сценах (фэнтэзи и городская), мы не выполняли пред-вычисления нужного процента, но мы применили упомянутый критерий остановки, и скорость работы TSM сравнима с SSM.
  
• Можно применить неквадратную текстуру (текстуру с различными длиной и шириной) для уменьшения aliasing в направлении Z.
  
• Иллюстрация 6 в статье показывает кривую над областью focus region на карте тени. Эта область уменьшается, когда камера и свет приближаются к ситуации dueling frusta. Это трудная ситуация для всех методов генерации карт тени. Кривая TSM (Иллюстрация 6 в статье) может быть использована, чтобы "адаптивно" постепенно менять размер текстуры (сохраняя при этом непрерывность), чтобы получить бОльшую текстуру для ситуаций dueling frusta и меньшую для других случаев (требует поддержки текстур размеров кроме степеней двойки).
  
• Использовать только vertex program для решения проблемы polygon offset. Однако у этог есть ограничения (как было описано выше и в статье).

Для вопросов и предложений пожалуйста свяжитесь tobiasmartin@t-online.de или tants@comp.nus.edu.sg. Вот часть исходного кода, которая может помочь в разработке по приведенному методу.