суббота, 31 октября 2009 г.

Column-major vs Row-major

Теперь жалею, что когда-то решил использовать column-major матрицы в своей математической библиотеке. Во-первых, теперь приходится транспонировать их при передаче в D3D-шейдер, чтобы умножение вектора на матрицу можно было записывать как mul(v, M) (Direct3D), а не mul(M, v) (OpenGL). В первом случае это ложится на четыре красивые dp4 инcтрукции, во-втором - на комбинацию mul/mad. Во-вторых, для column-major матриц приходится записывать перемножение по правилу "post-multiplying":

WorldViewProj = Proj * View * World
HPos = WorldViewProj * Pos

Я не араб, и меня этот способ начал со временем напрягать. Придётся переписывать всю матричную библиотеку - ошибки молодости :( Думаю через шаблоны сделать возможность указать, какая именно матрица собирается использоваться - row-major или column-major, и если идёт присваивание матриц с разным порядком записи, чтобы осуществлялось автоматическое траспонирование:

Mat4< column_major > m1 = GetWorldViewProj();
Mat4< row_major > m2( m1 ); // transpose

среда, 28 октября 2009 г.

Ray-BV Intersection

Сегодня весь день обдумывал, в какую сторону двигаться дальше.

[Здесь была неверная инфа :)].
Полдня лазания по Гуглу и отладки кода, и в шейдеры легли максимально эффективные функции, тестирующие пересечения. Пересечение с боксом просто как две копейки. А люди выдумывают какие-то извраты через координаты Плюкера...

Тест на пересечение со сферой - 6 слотов инструкций:



Тест на пересечение с эллипсом - 14 инcтрукций:



Branchless тест на пересечение с AABB - 14 инcтрукций (включая одно деление):



При этом мне не нужно само пересечение, а только булево значение, а значит, нет необходимости в квадратном корне.

вторник, 27 октября 2009 г.

HD 2400 is slow

В общем провёл я предварительные тесты на производительность. Конкретно для Radeon HD 2400 результаты неутешительные - скорость падает приблизительно линейно с возрастанием кол-ва треугольников... Если для квада из двух треугольников получается ~100 fps, то для кубика из 12 треугольников - 20 fps. Для двух кубиков - 11... Скорость мало зависит от площади кадра, которую покрывают кубики, и от типа кушаемых данных - вещественные или байты, хотя при вещественных скорость снижается.

Понятно, откуда линейное падение скорости - от линейного увеличения длины цикла. Я не считал точно, но на глазок длина цикла ~60 инструкций (весь шейдер - 120). 24 треугольника требуют выполнить ~1440 + 60 = ~1500 инструкций. Боттлнек в недостатке вычислительной способности при тупом переборе циклом. Можно взять high-end видеокарту, но это только отодвинет верхнюю границу, после которой симптомы будут аналогичны.

В принципе результаты не так уж плохи (в расчёте на high-end), таким способом можно рассчитывать на возможность играться с несколькими десятками треугольников (мне нужно даже меньше). К тому же на HD 2400 динамический бранчинг нифига не фурычит, а только просаживает скорость :( Хотя очевидно, что он должен давать прирост, например, при cull-инге back-faced треугольников. На нормальном железе с нормальным бранчингом должен быть заметный выигрыш.

Очевидный резерв для таких забав на растеризаторах - использование какого-нить иерархического дерева для геометрии, с обязательным условием, что динамический бранчинг работает эффективно. Это позволит снизить кол-во вычислений, связанных с тестом на пересечение луча и геометрии, правда, появятся затраты на обход дерева. Честно говоря, мне сложно представить, как это можно впихнуть в пиксельный шейдер, и стоит ли - очень быстро всё равно не получится.

Вот как-то так. Но всё же я надеялся на чуть лучшие результаты.

PS. Ну, а что можно сделать с несколькими треугольниками? Поверьте, в умелых руках - очень многое :)

