Пишем настоящий шум Перлина

от автора

По поисковому запросу шум перлина сразу попадается этот перевод на Хабре. Как справедливо заметили в комментариях к публикации, речь идёт вовсе не о шуме Перлина. Возможно, автор перевода и сам был не в курсе.

Чем выгодно отличается шум Перлина, легко заметить, если сравнить картинки.

Обычный шум (из той самой статьи):
image

Шум Перлина:
image

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

Рассмотрим двухмерный вариант. Пока напишем только класс-заготовку. Входящие данные — двухмерный вектор или два числа с плавающей точкой: x, y.

Возвращаемое значение — число от -1.0 до 1.0:

public class Perlin2d {   public float Noise(float x, float y)   {     throw new NotImplementedException();   } } 

Пару слов об интерполяции

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

Здесь в третьем условном квадрате точка в самом центре после интерполяции будет иметь значение 3:

image

Рассмотрим детальнее как там получается тройка. Координаты точки:
x:2.5, y:0.5

Целочисленные координаты точки (верхний левый угол квадрата):
x:2, y:0
получаются округлением в сторону 0 (функция floor).

Локальные координаты точки внутри квадрата получаются вычитанием:
x = 2.5 – 2 = 0.5, y = 0.5 – 0 = 0.5

Берём значение левого верхнего угла квадрата (1) и верхнего правого (2). Интерполируем верхнюю грань используя локальную координату x (0.5). Линейная интерполяция выглядит так:

static float Lerp(float a, float b, float t) { // return a * t + b * (1 - t); можно переписать с одним умножением (раскрыть скобки, взять в другие скобки):   return a + (b - a) * t; } 

Берём значение левого нижнего угла квадрата (2) и нижнего правого (7). Интерполируем нижнюю грань используя всё ту же локальную координату x (0.5).
Результаты:
верхняя: 1.5 нижняя: 4.5

Теперь осталась интерполяция верхней и нижней с использованием локальной координаты y (тоже 0.5):
1.5 * 0.5 + 4.5 * (1 – 0.5) = 3

Билинейная интерполяция самая простая но и результат не самый привлекательный.

Другие варианты интерполяции подразумевают модифицирование локальной координаты (параметра t) перед интерполяцией. Получаются более плавные переходы возле граничных значений (0 и 1).

image

В шуме Перлина задействован первый вариант, он даёт достаточно сильное искривление.

static float QunticCurve(float t) {   return t * t * t * (t * (t * 6 - 15) + 10); } ... // комбинирование с функцией линейной интерполяции: Lerp(a, b, QuinticCurve(t)) 

Главная идея и отличие шума Перлина

Всё очень просто:
1. В узлах сетки — псевдослучайные вектора (двухмерные для двухмерного шума, трехмерные для трехмерного и так далее), а не псевдослучайные числа.
2. Интерполируем между скалярными произведениями a) векторов от вершин квадрата до точки внутри квадрата (куба в трехмерном варианте) и b) псевдослучайных векторов (при описании шума Перлина их называют градиентными векторами).

В своём улучшенном варианте шума Кен Перлин использует всего 12 градиентных векторов. Для двухмерного варианта требуется всего 4 — по количеству граней (у квадрата их 4). Вектора направлены (условно из центра куба/квадрата) в сторону каждой из граней и не нормализованы.

Вот они:

{  1, 0 } { -1, 0 } {  0, 1 } {  0,-1 } 

image

