Каустика в реальном времени

Автор: Сергей Резник

Данная статья расскажет о создании и отрисовке каустики от прозрачных и полупрозрачных объектов.

Введение
Физическая основа
Предварительные этапы отрисовки
Раcчет каустики
  Расчет каустики при одном преломлении
  Расчет каустики при двух преломлениях
  Дисперсия света
Заключение
Ссылки

Введение

Каустика значительно улучшает восприятие прозрачных и полупрозраных объектов в сцене. Конечно, для игр этот эффект реализовывать пока что дорого (в вычислительном плане), но для демок и различных заставок (типа хранителей экранов) такой эффект может быть очень полезен. Реализовать его достаточно просто — рисование непосредственно каустики идет всего в три прохода (без учета финальной отрисовки).

Данный метод чем-то напоминает технику фотонных карт (photon mapping), поэтому в данном контексте можно говорить о фотонах. И в то же время, метод похож на теневые карты (shadow map). Общая идея метода такова — мы испускаем фотон (или луч света, если угодно) и далее вычисялем его искаженную преломлением на границе прозрачного объекта траекторию. После чего мы находим пересечение этой траектории (луча света) с геометрией. Затем мы сохраняем в текстуру положение этого фотона. Рассчитав таким образом положения n фотонов мы получим текстуру, содержащую «узор» каустики. Теперь нам остается наложить эту текстуру на геометрию (в этом состоит похожесть метода с теневыми картами)

Итак, для создания этого эффекта (в самом простом случае) нам потребуются следующие этапы отрисовки (проходы):
1) запись в текстуру положения и нормалей преломляющего объекта (далее refractive — преломляющий объект),
2) запись в текстуру положений объектов, на которые будет падать преломленный свет (далее receiver — принимающая геометрия),
3) непосредственно расчет каустики.

Три прохода идеально подходят для создания каустики от водной поверхности. Для отрисовки объектов желательно использовать еще один проход. Но об этом чуть позже.

Физическая основа

Прежде чем приступать к реализации данного эффекта, думаю, полезно будет разобраться, что же все-таки мы моделируем. Здесь я дам только краткое объяснение происходящего без каких-либо формул.

Итак, как известно лучи света при переходе из одной среды в другую (например из воздуха в стекло) преломляются по закону Снелла. При этом некоторые лучи, оклонившить от первоначальной траектории рассеиваются, а некоторые концентрируются на небольших площадях. При этом появляются яркие узоры на поверхности в точке концентрации преломленных лучей. То же самое происходит и при концетнтрации отраженных лучей. Эти яркие узоры и называются каустикой. Наблюдать каустику можно повсюду — отражения на дне фарфоровой чашки, световое пятно от линзы для выжигания, радуга и т.д.

сконцентрированные лучи света после прохождения прозрачной сферы

Предварительные этапы отрисовки

Давайте же приступим к созданию этого эффекта. В простейшем случае нам потребуется сделать два предварительных прохода отрисовки. Причем, если у вас есть уже теневые карты, то проход остается всего один. Теперь подробнее.

В первом проходе (они никак не связанны, и первым может быть любой, но я для определенности буду иметь ввиду первым тот, о котором речь пойдет дальше) нам потребуется получить позиции геометрии, на которую в последствии будет проецироваться каустика. В оригинальном документе, на основе которого я делал программу, предлагалось использовать float-текстуры для записи положения геометрии в world-space. Я же предлагаю ограничиться стандартной записью в глубину с позиции источника света (карты теней, привет!) геометрии и восстанавливать потом положение геометрии путем домножения на обратную матрицу. Так как нам понадобится всего несколько таких умножений на один «фотон», то сильно накладно это не выйдет. Еще одним дополнительным плюсом такого подхода будет скорость этого прохода — записывать будем только в глубину, не трогая цвет. Для этого нам понадобится отключить запись в цвет для framebuffer'а:

  glReadBuffer(GL_NONE);
  glDrawBuffer(GL_NONE);

Итак, первый проход — отрисовка с позиции источника света «принимающей» геометрии без записи в цвет.

Во втором проходе нам потребуется записать положение и нормали объектов, которые будут преломлять свет. В оригинальном документе предлагают делать два float rendertarget'a, но это слишком жирно. Мы также ограничимся записью в глубину для восстановления положения в world-space и одним простым RGBA8 rendertarget'ом для нормалей. Точно так же рисуем преломляющую геометрию с позиции источника света с примерно вот таким шейдером:

// VERTEX
#version 150
precision highp float;

uniform mat4 mModelViewProjection;

in vec4 Vertex;
in vec3 Normal;
out vec3 vNormalWS;

void main()
{
 vNormalWS = Normal;
 gl_Position = mModelViewProjection * Vertex;
}

Во фрагментный шейдер будем передавать показатель преломления материала. Это делается для того, чтобы каждый объект на сцене мог «быть сделан» из разных материалов, и в то же время обрабатывался одним проходом рассчета каустики.