PPS. Вот перед сном думал, как можно ускорить.
1) Я делаю рэйкастинг для всех пикселей фреймбуфера, очевидно что это худший случай. Рэйкастинг - это потому что простейшая реализация. Для рэйтрейсинга легче. Если у нас есть рефрактор, не занимающий всю площадь экрана, то растеризатором сначала рисуем front-faced треугольники, а потом выполняем рэйтрейсинг.
2) Early stencil rejection. Делаем рэйтрейсинг в малю-ю-юсенький render target, разблуриваем получившуюся маску, "натягиваем" её на фреймбуфер и помечаем в стенсиле. Затем делаем нормальный рэйтрейсинг, стенсил тест отбрасывает ненужные фрагменты вне маски.
3) Почему-то подумалось об oct-tree со сферами, описывающими кубики-узлы. Тест на пересечение со сферой тривиален.

понедельник, 26 октября 2009 г.

Raycaster vs Rasterizer

It has begun... Первая реализация рэйкастера на DirectX 11!
Два текстурированных треугольника, диффуз + спекуляр, Radeon HD 2400, окно 1024x768, без мультисэмплинга.

Растеризатор, FPS ~= 390:


Рэйкастинг, FPS ~= 114:


Растеризатор, FPS ~= 490:


Рэйкастер, FPS ~= 123:



Два (нет, три) отличия этого метода от растеризатора - нельзя напрямую использовать MIP-фильтрацию и не работает мультисэмплинг. Ну и скорость :) Шейдер для рэйкастинга + шейдинга занимает сейчас около 100 инструкций. Я ещё пооптимизирую что можно, потом буду мерять, во что это выливается.

Дальше - сжатие данных и тесты. Посмотрим, сколько треугольников удастся отобразить, прежде чем FPS на моём лоу-энде упадёт до отметки 20.

How to pack normal into just 2 bytes

Необходимость в жёсткой оптимизации по чтению из памяти привели меня к упаковке предрасчитанных данных и вершин. Хотел было написать, как мне удалось реализовать упаковку нормали в unorm2 с потерей одного бита точности для хранения z sign. Но, порыскав по гуглу, нашёл великолепный анализ методов сжатия нормалей от Aras Pranckevičius (я и не думал, что их столько, видать deferred техники здорово включают изобреталку:):

Compact Normal Storage for small G-Buffers

Все представленные алгоритмы упаковывают нормаль в unorm2 и восстанавливают потерянную информацию, используя квадратный корень. Метод простой упаковки xy с последующей реконструкцией z оказался наименее качественным (зато наиболее быстрым). Видимо из-за dot(n.xy, n.xy), который усиливает ошибку квантования (а для world нормалей придётся пожертвовать ещё одним битом точности). Любопытый обзор, вероятно лучший в Сети - для поклонников deferred shading-a будет весьма ценнен. Я же собираюсь использовать наиболее быстрый вариант, т. к. у меня нормали будут интерполироваться по треугольнику, и каждую надо распаковать... Сойдёт и удовлетворительное качество, главное - минимум чтений из памяти.

четверг, 22 октября 2009 г.

Woop's unit triangle intersection test

После курения пейпров и вникания в кудовские кернелы реализовал и этот алгоритм. Математика - это вещь! Диву даёшься, как только получается реализовать полную проверку такого пересечения в несколько инструкций.

Кстати я не знаю точно, кто именно автор алгоритма. S. Woop указан в качестве соавтора пейпра, в котором рассматривается математика алгоритма, это правда. Но там же идут отсылки к более ранней работе J. Arenberg 1988 года, и говорится что это просто его расширенная версия. Но везде (презентации, на форумах) упоминается именно Woop.

Моя идея с матрицей преобразования хороша, но не пригодна. Впрочем, преобразование из одной системы координат в другую матрицей афинного преобразования - задача не новая, алгоритм Woop-a использует эту же идею. Проблема в том, что кроме самих текстурных координат мне нужны ещё как минимум интерполированные позиция и нормаль (растеризатор-то теперь мне ничего интерполировать не будет), и искать её удобно через те же барицентрические координаты. И алгоритм Мюллера, и алгоритм Вупа делают проверку через поиск барицентрических координат, а раз они уже есть, то использовать их для интерполяции вершинных значений - наиболее разумное решение.

