Методы ускорения расчета освещенности с использованием квази- Монте Карло интегрирования

( Speed up methods for global illumination problem using quasi- Monte Carlo integration
Preprint, Inst. Appl. Math., the Russian Academy of Science)

Содержание

Содержание......................................................................................................................... 3

1. Введение.......................................................................................................................... 4

2. Метод распространения................................................................................................. 6

3. Применение квази- Монте Карло для метода распространения............................... 8

4. Когерентная трассировка фотонов.............................................................................. 11

5. Снижение размерности интегрирования квази- Монте Карло............................... 13

6. Трассировка вторичных фотонов................................................................................ 17

7. Сравнение псевдо- и квази- Монте Карло методов................................................. 18

8. Заключение.................................................................................................................... 20

Список литературы........................................................................................................... 21

1.    Введение

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

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

,               (1)

Рис.1.

где  - излучательность точки поверхности от прямого освещения, а интегральный оператор , с учетом функции видимости точек сцены , описывает взаимодействие света пришедшего от других освещенных поверхностей (Рис.1):

              (2)

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

                                  (3)

Рис.2.

Наряду с энергетическими формами уравнения рендеринга (1, 2) рассматривают его сопряженный вид, описывающий функцию потенциальной видимости  из положения наблюдения  в направлении  (Рис.2):

                  (4)

где  - прямая видимость в направлении наблюдения, а интегральный оператор  сопряжен с оператором :

      (5)

Принцип взаимности Гельмгольца (при условии его применимости) для подынтегральной функции сопряженного оператора:

допускает использование канонической формы функции двунаправленного отражения как для уравнения рендеринга, так и для уравнения видимости:

      (6)

Методы построения изображений на основе решения уравнения рендеринга для излучательности (2) носят название прямых методов, в то время как обратные методы для тех же целей используют уравнение потенциальной видимости (6). Отдельный класс занимают двунаправленные методы, комбинирующие прямые и обратные методы построения изображений. Общая классификация способов решения интегральных уравнений для любого из методов приведена в [2] и сводится к трем типам:

Обращение

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

Итерирование

Сходимость итерирования обусловлена поглощением энергии происходящей на каждом этапе взаимодействия света с поверхностью. Для представления аппроксимирующих функций , как правило, используются методы конечных элементов. Примером решений, основанных на итерировании, являются методы диффузной излучательности (diffuse radiosity) [3], разделенной полусферы (partitioned hemisphere) [4] и направленных распределений (directional distributions) [5]. Существенный недостаток итерационных методов проявляется на примере сцен со слабо поглощающими поверхностями [6], когда аппроксимация оператора   требует значительных вычислительных ресурсов.

Распространение

Решения, основанные на методе распространения, в отличие от предыдущего метода, не требуют промежуточного представления для искомой функции излучательности L или функции потенциальной видимости W.  В целом, метод распространения обеспечивает более гибкий и универсальный механизм взаимодействия света с точкой поверхности, анализу преимуществ которого посвящена работа [7].

2.    Метод распространения

Прямая трассировка лучей, использующая метод распространения для решения уравнения излучательности, носит название трассировки фотонов. Распространение фотонов в качестве бесконечного ряда  имеет простой интуитивный смысл (Рис.3):  вносит вклад прямого освещения, - вклад света, пришедшего из источника после первого отражения,  - вклад после двух переотражений и т.д. Нетрудно заметить,

Рис.3.

что член  имеет размерность 2d+2, что естественно сразу же приводит к высоким размерностям интегрирования.

Рис.4.

Обратную трассировку лучей, дающую решение уравнения потенциальной видимости W, называют трассировкой путей. Каждый трассируемый путь берет начало из камеры наблюдателя и вносит вклад в соответствующий пиксель экрана (Рис.4). Как сопряженный метод, трассировка путей, очевидно, имеет ту же размерность, что и трассировка фотонов.

Таб.1.

