Участник:Kst179/Metropolis light transport

Metropolis light transport — это алгоритм, применяющий один из вариантов метода Монте Карло — метод Метрополиса-Гастинга, Монте Карло на Марковских цепях (Metropolis-Hastings, Markov chain Monte Carlo, MCMC), для генерации изображений при помощи физического описания трехмерных сцен. Алгоритм используется в компьютерной графике, решает задачу глобального освещения (точное моделирование световых эффектов). Отличается от классических методов, таких как Bidirectional path tracing скоростью работы и точностью моделирования (возвращает несмещенное значение).

1 Свойства и структура алгоритма

1.1 Общее описание алгоритма

Основной объект с которыми работает трассировщик лучей это пути от источников света до камеры. Путь это ломанная, каждое звено которой имитирует направление движения фотона, ломанная начинается в источнике света, изменяет направление при попадании на какую либо поверхность, а последнее звено должно заканчиваться в камере. Обозначим путь как $\overline{x} = x_1, x_2, \dots, x_k$, где $k\gt 1, x_i$ — координата на поверхности некоторого объекта от которого отражается луч). Для каждого пути введем некоторую функцию $f(\overline{x})$, обозначающую яркость, с которой этот луч освещает пиксель матрицы камеры.

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

В MLT путь считается случайной величиной, с плотностью вероятности $p(\overline{x}) = \frac{1}{b}f(\overline{x}), b=const$, а генерация путей происходит согласно этому распределению. Тогда получается что более весомых путей мы генерируем больше, а пути проходящие через стены или на которых теряется слишком много яркости мы практически не рассматриваем. Также если генерировать пути согласно вышеописанному распределению то выражение для подсчета окончательного изображения $m_j = \int _{\Omega} w_j(\overline{x})f(\overline{x})d\overline{x}$, (где $m_j$ — значение $j$-го пикселя, $\Omega$ — пространство всевозможных путей, $w_j(\overline{x})$ — индикатор того, что $\overline{x}$ попадает в $j$-й пиксель) можно заменить на следующую несмещенную оценку:

$m_j = \int _{\Omega} w_j(\overline{x})f(\overline{x})d\overline{x} = \int _{\Omega} w_j(\overline{x})\frac{f(\overline{x})}{p(\overline{x})} p(\overline{x})d\overline{x} = \mathbb{E} w_j(\overline{x})\frac{f(\overline{x})}{p(\overline{x})} \approx \frac{1}{N} \sum _{i=0} ^{N} w_j(\overline{x_i})\frac{f(\overline{x_i})}{p(\overline{x_i})} = \frac{b}{N} \sum _{i=0} ^{N} w_j(\overline{x_i})$.

То есть после семплирования путей $\overline{x_1}, \overline{x_2}, \dots, \overline{x_N}$ получаем что яркость j-го пикселя пропорциональна количеству приходящих в этот пиксель лучей а от их вида никак не зависит (но для этого необходимо чтобы пути строго соответствовали своему распределению).

1.2 Математическое описание алгоритма

Для семплирования путей из, вообще говоря, неизвестного нам распределения $p(\overline{x})$ используется метод MCMC. Построим следующую марковскую цепь: введем некоторую вероятность $K(\overline{y}|\overline{x})$ перехода $\overline{x} \rightarrow \overline{y}$. Выберем произвольный путь $\overline{x_0}$ и случайным блужданием будем переходить от $\overline{x_{i-1}}$ к $\overline{x_i}$. Пусть $\overline{x_0}$ распределено с некоторой плотностью вероятности $p_0(\overline{x})$, тогда $x_i$ будет распределено как $p_i(\overline{x}) = \int_{\Omega} K(\overline{x}|\overline{y})p_{i-1}(\overline{y})d\overline{y}$. При некоторых условиях на $K(\overline{x}|\overline{y}), ~\underset{i\rightarrow\infty}{lim}p_i(\overline{x}) = p(\overline{x})$, откуда следует что при правильном выборе $K(\overline{x}|\overline{y})$, и большом количестве переходов, независимо от выбора $\overline{x_0}$, все $\overline{x_i}$ будут распределены согласно вероятности $p(\overline{x})$.

Пусть мы находимся в состоянии (пути) $\overline{x}$, сделаем над путем некоторую мутацию (небольшое изменение) и получим некоторый другой путь $\overline{y}$. На пространстве всевозможных мутаций над $\overline{x}$ введем вероятность $T(\overline{y}| \overline{x})$ того что при мутации $\overline{x}$ получится $\overline{y}$. Далее мы можем принять эту мутацию или отклонить с некоторой вероятностью $a(\overline{y}|\overline{x})$. Если мутация принимается, то совершается переход по марковской цепи $\overline{x} \rightarrow \overline{y}$, иначе $\overline{x}\rightarrow\overline{x}$, т.е. мы остаемся в том же состоянии. Таким образом мы сконструировали $K(\overline{y}|\overline{x}) = a(\overline{y}|\overline{x})T(\overline{y}|\overline{x})$, при $\overline{x} \neq \overline{y}, K(\overline{x}|\overline{x}) = \int_{\Omega}(1-a(\overline{y}|\overline{x}))T(\overline{y}|\overline{x})d\overline{y}$, где $T(\overline{y}|\overline{x})$ определяется тем, как именно мы производим мутации над путями, будет обсуждено чуть позже, а $a(\overline{y}|\overline{x})$ выбирается таким, чтобы цепь сходилась к стационарному состоянию, и скорость сходимости была как можно выше. Согласно принципу детального равновесия Марковская цепь сходится к стационарному распределению $p(\overline{x})$ при $p(\overline{x})K(\overline{y}|\overline{x}) = p(\overline{y})K(\overline{x}|\overline{y})$, или, что тоже самое: $\frac{1}{b}f(\overline{x})a(\overline{y}|\overline{x})T(\overline{y}|\overline{x}) = \frac{1}{b}f(\overline{y})a(\overline{y}|\overline{x})T(\overline{x}|\overline{y})$.

Откуда получаем $\frac{a(\overline{y}|\overline{x})}{a(\overline{x}|\overline{y})} = \frac{f(\overline{y})T(\overline{x}|\overline{y})}{f(\overline{x})T(\overline{y}|\overline{x})}$. Так как скорость сходимости необходима как можно большая, то выберем $a(\overline{y}|\overline{x})$ и $a(\overline{x}|\overline{y})$ максимально возможным: $a(\overline{y}|\overline{x}) = min\left\{ 1, \frac{f(\overline{y})T(\overline{x}|\overline{y})}{f(\overline{x})T(\overline{y}|\overline{x})} \right\}$.

Теперь определим функцию $f(\overline{x})$. Функция $f(\overline{x})$ складывается из нескольких промежуточных значений, каждое из которых отображает какая часть освещенности падает при прохождении звена пути/отражения луча от поверхности. Рассмотрим одно звено пути $x' \rightarrow x''$. Пусть x' находится на поверхности объекта излучающего свет. Тогда определим $L_e(x'\rightarrow x'')$ как освещенность точки $x''$ при излучении света из точки $x'$. Для точечного источника $L_e(x'\rightarrow x'') = \frac{I}{||x' - x''||^2}$, где I — интенсивность источника, для излучающей поверхности с нормалью $\overline{n}$ в точке $x'$, $L_e(x' \rightarrow x'') = I\frac{|\lt \overline{n}, x''-x'\gt |}{||x'-x''||^3}$ и т.д. для каждого типа источника функция будет разной.