Пока я фигачил тестовые реализации для этого тормоза CPU, в голове крутились мысли, как это должно быть реализовано на GPU класса DX10/11. Во-первых, зачем таскать текстурные координаты (а в перспективе - и нормали) вместе c precomputed данными? Ведь для поиска самого пересечения они не нужны, а требуются лишь тогда, когда найдено пересечение с треугольником и его индекс известен. Нужно отделить данные для интерполяции от precomputed данных для поиска пересечения, разнести их по двум буферам. Тогда буфер precomputed данных ужмётся до минимума, а значит, и фетчинг данных при линейном поиске по массиву пойдёт быстрее. А когда треугольник найден, то по индексу делаем выборку атрибутов вершин из второго буфера. Хотя тут и есть скачок по памяти, но он всего один на пиксель.

Тогда структура для алгоритма Мюллера должна быть такой:
struct prec_tri
{
float3 v0;
float3 e0;
float3 e1;
};

Тут к сожалению, в cbuffer пропадает впустую 3 вещественных, а это 12 байт из 48. Один из вариантов - использовать Buffer < float3 > и Load(). Для алгоритма Вупа достаточно матрицы 4x3 (четыре столбца по 3 элемента - поворот и перенос). Это те же 48 байт. А вот если мне нужна позиция в точке пересечения (для L, V), то не обязательно тащить вершины - можно попробовать найти по O + t * D, хотя неясно что будет с точностью и визуальным качеством. Но это так, мысли наперёд.

А вот вопрос, как именно фетчить: через cbuffer, tbuffer или Load() - очень интересен. Судя по описаниям в SDK, tbuffer оптимизирован для случайного доступа, cbuffer - для последовательного. Может оказаться, что это "шило на мыло". Можно упереться не в память, а в расчёты, kd-tree то не планирую, цикл соответственно будет только расти с ростом кол-ва треугольников. Да и моя HD 2400 очень тормозная железяка, грех на ней рэйтресингом баловаться, а может заниматься оптимизациями - в самый раз? :). Думаю, не купить ли на последние деньги HD 4770 или 4850, в районе 100$ сейчас...

Кстати, привет от компайлера шейдеров! Что-то вроде такого он пережёвывает секунд 10:

#define MAX_TRIANGLES 200

[loop]
for (int i = 0; i less than MAX_TRIANGLES; ++i)
{
prec_tri t = g_buf[ i ];
// Дальше куча расчётов.
}

Если взять, скажем, 500 - думает секунд 40. Ассемблерный листинг показывает, что цикл не разворачивается. Может компайлер проверяет валидность i для всех итераций цикла? Пробовал разные флаги подсовывать, убирать оптимизацию - бесполезно. Единственное, что помогло - замена compile-time MAX_TRIANGLES на значение из буфера констант. Тормоза пропадают. Какой из этого вывод? Незнание - это блаженство (с) The Matrix.

вторник, 20 октября 2009 г.

GPU Ray-Triangle Intersection

Общепризнанным стандартом здесь является алгоритм Мюллера-Трумбора:

Fast, Minimum Storage Ray/Triangle Intersection.

Написан давно (1997 г), но хорош и подходит для GPU. Принцип работы - поиск барицентрических координат пересечения u и v, имея которые, можно легко найти точку пересечения или её текстурные координаты:
w = 1 - u - v
tc = tc0 * w + tc1 * u + tc2 * v
В CPU версии на вход подаются начало луча, его направление и три вершины треугольника. На выходе - скаляр t на луче (наименьший означает ближайший треугольник), и uv. После этого, проверив (t < min_t), надо найти собственно сами текстурные координаты точки пересечения, для этого нужны текстурные координаты вершин треугольника.