При реализации метода распространения, особой популярностью пользуется метод Монте Карло, характерным свойством которого является независимость от размерности интегрирования. Погрешность метода представляет собой случайный шум, не зависящий от количества переотражений и преломлений, претерпеваемых трассируемыми фотонами (путями), а погрешность оценивается как O(N^-1/2), где N – количество трассируемых лучей, что специфично для любого метода Монте Карло [8]. Таблица (Таб.1) содержит изображения, полученные методом Монте Карло для первых n членов ряда фотонов при n = 0, n = 1, n = 2, n = 4.

3.    Применение квази- Монте Карло для метода распространения

При интегрировании с помощью метода Монте Карло, необходимо придумать некоторую случайную величину h, математическое ожидание которой равно значению искомого интеграла. Тогда это значение может быть оценено с любой наперед заданной точностью с помощью статистического среднего величины h  по N  испытаниям, где N достаточно велико:

 ,                               (7)

где h_i – значение случайной величины h  из области интегрирования W в i-м испытании, p_h(x) - плотность распределения h  в W .

Известно, что случайность точек h_i не является необходимым условием сходимости суммы в правой части (7) к искомому интегралу. Для этого достаточно, чтобы эти точки были распределены с плотностью . Построен целый ряд неслучайных последовательностей точек, позволяющий производить интегрирование в соответствии с (7). Такие последовательности получили название квазислучайных. Обычно строятся равномерно распределенные последовательности точек, а последовательности с требуемым распределением  получают из них путем некоторой трансформации. Примером квазислучайных последовательностей точек могут служить ЛП_τ-последовательности, построенные в работе [9]. Обычно демонстрация большей равномерности квазислучайных последовательностей точек по сравнению с псевдослучайной производится путем проектирования точек на единичный квадрат. Квадрат внизу равномерно разбит на 64 подквадрата, и на него нанесены 64 псевдослучайных точки (Рис.5) и 64 квазислучайных точки ЛП_τ-последовательности (Рис.6). Видно, что в каждый подквадрат попало ровно по одной квазислучайной точке, в то время как для псевдослучайных точек равномерное заполнение  подквадратов не выполняется:

Известно, что алгоритмы квази- Монте Карло, использующие квазислучайные последовательности точек вместо псевдослучайных чисел, сходятся к искомому результату, по крайней мере не медленнее, а обычно быстрее, чем соответствующие алгоритмы Монте Карло.

Теоретическая скорость сходимости интегрирования Монте Карло с использованием квазислучайных последовательностей оценивается неравенством Коксма-Хлавка, которое утверждает, что для интегрирования Монте Карло функции, определенной в единичном d-мерном кубе R^d с использованием N точек, ошибка интеграции E_N(f) ограничена сверху:

                               (8)

где D_N - невязка используемой последовательности, а V(f) – вариация функции в смысле Харди-Крауза. В отличие от интегрирования Монте Карло с использованием псевдослучайной последовательности, скорость сходимости интегрирования с использованием квазислучайной последовательности зависит от кратности интеграла. Это объясняется увеличением величины невязки квазислучайных последовательностей с ростом числа измерений. При достаточно большом количестве лучей, ошибка Монте Карло интегрирования с d-мерной квазислучайной последовательностью оценивается как . Для заранее известного числа N  трассируемых лучей в качестве квази-последовательности может использоваться набор точек Хаммерсли с чуть более улучшенной оценкой ошибки: .

Для небольших значений N невязка интегрирования с использованием квази-последовательности больше и практически равна ошибке псевдослучайной последовательности, т.е. O(N^-1/2).

Для увеличения производительности алгоритма Монте Карло трассировки лучей, обычно, достаточно просто формально заменить датчик псевдослучайных чисел на генератор квазислучайных точек. Однако, размерность метода распространения, под которой понимается максимальная размерность точек h_i, необходимая для трассировки одного фотона или обратной трассировки пути, способна нивелировать все преимущества квазислучайности. В частности, в работе [10] показано, что при размерностях интегрирования, больших 12, сходимость квазислучайной последовательности практически равна псевдослучайной последовательности.

Заметное влияние на результаты квази- Монте Карло интегрирования оказывает и характер интегрируемой функции. Очевидно, что если функция d-мерного аргумента x существенно зависит только от координат с младшими индексами, тогда, по сути, интегрирование производится не в d-мерном кубе, а в кубе меньшей размерности. Функция может существенно зависеть и от всех координат, но зависимость от старших координат может быть существенно меньшей, чем от младших. Именно так происходит при решении уравнения рендеринга. Свет поглощается в процессе переотражений, поэтому члены ряда фотонов с большими степенями имеют меньшее влияние на значение интегральной суммы. Чем меньше эта зависимость, тем больший выигрыш можно надеяться получить от замены псевдослучайных чисел квазислучайными точками. Отсутствие гладкости интегрируемой функции, что вполне типично для уравнения рендеринга, напротив, ухудшает сходимость квази- Монте Карло и приближает ее к сходимости Монте Карло.

4.    Когерентная трассировка фотонов

Рассмотрим квазислучайную последовательность, получаемую обращением порядкого номера числа последовательности. Функция обращения может быть представлена в виде:

          ,                            (9)

где  – базис последовательности,  - j-цифра числа i, такая что .

Значения  для некоторых значений i  приведены в таблице 2. Как ясно видно из таблицы, чтобы получить значение функции обращения , необходимо написать i справа налево после знака десятичной дроби. Аналогичная процедура имеет место для произвольного базиса последовательности b, при условии, что i и  выражены для того же базиса. Погрешность интегрирования с использованием последовательности обращения  при любом базисе оценивается как .

Для трассировки лучей необходимо ассоциировать индекс каждого трассируемого луча с последовательностью квазислучайных чисел, рассматривая, таким образом, многомерную последовательность. Первое квазислучайное число рассматриваемой последовательности может, например, использоваться для вероятностного выбора точечного источника света из общего набора источников, второе и третье - для генерации направления луча, четвертое – для выбора события в точке пересечения луча с поверхностью, следующие два – для генерации нового направления луча в соответствии с выбранным событием и т.д. Многомерная последовательность Халтона строится с использованием простых чисел в качестве базисов последовательности обращения:

         (10)

Последовательность обращения обладает свойством периодичности. В самом деле, из определения последовательности следует, что

                             (11)

(n, k – неотрицательные целые)

Свойство периодичности иллюстрируется таблицей для , из которой видно, что  и  имеют близкие значения (i = 123). Свойство периодичности справедливо и для многомерной последовательности Халтона. Заметим, что если  период для , а  для , то  является периодом для двумерной последовательности {,}. Аналогично, для общего случая многомерной последовательности Халтона, … является периодом многомерной последовательности {,, ... }. Свойство периодичности многомерной последовательности Халтона было, в частности, использовано в [11] для создания высокоэффективной трассировки пространственно-когерентных фотонов в динамически изменяющихся сценах. Набор квазислучайных фотонов на этапе генерации разбивается на пространственно-когерентные группы, каждая из которых занимает свою область в многомерном пространстве фотонов. На практике это означает, что когерентные фотоны из одной группы, как правило, излучаются одним и тем же источником света в близких к друг другу направлениях, часть из них имеют близко расположенные друг к другу точки пересечения с некоторой поверхностью сцены, взаимодействие с которыми дает похожий результат. Подобная когерентность фотонов, сама по себе, демонстрирует значительное преимущество при кэшировании данных в ходе трассировки лучей.

5.    Снижение размерности интегрирования квази- Монте Карло

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

В процессе трассировки ряд событий (например, выбор источника света) относится к выбору одного исхода из ряда возможных и имеет дискретную плотность распределения. Известно, что при выборе дискретных событий можно использовать одно и то же псевдослучайное значение несколько раз, применяя каждый раз процедуру перенормировки. Этот метод, первоначально предложенный для экономии псевдослучайных чисел в целях уменьшения времени работы генератора псевдослучайных чисел [12], носит название модифицированного метода суперпозиции. Формально, если дискретное событие зависит от того, попадает или нет (квази-) случайное число s в полуинтервал [a, b), то перенормировка означает, что следующее случайное число s’ может быть получено без увеличения конструктивной размерности алгоритма:

s’ = (s – a) / (b – a).

Традиционно, в алгоритме Монте Карло, для моделирования источника света с заданным фотометрическим распределением, или рассеивающих свойств поверхности, заданными функцией двунаправленного отражения, применяется метод отбора-отказа, предложенный фон Нейманом. В основе этого метода лежит простое геометрическое соображение, применяемое при вычислении определенных интегралов методом проб и ошибок [8].

Рис.7.

Фотометрическое распределение источника, равно как и функция двунаправленного отражения, задается на сетке в сферической системе координат (Рис.7). Для определения значения функции по направлению (φ, θ), находится ячейка сетки, которой принадлежит это направление. Искомое значение вычисляется билинейной интерполяцией значений A₁, A₂, B₁, B₂ в углах ячейки. Например, сила света внутри изображенной на рисунке ячейки фотометрической сетки определяется по формуле:

Идея метода отбора-отказа состоит в том, что в соответствии с равномерным распределением генерируются пары углов (φ, θ) внутри ячейки до тех пор, пока очередная пара не удовлетворит условию, при котором выбор (φ, θ) пропорционален плотности фотонов I(φ, θ). Применение метода отбора-отказа в алгоритмах квази-Монте Карло ведет к потере производительности по двум причинам: возрастает конструктивная размерность алгоритма, и увеличивается «разрывность» задачи [12]. Поэтому, при переходе на квазислучайные последовательности точек, необходимо заменить его методом обратных функций.

Пусть выбрана ячейка C_{k ,} вероятность выбора которой пропорциональна потоку через нее. Направление излучения фотона внутри C_k обладает плотностью вероятности, пропорциональной I(φ, θ):

Условные плотности вероятности требуемых углов:

.

Соответствующие p₁ и p₂ функции распределения будут равны:

.

Тогда углы (φ, θ) излучения фотона определяются из уравнений:

,   ,

откуда

Необходимое распределение находится численным решением трансцендентного уравнения для η₁, например методом секущей.

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

Заметим, что изложенный алгоритм имеет смысл только в случае использования квази- Монте Карло интегрирования, т.к. улучшает характер сходимости. В случае алгоритмов Монте Карло, сходимость которых не зависит от конструктивной размерности, замена метода отказов на метод обратных функций, в данном случае, существенного выигрыша не дает.

Значительно снизить размерность интегрирования d и улучшить оценку ошибки квазислучайной последовательности позволяет комбинирование прямой и обратной трассировок лучей в виде двунаправленных методов распространения, нашедших реализацию в современных алгоритмах моделирования распространения света [15, 16, 17].

Двунаправленные методы распространения сочетают преимущества прямой и обратной трассировок лучей, значительно повышая эффективность вычислений. В таких методах, перенос световой энергии между трассируемыми фотонами и обратными путями осуществляется посредством карты освещенности, впервые предложенной в работе [18]. Карта освещенности регистрирует энергию фотонов на сетке, покрывающей диффузные грани геометрических примитивов в сцене. В любой момент времени трассировки фотонов, освещенность в j-элементе сетки карты вычисляется как:

,

где m_j – количество фотонов зарегистрированных на j–элементе, Ф – энергия одного фотона, S_j – площадь j–элемента сетки. При обратной трассировке путей, попадание в элемент карты освещенности формирует значение искомой светимости на соответствующем пикселе экрана без дальнейшего распространения. Для квази- Монте Карло интегрирования важно, что метод распространения для трассировки обратных путей прерывается на низких размерностях, по мере попадания трассируемых путей в элементы карты освещенности. В свою очередь, элементы карты освещенности, вычисленные трассировкой фотонов с высокими размерностями и имеющими худшую оценку точности квази- Монет Карло интегрирования, могут быть проигнорированы при обратной трассировке, тем самым обеспечивая лучшую сбалансированность конструктивных размерностей для прямого и обратного методов распространения. При идеальном балансе, размерность интегрирования составляет d/2, а следовательно, ошибка интегрирования с использованием квазислучайной последовательности .