Итак, каждому узлу сетки соответствует один из четырёх векторов. Пусть вектор у нас будет массивом float-ов.

    float[] GetPseudoRandomGradientVector(int x, int y)     {         int v = // псевдо-случайное число от 0 до 3 которое всегда неизменно при данных x и y          switch (v)         {             case 0:  return new float[]{  1, 0 };             case 1:  return new float[]{ -1, 0 };             case 2:  return new float[]{  0, 1 };             default: return new float[]{  0,-1 };         }     } 

Реализация

Нам понадобится скалярное произведение векторов:

    static float Dot(float[] a, float[] b)     {         return a[0] * b[0] + a[1] * b[1];     } 

Главный метод:

    public float Noise(float fx, float fy)     {         // сразу находим координаты левой верхней вершины квадрата         int left = (int)System.Math.Floor(fx);         int top  = (int)System.Math.Floor(fy);          // а теперь локальные координаты точки внутри квадрата         float pointInQuadX = fx - left;         float pointInQuadY = fy - top;          // извлекаем градиентные векторы для всех вершин квадрата:         float[] topLeftGradient     = GetPseudoRandomGradientVector(left,   top  );         float[] topRightGradient    = GetPseudoRandomGradientVector(left+1, top  );         float[] bottomLeftGradient  = GetPseudoRandomGradientVector(left,   top+1);         float[] bottomRightGradient = GetPseudoRandomGradientVector(left+1, top+1);          // вектора от вершин квадрата до точки внутри квадрата:         float[] distanceToTopLeft     = new float[]{ pointInQuadX,   pointInQuadY   };         float[] distanceToTopRight    = new float[]{ pointInQuadX-1, pointInQuadY   };         float[] distanceToBottomLeft  = new float[]{ pointInQuadX,   pointInQuadY-1 };         float[] distanceToBottomRight = new float[]{ pointInQuadX-1, pointInQuadY-1 };          // считаем скалярные произведения между которыми будем интерполировать /*  tx1--tx2   |    |  bx1--bx2 */         float tx1 = Dot(distanceToTopLeft,     topLeftGradient);         float tx2 = Dot(distanceToTopRight,    topRightGradient);         float bx1 = Dot(distanceToBottomLeft,  bottomLeftGradient);         float bx2 = Dot(distanceToBottomRight, bottomRightGradient);          // готовим параметры интерполяции, чтобы она не была линейной:         pointInQuadX = QunticCurve(pointInQuadX);         pointInQuadY = QunticCurve(pointInQuadY);          // собственно, интерполяция:         float tx = Lerp(tx1, tx2, pointInQuadX);         float bx = Lerp(bx1, bx2, pointInQuadX);         float tb = Lerp(tx, bx, pointInQuadY);          // возвращаем результат:         return tb;     } 

В качестве бонуса:

мультиоктавный шум

    public float Noise(float fx, float fy, int octaves, float persistence = 0.5f)     {         float amplitude = 1; // сила применения шума к общей картине, будет уменьшаться с "мельчанием" шума         // как сильно уменьшаться - регулирует persistence         float max = 0; // необходимо для нормализации результата         float result = 0; // накопитель результата          while (octaves-- > 0)         {             max += amplitude;             result += Noise(fx, fy) * amplitude;             amplitude *= persistence;             fx *= 2; // удваиваем частоту шума (делаем его более мелким) с каждой октавой             fy *= 2;         }          return result/max;     } 

И последнее — использование таблицы со случайными числами. В коде Кена Перлина такая таблица прописана вручную и достаются оттуда значения совсем по-другому. Здесь можно экспериментировать и от этого немало зависит равномерность шума и отсутствие в нём явных паттернов.

Я сделал

так

class Perlin2D {     byte[] permutationTable;      public Perlin2D(int seed = 0)     {         var rand = new System.Random(seed);         permutationTable = new byte[1023];         rand.NextBytes(permutationTable); // заполняем случайными байтами     }      private float[] GetPseudoRandomGradientVector(int x, int y)     { // хэш-функция с Простыми числами, обрезкой результата до размера массива со случайными байтами         int v = (int)(((x * 1836311903) ^ (y * 2971215073) + 4807526976) & 1023);         v = permutationTable[v]&3;          switch (v)         {             ... 

& 3 здесь обрезает любое int32 число до 3, читайте об операции AND на википедии
Операция типа % 3 тоже сработала бы, но намного медленней.

Исходный код целиком (без комментариев)

class Perlin2D {     byte[] permutationTable;      public Perlin2D(int seed = 0)     {         var rand = new System.Random(seed);         permutationTable = new byte[1023];         rand.NextBytes(permutationTable);     }      private float[] GetPseudoRandomGradientVector(int x, int y)     {         int v = (int)(((x * 1836311903) ^ (y * 2971215073) + 4807526976) & 1023);         v = permutationTable[v]&3;          switch (v)         {             case 0:  return new float[]{  1, 0 };             case 1:  return new float[]{ -1, 0 };             case 2:  return new float[]{  0, 1 };             default: return new float[]{  0,-1 };         }     }      static float QunticCurve(float t)     {         return t * t * t * (t * (t * 6 - 15) + 10);     }      static float Lerp(float a, float b, float t)     {         return a + (b - a) * t;     }      static float Dot(float[] a, float[] b)     {         return a[0] * b[0] + a[1] * b[1];     }      public float Noise(float fx, float fy)     {         int left = (int)System.Math.Floor(fx);         int top  = (int)System.Math.Floor(fy);         float pointInQuadX = fx - left;         float pointInQuadY = fy - top;          float[] topLeftGradient     = GetPseudoRandomGradientVector(left,   top  );         float[] topRightGradient    = GetPseudoRandomGradientVector(left+1, top  );         float[] bottomLeftGradient  = GetPseudoRandomGradientVector(left,   top+1);         float[] bottomRightGradient = GetPseudoRandomGradientVector(left+1, top+1);          float[] distanceToTopLeft     = new float[]{ pointInQuadX,   pointInQuadY   };         float[] distanceToTopRight    = new float[]{ pointInQuadX-1, pointInQuadY   };         float[] distanceToBottomLeft  = new float[]{ pointInQuadX,   pointInQuadY-1 };         float[] distanceToBottomRight = new float[]{ pointInQuadX-1, pointInQuadY-1 };          float tx1 = Dot(distanceToTopLeft,     topLeftGradient);         float tx2 = Dot(distanceToTopRight,    topRightGradient);         float bx1 = Dot(distanceToBottomLeft,  bottomLeftGradient);         float bx2 = Dot(distanceToBottomRight, bottomRightGradient);          pointInQuadX = QunticCurve(pointInQuadX);         pointInQuadY = QunticCurve(pointInQuadY);          float tx = Lerp(tx1, tx2, pointInQuadX);         float bx = Lerp(bx1, bx2, pointInQuadX);         float tb = Lerp(tx, bx, pointInQuadY);          return tb;     }      public float Noise(float fx, float fy, int octaves, float persistence = 0.5f)     {         float amplitude = 1;         float max = 0;         float result = 0;          while (octaves-- > 0)         {             max += amplitude;             result += Noise(fx, fy) * amplitude;             amplitude *= persistence;             fx *= 2;             fy *= 2;         }          return result/max;     } } 

Результат:
image

ссылка на оригинал статьи http://habrahabr.ru/post/265775/


Комментарии

Добавить комментарий

Ваш адрес email не будет опубликован. Обязательные поля помечены *