Одними из факторов, серьёзно ударяющих по производительности RTI, являются фетчинг данных для треугольника и загруженность регистров (регистрового файла GPU):

Основная трудность трассировки лучей на GPU

И хотя это становится проблемой при больших массивах трассируемых данных, это своего рода указатель и для ограниченного случая. Что касается фетчинга данных, то статические треугольники должны браться из буфера констант (идём тупым перебором по всему массиву, упаси господи юзать какие-то kd-tree). В базовом варианте нам нужно:
struct tri
{
float3 v0;
float3 v1
float3 v2;
float2 tc0;
float2 tc1;
float2 tc2;
};
Для буферов констант важен alignment, поэтому предъявляется требование к выравниванию вектора по границе, кратной 16 байт (float4) (NOTE: Более точно, данные вектора не должны пересекать границу, кратную 16 байтам). Поэтому размер структуры - sizeof(float4) * 5 = 80 байт (два float2 группируются в один float4). Много, хотя хорошо укладываемся в 16-byte alignment. Попробуем сжать.

Как видно из кода алгоритма (ссылку приводил в начале), из трёх вершин треугольника используется только один, остальные нужны для нахождения векторов рёбер треугольника:
edge1 = v1 - v0
edge2 = v2 - v0
В структуру вместо вершин! Заодно осводождаем два слота инструкций. Третьи текстурные координаты распихиваем по w компонентам векторов рёбер. Итого получаем:
struct prec_tri
{
float3 v0;
float4 e0;
float4 e1;
float4 tc0tc1;
};
Укладываемся в 64 байта.

В героической попытке придумать алгоритм покороче и избавиться от поиска барицентрических координат, я решил искать текстурные координаты через матрицу трансформации. Представим, что e0 и e1 формируют не ортонормированный базис в object space, тогда:
t0 = tc1 - tc0
t1 = tc2 - tc0
формируют базис в пространстве текстуры. Третий орт не нужен, т. к. текстурная плоскость совпадает с плоскостью треугольника. Тогда:
Muv = Mobj * Mx
откуда:
Mx = Mobj^-1 * Muv
(запись для row-major матриц). Из-за того, что наш треугольник представлен в object space, нужно переводить точку пересечения в локальные координаты треугольника, а потом переводить в локальные текстурные координаты:
tc = (p - v0).xy * Mx + tc0
Достаточно оперировать матрицами второго порядка, т. к. фактически из всех афинных преобразований используются поворот и масштабирование на плоскости. Но это не всё. Т. к. оси Oy и Ot в Direct3D направлены в разные стороны, нужна матрица преобразования координаты t:
t = t * (-1) + 1
Кроме того, можно не хранить текстурные tc0, а записать их как translation в матрицу 2x3.
Итого:
Muv * Mflip = Mobj * Mx
Mx = Mobj^-1 * Muv * Mflip
tc = ((p - v0).xy, 1) * Mx
Естественно, что матрица предрассчитывается и заливается в буфер констант. Матрица 2x3 займёт 2 * 4 * 4 = 32 байта (NOTE: Впрочем, матрица 2x2 укладывается во float4, и если tc0 дополнить какими-то полезными данными до float4, матрицу 2x3 использовать нет смысла). Я протестировал этот метод в случае, если известны координаты пересечения луча с треугольником - работает отлично. Единственное, что ускользнуло от меня - вместо матрицы флипа t = t * (-1) + 1 удаётся использовать t = t * (-1). Загадка :)

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

Одной из самых перспективных альтернатив алгоритму Мюллера является Woop’s unit triangle test. Вот что удалось накопать по этой теме:

Realtime Ray Tracing of Dynamic Scenes on an FPGA Chip
RPU: A Programmable Ray Processing Unit for Realtime Ray Tracing
Пост whiteambit на ompf.org

Сейчас копаю эту тему. По идее, в unit space отпадает необходимость в (p - v0). Не знаю, что получится, но поиск текстурных координат пересечения по двум матрицам - это жёстко :)