// FRAGMENT
#version 150
precision highp float;

uniform float fRefractionIndex;

in vec3 vNormalWS;
out vec4 FragColor;

void main()
{
 FragColor = vec4(vec3(0.5) + 0.5 * vNormalWS, 1.0 / fRefractionIndex);
}

Раcчет каустики

Итак, после первых двух проходов у нас уже есть достаточно данных для того, чтобы перейти непосредственно к расчету карты каустики.

Основная идея алгоритма заключается в том, что все расчеты происходят в вершинном шейдере. Как говорилось раньше, метод похож на photon mapping, следовательно нам нужно будет создать эти фотоны. Для этого мы создаем вершинный буфер, который состоит только из позиций вершин (т.е. без нормалей и текстурных координат), расположенных в квадрате [0..1]x[0..1]. Этот буффер мы будем рисовать точками (GL_POINTS). Создание такого буфера выглядит примерно так:

 // photonMapSize - разрешение конечной текстуры каустики
 int numPhotons = sqr(photonMapSize); 

 vec2* photons = new vec2[numPhotons];
 Index* indices = new Index[numPhotons];
 int k = 0;
 for (int i = 0; i < photonMapSize; ++i)
  for (int j = 0; j < photonMapSize; ++j)
  {
   photons[k] = vec2( float(j) / (photonMapSize), float(i) / (photonMapSize) );
   indices[k] = k;
   k++;
  }
 
 // процедуру создания вершинного буфера считать псевдокодом
 photonBuffer = createVertexBuffer(numPhotons, photons, numPhotons, indices, GL_POINTS);

Далее в вершинном шейдере нам нужно будет определить траекторию луча света и найти пересечение этого луча с принимающей геометрией. Причем пересечения находятся только в image-space. За счет этого мы получаем большую производительность (нет необходимости «знать» принимающую геометрию повершинно и находить пересечения луч—треугольник).

Общий алгоритм нахождения пересечений в images-space такова:

1) находим начальную точку и направление вектора (луча света), для начала берем единичную длинну вектора,
2) проецируем этот вектор на текстуру (в квадрат [0..1]x[0..1]) — получаем из этого вектора текстурные координаты,
3) по полученным текстурным координатам читаем из предварительно заготовленной текстуры (из первого прохода) глубину и восстанавливаем положение прочитанной точки в world-space,
4) находим расстояние между начальной и прочитанной точкой,
5) если расстояние больше заданного, то домножаем начальный вектор на это расстояние. Переходим к пункту два,
6) если расстояние небольшое — значит мы нашли пересечение.

Обычно для нахождения пересечения хватает всего двух итераций.

Расчет каустики при одном преломлении

Таким образом, алгоритм расчета карты каустики для одного преломления выглядит так:

• включаем аддитивное смешивание (BlendFunc(GL_ONE, GL_ONE)),
• рисуем созданный ранее буфер в текстуру с описанным далее шейдером,
• на входе у нас три текстуры (глубина принимающей геометрии, глубина и нормали преломляющей геометрии).

Для каждой точки (вершины) из вершинного буфера выполняем следующий алгоритм в вершинном шейдере:

1) читаем из текстуры нормаль (и вместе с ней и показатель преломления),
2) читаем и восстанавливаем положение точки на преломляющей геометрии,
3) находим направление луча света,
4) с помощью встроенной функции refract находим направление преломленного вектора (луча света),
5) полагаем начальной точкой луча света — восставновленное положение преломляющей геометрии,
6) находим пересечение луча света с принимающей геометрией,
7) проецируем найденную точку пересечения на текстуру.

Теперь немного шейдерного кода:

• восстановление положения в world-space по экранным координатам и глубине:

vec3 restoreWorldPosition( vec2 texCoords, float depth )
{
 vec4 vWorldSpaceVertex = mInverseModelViewProjection * vec4(texCoords, depth, 1.0);
 return vWorldSpaceVertex.xyz / vWorldSpaceVertex.w;
}

• получение положения принимающей геометрии по текстурным координатам:

vec3 sampleReceiverPosition(vec2 texCoords) 
{
 float fSampledDepth = texture(receiver_depth, vec2(0.5) + 0.5 * texCoords).x;
 return restoreWorldPosition(texCoords, 2.0 * fSampledDepth - 1.0 );
}

• получение положения преломляющей геометрии по текстурным координатам:

vec3 sampleRefractivePosition(vec2 texCoords)
{
 float fSampledDepth = texture(refractive_depth, texCoords).x;
 return restoreWorldPosition( 2.0 * texCoords - vec2(1.0), 2.0 * fSampledDepth - 1.0 );
}

и теперь непосредственно сам шейдер (без описанных выше ф-ций):

#version 150
precision highp float;

uniform sampler2D refractive_normals;
uniform sampler2D refractive_depth;
uniform sampler2D receiver_depth;

uniform vec3 vCamera;
uniform mat4 mModelViewProjection;
uniform mat4 mInverseModelViewProjection;

