Данный метод является, возможно, наиболее известным способом процедурной генерации текстур на данный момент.
Он описан в широко известной статье Хьюго Элиаса (Hugo Elias) и её различныхпереводах, где он (ошибочно) называется шумом Перлина. Кроме этого статья содержит несколько недостатков:
1. Рекомендация использования косинусоидальной интерполяции.
2. Двумерная функция Noise1 использует одномерную функцию с индексом n = x + y * 57. Т. к. 57 - составное число, то это приводит к заметным визуальным артефактам.
3. Используется подход через шумовую функцию, без упоминания о том, что для заполнения текстуры (что является основным упоминаемым в статье использованием) может быть достаточно ГСЧ, что улучшает производительность. Шумовая функция позволяет, впрочем, не хранить текстуру явно.
Также использованная функция шума на решётке (Noise1) при прямолинейном переводе с Pascal на C приводит к коду, который содержит undefined behavior (signed arithmetic overflow):
+ Например
− Скрыть
double IntegerNoise (int n)
{
n = (n >>13) ^ n;
int nn = (n * (n * n * 60493 + 19990303) + 1376312589) & 0x7fffffff;
return1.0 - ((double)nn / 1073741824.0);
}
(пример взят из библиотеки libnoise)
что, впрочем, не является виной автора.
Вкратце, суть метода:
1. В узлах решётки (точках с целочисленными координатами) детерминированным образом задаются псевдослучайные числа в диапазоне [-1;+1].
2. Определяется функция от вещественных координат, как интерполяция между значениями функции из п.1, в зависимости от дробной части координат.
3. Шумовая функция определяется как сумма заданного числа октав функции из п.2 с заданным коэффициентом a.
Метод может быть использован для произвольного количества измерений.
Наиболее частыми функциями интерполяциями являются: линейная (или би-, трилинейная, и т. д. для пространств более высокой размерности), кубическая (Catmull–Rom spline), S-curve (3-й и 5-й степени).
Их формулы для одномерного случая:
+ Показать
− Скрыть
Линейная:
\(x=x_0 (1-t)+x_1 t\)
обеспечивает непрерывность функции.
Кубическая (Катмулл-Ром):
\(\begin{align*}
x=\frac{1}{2}( & (-t^3+2 t^2-t) & x_{-1} & + \\
& (3 t^3-5 t^2+2) & x_{0} & + \\
& (-3 t^3+4 t^2+t) & x_{1} & + \\
& (t^3-t^2) & x_{2} & )
\end{align*}\)
обеспечивает непрерывность функции и её 1-й производной.
Кубическая S-curve:
\(x=x_0 (1-(3 t^2-2 t^3))+x_1 (3 t^2-2 t^3)\)
обеспечивает непрерывность функции и её 1-й производной.
S-curve 5-й степени:
\(x=x_0 (1-(6 t^5-15 t^4+10 t^3))+x_1 (6 t^5-15 t^4+10 t^3)\)
обеспечивает непрерывность функции и её 1-й и 2-й производной.
Примечание:
В отличие от остальных приведённых функций интерполяции, Catmull–Rom spline:
1. Требует 4 контрольные точки, а не 2.
2. Может давать значения, выходящие за диапазон входных данных (например, последовательность [0,1,1,0] в центре даёт 1.25).
Пример реализации для двумерного случая и линейной интерполяции (диапазон значений [-1;+1]):
Пример демонстрационный и является сильно неоптимальным по производительности (т. к. оптимизировалась в нём понятность).
Замер времени его работы на машинах ideone дал следующие результаты (не стоит на них полагаться, т. к. качество экспериментальной методики вызывает вопросы):
+ Показать
− Скрыть
octaves_of_value_noise:
octaves
tacts per call
1
337
2
613
3
894
4
1171
5
1466
grid_value_noise: 32 tacts.
Для получения текстур предлагается использовать самплинг данной функции. Следует заметить, что в таком виде функция даёт текстуры, не являющиеся бесшовными, что впрочем легко исправляется (добавлением модуля в grid_value_noise).
Однако если целью является однократное получение статической текстуры с линейной интерполяцией, то возможности шумовой функции (которая позволяет, например, элементарно получать масштабированные/смещённые/повёрнутые версии одной и той же текстуры) являются избыточными, и тот же подход можно реализовать эффективнее: генерируем серию текстур, размера size/2^(octaves-1), size/2^(octaves-2), ..., size - первая заполнена белым шумом, каждая следующая - суммой белого шума и upscale-версии предыдущей текстуры с весами, соответствующими параметру a.
Это приводит к примерно следующей реализации (диапазон [0;1]):
http://ideone.com/2VP3jg
Эта реализация генерирует квадратную power-of-two бесшовную текстуру, не использует выделение дополнительной памяти (идея layout'а взята у TarasB), и её быстродействие слабо зависит от количества октав. Несмотря на не совсем линейный порядок обхода (интуитивно могущий приводить к дополнительным кеш-промахам), время работы в пересчёте на 1 тексель не очень сильно (менее, чем вдвое) зависит от размера текстуры.
Как легко видеть, код примерно в 4 раза быстрее предыдущей реализации, даже для 1 октавы.
Впрочем, обобщение этого способа на другие функции интерполяции не является очевидным.
Рассматриваемый в этом пункте шум имеет некоторые недостатки. В частности, особенно при низком количестве октав, он обладает хорошо заметной анизотропией.
Качество получившейся текстуры зависит также от выбранной функции интерполяции. В зависимости от цели использования текстуры, подходящая функция может быть разной. В частности, использование текстуры как карты высот обычно более требовательно к функции интерполяции.
Статья http://libnoise.sourceforge.net/noisegen/index.html содержит хорошее изложение данной темы. Примеры оттуда (во всех примерах используется только 1 октава, для большей наглядности артефактов):
+ Показать
− Скрыть
Интерполяция 0-го порядка.
Билинейная (т. е. 1-го порядка) интерполяция. Справа - полученная по ней карта высот.
Интерполяция S-curve 3-го порядка. Справа - полученная по ней карта высот.
Интерполяция S-curve 5-го порядка. Справа - полученная по ней карта высот.
Для сравнения текстуры градиентного шума (тоже 1 октава):
+ Показать
− Скрыть
Градиентный шум. Справа - полученная по нему карта высот.
Шум Перлина
Шум Перлина разработан Кеном Перлином в 1983 г. и является примером градиентного шума.
Алгоритм реализует шумовую функцию в произвольном количестве измерений. Исходная реализация была 3-мерной. Её код доступен (для ознакомительных целей; он всё ещё защищён копирайтом) на странице Кена Перлина:
http://mrl.nyu.edu/~perlin/doc/oscar.html#noise
Описание алгоритма и некоторое количество сопутствующей информации можно найти в слайдах Перлина "Making noise" (текст справа может быть не очень удобно читать):
http://www.noisemachine.com/talk1/index.html
Алгоритм нахождения значения функции в заданной точке следующий:
1. В узлах решётки (точках с целочисленными координатами) детерминированным образом задаются псевдослучайные n-мерные единичные (или, по крайней мере, одинаковые по длине) векторы (называемые градиентами, о причинах - см. далее).
2. Для точки находится целочисленная ячейка (единичный гиперкуб с целочисленными вершинами), которая её содержит.
3. Для каждого угла ячейки строится вектор из этого угла в точку, и находится его скалярное произведение с градиентом в этом угле.
4. Полученные скалярные произведения интерполируют (линейная интерполяция с весами, соответствующим S-curve от дробных частей координат: важно, чтобы 1-е производные от функции интерполяции были равны 0 на краях). Такая интерполяция сепарабельна. С практической точки зрения, вычисляются n значений S-curve (по 1 для каждой координаты) и 2^{n}-1 раз проводится линейная интерполяция (2^{n-1} по первой координате, затем 2^{n-2} - по 2-й для результатов предыдущего шага, и т. д.).
+ Пример для 3D
− Скрыть
Пусть имеются 8 посчитанных значений скалярного произведения в углах куба:
\(s_{000},\,s_{100},\,s_{010},\,s_{110},\,s_{001},\,s_{101},\,s_{011},\,s_{111}\)
Сначала проинтерполируем по x:
\(s_{*00}=s_{000} (1 - t_0)+s_{100} t_0\) \(s_{*10}=s_{010} (1 - t_0)+s_{110} t_0\) \(s_{*01}=s_{001} (1 - t_0)+s_{101} t_0\) \(s_{*11}=s_{011} (1 - t_0)+s_{111} t_0\)
где
\(t_0=3 x^2-2 x^3\)
Дальше - по y:
\(s_{**0}=s_{*00} (1-t_1)+s_{*10} t_1\) \(s_{**1}=s_{*01} (1-t_1)+s_{*11} t_1\)
где
\(t_1=3 y^2-2 y^3\)
И окончательно, по z:
\(s_{***}=s_{**0} (1-t_2)+s_{**1} t_2\)
где
\(t_2=3 z^2-2 z^3\)
Легко убедиться, что интерполяция сепарабельна, т. е. результат не зависит (с точностью до погрешностей округления при работе с числами с плавающей точкой) от порядка, в котором мы обрабатываем оси.
То же самое, более формально:
1. \(\vec{g}(\vec{k}): \vec{k} \in \mathbb{Z}^n, \vec{g}(\vec{k}) \in \mathbb{R}^n, \vec{g}(\vec{k}) \cdot \vec{g}(\vec{k}) = 1\) — градиенты в узлах решётки.
2. \(\vec{p}=\vec{k}+\vec{r}, \vec{k} \in \mathbb{Z}^n, \vec{r} \in [0;1]^n\) — целочисленная и дробная части координат.
3. \(s_{\vec{i}}=s_{i_0,i_1,\dots,i_{n-1}}=\left(\vec{g}(\vec{k}+\vec{i}) \right) \cdot (\vec{r}-\vec{i}), i_q \in \{0, 1\}\) — скалярные произведения в углах ячейки.
4.
\(\mathrm{Noise}(\vec{p})=
\sum_{\vec{i} \in \{0,1\}^n}
s_{\vec{i}}
\left(
\prod_{q=0}^{n-1}
1-(3 r_q^2 - 2 r_q^3)
\right)\) — интерполяция.
Легко показать, что если условие равенства 0 производной на краях выполнено, то градиент получившейся скалярной функции в целочисленных точках будет равен заданным там векторам.
Данная функция несколько медленнее интерполяции численного шума, рассмотренной в предыдущем пункте, но даёт более качественный шум. Несмотря на то, что шум Перлина не является изотропным, его анизотропия гораздо менее ярко выражена, чем у численного шума на декартовой сетке.
Демонстрация:
+ Показать
− Скрыть
Слева - шум Перлина с осями "вправо/вверх", справа - его поворот на 37°.
Public Domain реализацию 3-мерного шума Перлина (одним файлом) можно найти в составе коллекции библиотек stb:
https://github.com/nothings/stb
Она реализует несколько другую модификацию алгоритма. В stb-реализации случайный вектор - одна из 12 возможных перестановок {±1,±1,0}, выбираемой с помощью хеш-функции от (целочисленных) координат (revised Perlin noise, http://mrl.nyu.edu/~perlin/noise/ ) и процедура выбора тоже несколько отличается.
Также шум Перлина реализован в составе библиотеки libnoise.
Как легко видеть алгоритм имеет сложность O(2^n), где n - размерность пространства. Для высоких размерностей его быстродействие получается очень низким. Для борьбы с этим Кен Перлин в 2001 г. предложил симплекс-шум, который разбивает пространство на симплексы, а не на гиперкубы, что даёт возможность снизить сложность алгоритма до O(n^2). Симплекс-шум подпадает под действие патента (по крайней мере в США), но есть альтернативные открытые реализации (напр. OpenSimplex noise). Статья Перлина с описанием алгоритма доступна по ссылке:
http://www.csee.umbc.edu/~olano/s2002c36/ch02.pdf
Sparse Convolution Noise
Этот вид шума был предложен Джоном Питером Льюисом в 1989 г. в статье Algorithms for Solid Noise Synthesis.
Он представляет собой свёртку некоторым образом выбранного ядра (англ. kernel) \(h\) с Пуассоновским шумовым процессом:
\(N(\vec{r})=\left[ h * \gamma \right] (\vec{r})\)
другими словами:
\(N(\vec{r})=\int_{\mathbb{R}^n} h(\vec{r}') \gamma(\vec{r}-\vec{r}') dr'^n\)
Пуассоновский шумовой процесс - это сумма (обычно - бесконечного числа) точечных импульсов в случайно выбранных равномерно распределённых точках со случайными нескореллированными амплитудами:
\(\gamma(\vec{r})=\sum_i w_i \delta(\vec{r}-\vec{r}_i)\),
где \(\delta(\vec{x})\) - дельта-функция Дирака.
Название отражает тот факт, что количество импульсов в отдельно выбранном элементе пространства описывается распределением Пуассона.
Из этого следует, что:
\(N(\vec{r})=\sum_i w_i h(\vec{r}-\vec{r}_i)\)
Эта сумма содержит, в общем случае, бесконечное количество слагаемых.
Для того, чтобы её вычисление занимало конечное время, выбирают ядро, являющееся финитной функцией (т. е. функцией с ограниченным носителем; примечание: в \(\mathbb{R}^n\) понятия ограниченный носитель и компактный носитель тожественны), а следовательно - обращающееся в 0 на некотором расстоянии от начала координат (радиусе ядра).
Льюис использовал ядро
\(h(\vec{r})=
\left\{\begin{matrix}
\frac{1}{2}+\frac{1}{2} \cos(\pi |r|) & , |r|<1 \\
0 & , |r| \geq 1
\end{matrix}\right.\)
В статье Procedural Noise using Sparse Gabor Convolution бесконечное ядро обрезается по ~ 5% максимума. Функция при этом получается разрывной, что, по утверждениям авторов, не приводит к визуальным артефактам даже для производных. При желании разрыв можно убрать использованием оконной функции.
Для ядра радиуса \(R\) алгоритм разбивает пространство на ячейки ребром \(R\). Тогда значение шумовой функции в точке зависит только от импульсов, содержащихся в ячейке, содержащей эту точку, и её непосредственных соседях.
Импульсы в ячейке генерируются ГСЧ с зерном (начальным состоянием) являющимся хеш-функцией от координат ячейки. Хранить их явным образом не требуется.
Это приводит к примерно следующей реализации для 2-мерного случая (оптимизировалась понятность, а не эффективность):
+ Показать
− Скрыть
// Sparse convolution noise function.// Repeats itself with period wrap_around in both x and y.float sparse_convolution_noise(float x,float y)
{
using distr=std::uniform_real_distribution<float>;
float ret=0.0f;
int32_t
ix=int32_t(std::floor(x/kernel_radius)),
iy=int32_t(std::floor(y/kernel_radius));
for(int32_t cell_y=iy-1;cell_y<=iy+1;++cell_y)for(int32_t cell_x=ix-1;cell_x<=ix+1;++cell_x)
{
uint32_t cell_hash=
(uint32_t(cell_y)%wrap_around)*wrap_around+
(uint32_t(cell_x)%wrap_around);
std::mt19937 rng(cell_hash);
// The number of iterations of the following loop// has Poisson distribution.// Cf. Knuth's "Algorithm Q".// (Donald E. Knuth, The Art of Computer Programming, Volume 2,// "Seminumerical Algorithms", Addison-Wesley (1969),// page 117.float
p=1.0f,
L=std::exp(-lambda);
while(true)
{
p*=distr(0.0f,1.0f)(rng);
if(p<=L)break;
float
impulse_x=(float(cell_x)+distr(0.0f,1.0f)(rng))*kernel_radius,
impulse_y=(float(cell_y)+distr(0.0f,1.0f)(rng))*kernel_radius,
impulse_w=distr(-1.0f,+1.0f)(rng);
ret+=impulse_w*kernel(x-impulse_x,y-impulse_y);
}
}
return ret;
}
Результат работы для lambda=16 с Льюисовским ядром:
+ Показать
− Скрыть
Алгоритм производит в среднем \(3^{n} \lambda\) обращений к функции kernel, где n - размерность пространства.
Ещё одну версию кода можно найти в конце статьи Procedural Noise using Sparse Gabor Convolution .
Как можно догадаться из изложенного выше, выбор ядра сильно влияет на характер шума.
Популярным является Габоровское ядро:
\(g(\vec{r})=e^{-\pi a^2 r^2} \cos(2 \pi F_0 \vec{k} \cdot \vec{r}+\varphi_0)\)
Где \(a\) - "острота" (обратная ширина) гауссовской огибающей, \(F_0\), \(\vec{k}\) и \(\varphi_0\) - произвольным образом выбранные частота, единичный волновой вектор и начальная фаза гармонической осцилляции.
Эта функция была рассмотрена Денешем Габором в 1946 г. в статье Theory of communication.
Её преимущество состоит в том, что она (как показал Габор) имеет минимально возможное произведение неопределённостей (в терминах квантовой механики) в координатном и частотном пространстве:
\(\Delta x \cdot \Delta f=\frac{1}{2}\)
Что позволяет получить одновременно малый размер в координатном пространстве (ради малого радиуса ядра и быстрого вычисления) и малый размер в частотном пространстве (что позволяет хорошо контролировать спектр).
В статье Procedural Noise using Sparse Gabor Convolution для каждого импульса используется своё Габоровское ядро, т. е. каждому импульсу, кроме координат и веса сопоставлены параметры \(a\), \(F_0\) и \(\vec{k}\) связанного с ним ядра.