6.    Трассировка вторичных фотонов

Скорость сходимости Монте Карло интегрирования с использованием псевдопоследовательностей, а тем более квазипоследовательностей случайных чисел, зависит от гладкости интегрируемой функции. Наибольший вклад в нарушение гладкости вносит компонента прямого освещения (резкие тени, высокий градиент освещенности), в более узком смысле – прямое освещение точечными источниками света. Вычисление прямого освещения от точечных источников не представляет сложности, поэтому, логичным шагом представляется исключение точечных источников из стохастического алгоритма вообще. В работе [19] предлагается метод замены точечных источников (“first-shot” метод).

Представим искомую излучательность в виде:

L,

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

.

Заменяя светимость точечных источников на светимость элементов сетки поверхностей, ими освещенных (), получим новый набор источников света: . Таким образом, численное интегрирование применяется для модифицированной задачи вычисления  с более гладкой интегрируемой функцией. Для диффузных поверхностей хранение  требует одного числа на элемент сетки. На практике, однако, существуют поверхности со свойствами, отличными от диффузных. Для них  уже является функцией, зависящей от направления. В работе [19] предлагается хранить освещенность, полученную каждым элементом сетки от каждого точечного источника отдельно, тогда искомая  будет вычислена на основе хранимых данных «по требованию» для заданного направления.

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

7.    Сравнение псевдо- и квази- Монте Карло методов

Сравнение двух методов Монте Карло интегрирования проведено для сцен с аналитически вычисленными значениями полной излучательности в тестовых точках.

Сцена: октант сферы (Рис.8). Вторичная освещенность может быть вычислена точно при условии равенства фактора формы для всех пар точек. Сцена представляет собой 1/8 часть диффузной единичной сферы отсеченной тремя зеркальными координатными плоскостями. Обобщенный фактор формы характеризующий перенос световой энергии между диффузными и диффузно-зеркально-диффузными поверхностями одинаков для всех пар точек рассматриваемой части сферы. Теоретическое значение вторичной освещенности в центре октанта сферы в точке с координатами (3^-1/2, 3^-1/2, 3^-1/2) составляет 666.649 люкс, а полная освещенность равна 1353.247 люкс [21]. В таблице 3 приведены результаты моделирования освещенности трассировкой фотонов:

Таб.3.

Сцена: диффузный куб (Рис.9). Сцена представляет собой куб 10м´10м´10м, в центр которого помещен точечный источник света с силой света 50000 кд создающий прямое освещение 2000 люкс на ближайших точках внутренних поверхностей граней куба. Внутри куба свет распространяется диффузно, без выхода наружу. Отражательная способность диффузных граней куба 2/3. Таким образом,  большая компонента освещенности образуется за счет переотраженного света. Цвет источника и поверхностей куба белый.

Теоретические значения светимости в 6 точках на грани куба приведены в таблице 4 [21]. Исходя из симметрии сцены, для нахождения отклонения между моделированными и теоретическими значениями мы используем сетку 5´5 для каждой внутренней поверхности граней куба. Таким образом, отклонение (ошибка) теста вычисляется в L₂ норме по 150 известным точкам.

В таблице 5 приведены результаты моделирования освещенности трассировкой фотонов (квази-, псевдо-) и двунаправленной трассировкой лучей (квази²-, псевдо²-):

Таб.5.

Для сцены куб конструктивная размерность трассировки фотонов составила 6, а обратной трассировки лучей - 2.

6. Заключение

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

Практические результаты по сравнению квази- и псевдо- интегрирования для задач освещения на несложных тестовых сценах были проведены Келлером в работах [13, 14]. Основное преимущество использования квазислучайной последовательности было продемонстрировано трассировкой фотонов в пределах размерности d = 4, что соответствует выбору стартовой позиции на поверхности источника света и начальному направлению фотона. В данной работе приведено сравнение квази- и псевдо- интегрирования уравнения рендеринга для более высоких размерностей как с учетом диффузного, так и диффузно-зеркального переноса световой энергии. Для октанта сферы с зеркальными и диффузными поверхностями значительное преимущество использования квазислучайных последовательностей было получено на 3-ей минуте трассировки фотонов. Для диффузного куба на всех этапах вычислений многократное преимущество использования квазислучайных последовательностей было достигнуто с помощью двунаправленной трассировкой лучей.

Список литературы

[1]  J.T.Kajiya. The rendering equation. In Computer Graphics (SIGGRAPH'86 Proceedings), pages 143-150, 1986.

[2]  L.Szirmay-Karlos. Monte-carlo methods in global illumination. Institute of Computer Graphics, Vienna University of Technology, 1999.

[3]  F.Sillion, C.Puech. Radiosity and Global Illumination. Morgan Kaufmann Publishers, Inc., San Francisco, 1994.

[4]  D.S.Immel, M.F.Cohen, D.P.Greenberg. A radiosity method for non-diffuse environments. Computer Graphics (SIGGRAPH’91 Proceedings), 25(4):175-186, 1991.

[5]  F.X.Sillion, J.R.Arvo, S.H.Westin, D.P.Greenberg. A global illumination solution for general reflectance distributions. Computer Graphics (SIGGRAPH’91 Proceedings), 25(4):187-198, 1991.

[6]  L.Szirmay-Karlos, T.Foris, L.Neumann, B.Csebfalvi. An analysis to quasi- Monte Carlo integration applied to the transillumination radiosity method. Computer Graphics Forum (Eurographics’97), 16(3):271-281,1997.

[7]  A.Khodulev. Comparison of two Methods of Global Illumination Analysis. Technical report of Keldysh Inst., 1996.

[8]  Х.Гулд, Я.Тобочник. Компьютерное моделирование в физике. Часть 2. М.: Мир, стр. 39-42, 1990.

[9]  И.М.Соболь. Многомерные квадратурные формулы и функции Хаара. М.: Наука, 1969.

[10]B.Fox, P.Bratley, H.Niederreiter. Implementation and test of low discrepancy sequences, 1992.

[11]K.A.Dmitriev, S.Brabec, K.Myszkowski, H.P.Seidel. Interactive Global Illumination Using Selective Photon Tracing. 13-th Eutographics Workshop on Rendering, 2002.

[12]И.М.Соболь. Численные методы Монте-Карло. М.: Наука, 1973.

[13]A.Keller. Quasi-Monte Carlo radiosity. Eurographics Rendering Workshop 1996 Proceedings. Also in Rendering Techniques’96, Springer-Verlag, New York, 1996.

[14]A.Keller. Instant radiosity. SIGGRAPH’97 Proceedings, Addison-Wesley, pp. 49-56, 1997.

[15]E.P.Lafortune, Y.D.Willems. Bi-directional path tracing. Computer Graphics Proceedings, Alvor, Portugal, pp.145-153, 1993.

[16]E.Veach, L.J.Guibas. Optimally combining sampling techniques for Monte Carlo rendering. SIGGRAPH 95 Proceedings. Addison-Wesley, pp.419-428, 1995.

[17]S.N.Pattanaik, S.P.Mudur. Adjoint equations and random walks for illumination computation, ACM Transactions on Graphics 14: 77-102, 1995.

[18]P.S.Heckbert. Adaptive radiosity textures for bidirectional ray tracing. Computer Graphics (SIGGRAPH'90 Proceedings), volume 24, pages 145-154, August 1990.

[19]F.Castro, R.Martinez, M.Sbert. Quasi Monte-Carlo and extended first-shot improvements to the multi-path method. In Spring Conference on Computer Graphics’98, pages 91-102, 1998.

[20]Э.А.Копылов. Эффективный метод расчета освещенности с использованием функции распределения вероятности для генерации лучей, Графикон’2002, Труды конференции, стр. 225-229.

[21]E.A.Kopylov, A.B.Khodulev, V.L.Volevich. The comparison of Illumination Maps Technique in Computer Graphics Software, Graphicon’98, pp.146-152, 1998.