out vec4 cLightColor;
in vec2 Vertex;

vec3 restoreWorldPosition( vec2 texCoords, float depth );

vec3 sampleReceiverPosition(vec2 texCoords);
vec3 sampleRefractivePosition(vec2 texCoords);

// ф-ция, которая проецирует вектор и возвращает
// восстановленное положение принимающей геометрии
vec3 project(vec3 position) 
{
 vec4 vProjected = mModelViewProjection * vec4(position, 1.0);
 return sampleReceiverPosition(vProjected.xy / vProjected.w);
}

// нахождение пересечения:
vec3 estimateIntersection(vec3 startPoint, vec3 ray)
{
 return project( startPoint + ray * distance( startPoint, project(startPoint + ray)));
}

void main()
{
 // получение нормали
 vec4 vSampledNormal = texture(refractive_normals, Vertex);
 // получение положения преломляющей геометрии
 vec3 vRefractivePosition = sampleRefractivePosition(Vertex);
 // приведение нормали к промежутку [-1..1]
 vec3 vNormal = 2.0 * vSampledNormal.xyz - vec3(1.0);

 // нахождение направления луча света
 vec3 vLightVector = normalize(vRefractivePosition - vCamera);

 float eta = vSampledNormal.a;
 vLightVector = refract(vLightVector, vNormal, eta);
 vec3 vIntersectedVertex = estimateIntersection( vRefractivePosition, vLightVector );

 cLightColor = vec4(1.0);
 gl_Position = mModelViewProjection * vec4(vIntersectedVertex, 1.0);
}

Фрагментный шейдер в данном случае простой до безобразия:

#version 150
precision highp float;

in vec4 cLightColor;
out vec4 FragColor;

void main()
{
 FragColor = cLightColor;
}

Расчет каустики при двух преломлениях

Для моделирование перехода лучей света из одной среды в другую (например воздух-вода) нам достаточно было расчитать одно преломление, но для моделирования прохождения света через прозрачные объекты нам необходимо просчитывать два преломления: на «входе» в объект и на «выходе» их него. Суть метода при этом не меняется. Просто теперь нам нужно сделать два предварительных прохода для преломляющей геометрии вместо одного. В первом проходе мы запишем в текстуру передние грани объектов, во втором — задние. Делается это легко, для первого прохода мы выставляем отсечение задних граней:

glCullFace(GL_BACK);
glEnable(GL_CULL_FACE);

для второго прохода — отсечение передних граней:

glCullFace(GL_FRONT);
glEnable(GL_CULL_FACE);

Не забудьте потом выключить или вернуть отсечение задних граней ;)

Также, следует учесть такое явление, как полное внутреннее отражение. При углах, превышающих arcsin(n2 / n1), где n1 мы полагаем единице (вакуум или воздух), а n2 — показатель преломления нашего материала, луч света полностью отражается внутри материала. Чтобы не возиться с углами, и использовать только стандартные средства, сделаем после второго преломления следующую проверку (встроенная ф-я refract возвращает vec3(0.0) в случае полного внутреннего отражения):

 if (dot(vLightVector, vLightVector) == 0.0)
 {
  cLightColor = vec4(0.0);
  gl_Position = vec4(Vertex, 0.0, 1.0);
  return;
 }

Итак, теперь у нас вместо трех текстур на входе будет пятьЖ глубина принимающей геометрии, а так же глубина и нормали передних и задних граней преломляющей геометрии. Алгоритм расчета в шейдере дополнится нахождением пересечения лучей света сперва с задними гранями геометрии, а затем уже с принимающей геометрией.

Дисперсия света

Как известно абсолютный показатель преломления вещества зависит от длины волны света. То есть свет с большей длиной волны (например красный) преломляется меньше, чем свет с меньшей длиной волны (например фиолетовый). Из-за этого в узоре каустики можно иногда наблюдать весь спектр (радугу). Для определения показателя преломления существует формула Коши, вида:

 n = a + b/L^2

где n - показатель преломления, a и b — эмпирически подобранные константы, L — длина волны света. Так как нам необходимо обойтись одной величиной для задания показателья преломления, то я предлагаю использовать следующую формулу:

n = n0 * (1.0 + (1.0/n0 - 1.0) * L^2);

где n0 — прочитанное их текстуры значение показателя преломления. При n0 = 1.0 мы получим n = n0 = 1.0, при увеличении n0 зависимость показателя преломления n от длины волны будет линейно возрастать.

Дисперсия света при проходждении через прозрачный тор

Для того чтобы смоделировать дисперсию света, нам придется немного изменить алгоритм отрисовки. Если сейчас мы просто устанавливали константное значение цвета для точек (фотонов) и один раз отрисовывали вершинный буфер, то теперь нам нужно будет задавать цвет и длину волны и как минимум три раза выводить вершинный буфер на отрисовку (для R, G, B цветов соответственно). Раньше у нас было:

 photonBuffer->renderAll();

а теперь стало: