src/noise.cpp

   1 #include "noise.hpp"
   2
   3 #include <cmath>
   4
   5
   6 namespace {
   7
   8 constexpr float one_third = 1.0f/3.0f;
   9 constexpr float one_sixth = 1.0f/6.0f;
  10
  11 }
  12
  13 namespace blank {
  14
  15 GaloisLFSR::GaloisLFSR(std::uint64_t seed) noexcept
  16 : state(seed) {
  17
  18 }
  19
  20 bool GaloisLFSR::operator ()() noexcept {
  21         bool result = state & 1;
  22         state >>= 1;
  23         if (result) {
  24                 state |= 0x8000000000000000;
  25                 state ^= mask;
  26         } else {
  27                 state &= 0x7FFFFFFFFFFFFFFF;
  28         }
  29         return result;
  30 }
  31
  32
  33 SimplexNoise::SimplexNoise(unsigned int seed) noexcept
  34 : grad({
  35         {  1.0f,  1.0f,  0.0f },
  36         { -1.0f,  1.0f,  0.0f },
  37         {  1.0f, -1.0f,  0.0f },
  38         { -1.0f, -1.0f,  0.0f },
  39         {  1.0f,  0.0f,  1.0f },
  40         { -1.0f,  0.0f,  1.0f },
  41         {  1.0f,  0.0f, -1.0f },
  42         { -1.0f,  0.0f, -1.0f },
  43         {  0.0f,  1.0f,  1.0f },
  44         {  0.0f, -1.0f,  1.0f },
  45         {  0.0f,  1.0f, -1.0f },
  46         {  0.0f, -1.0f, -1.0f },
  47 }) {
  48         GaloisLFSR random(seed ^ 0x0123456789ACBDEF);
  49         unsigned char value;
  50         for (size_t i = 0; i < 256; ++i) {
  51                 perm[i] = random(value);
  52                 perm[i] &= 0xFF;
  53                 perm[i + 256] = perm[i];
  54                 perm12[i] = perm[i] % 12;
  55                 perm12[i + 256] = perm12[i];
  56         }
  57 }
  58
  59
  60 float SimplexNoise::operator ()(const glm::vec3 &in) const noexcept {
  61         float skew = (in.x + in.y + in.z) * one_third;
  62
  63         glm::vec3 skewed(glm::floor(in + skew));
  64         float tr = (skewed.x + skewed.y + skewed.z) * one_sixth;
  65
  66         glm::vec3 unskewed(skewed - tr);
  67         glm::vec3 relative(in - unskewed);
  68
  69         glm::vec3 second, third;
  70
  71         if (relative.x >= relative.y) {
  72                 if (relative.y >= relative.z) {
  73                         second = { 1, 0, 0 };
  74                         third = { 1, 1, 0 };
  75                 } else if (relative.x >= relative.z) {
  76                         second = { 1, 0, 0 };
  77                         third = { 1, 0, 1 };
  78                 } else {
  79                         second = { 0, 0, 1 };
  80                         third = { 1, 0, 1 };
  81                 }
  82         } else if (relative.y < relative.z) {
  83                 second = { 0, 0, 1 };
  84                 third = { 0, 1, 1 };
  85         } else if (relative.x < relative.z) {
  86                 second = { 0, 1, 0 };
  87                 third = { 0, 1, 1 };
  88         } else {
  89                 second = { 0, 1, 0 };
  90                 third = { 1, 1, 0 };
  91         }
  92
  93         glm::vec3 offset[4] = {
  94                 in - unskewed,
  95                 relative - second + one_sixth,
  96                 relative - third + one_third,
  97                 relative - 0.5f,
  98         };
  99
 100         int index[3] = {
 101                 (int)(skewed.x) & 0xFF,
 102                 (int)(skewed.y) & 0xFF,
 103                 (int)(skewed.z) & 0xFF,
 104         };
 105
 106         float n = 0.0f;
 107
 108         // 0
 109         float t = 0.6f - dot(offset[0], offset[0]);
 110         if (t > 0.0f) {
 111                 t *= t;
 112                 int corner = Perm12(index[0] + Perm(index[1] + Perm(index[2])));
 113                 n += t * t * dot(Grad(corner), offset[0]);
 114         }
 115
 116         // 1
 117         t = 0.6f - dot(offset[1], offset[1]);
 118         if (t > 0.0f) {
 119                 t *= t;
 120                 int corner = Perm12(index[0] + int(second.x) + Perm(index[1] + int(second.y) + Perm(index[2] + int(second.z))));
 121                 n += t * t * dot(Grad(corner), offset[1]);
 122         }
 123
 124         // 2
 125         t = 0.6f - dot(offset[2], offset[2]);
 126         if (t > 0.0f) {
 127                 t *= t;
 128                 int corner = Perm12(index[0] + int(third.x) + Perm(index[1] + int(third.y) + Perm(index[2] + int(third.z))));
 129                 n += t * t * dot(Grad(corner), offset[2]);
 130         }
 131
 132         // 3
 133         t = 0.6f - dot(offset[3], offset[3]);
 134         if (t > 0.0f) {
 135                 t *= t;
 136                 int corner = Perm12(index[0] + 1 + Perm(index[1] + 1 + Perm(index[2] + 1)));
 137                 n += t * t * dot(Grad(corner), offset[3]);
 138         }
 139
 140         return 32.0f * n;
 141 }
 142
 143
 144 int SimplexNoise::Perm(int idx) const noexcept {
 145         return perm[idx];
 146 }
 147
 148 int SimplexNoise::Perm12(int idx) const noexcept {
 149         return perm12[idx];
 150 }
 151
 152 const glm::vec3 &SimplexNoise::Grad(int idx) const noexcept {
 153         return grad[idx];
 154 }
 155
 156
 157 WorleyNoise::WorleyNoise(unsigned int seed) noexcept
 158 : seed(seed)
 159 , num_points(8) {
 160
 161 }
 162
 163 float WorleyNoise::operator ()(const glm::vec3 &in) const noexcept {
 164         glm::vec3 center = floor(in);
 165
 166         float closest = 1.0f;  // cannot be farther away than 1.0
 167
 168         for (int z = -1; z <= 1; ++z) {
 169                 for (int y = -1; y <= 1; ++y) {
 170                         for (int x = -1; x <= 1; ++x) {
 171                                 glm::vec3 cube(center.x + x, center.y + y, center.z + z);
 172                                 unsigned int cube_rand =
 173                                         (unsigned(cube.x) * 130223) ^
 174                                         (unsigned(cube.y) * 159899) ^
 175                                         (unsigned(cube.z) * 190717) ^
 176                                         seed;
 177
 178                                 for (int i = 0; i < num_points; ++i) {
 179                                         glm::vec3 point(cube);
 180                                         cube_rand = 190667 * cube_rand + 109807;
 181                                         point.x += float(cube_rand % 262144) / 262144.0f;
 182                                         cube_rand = 135899 * cube_rand + 189169;
 183                                         point.y += float(cube_rand % 262144) / 262144.0f;
 184                                         cube_rand = 159739 * cube_rand + 112139;
 185                                         point.z += float(cube_rand % 262144) / 262144.0f;
 186
 187                                         glm::vec3 diff(in - point);
 188                                         float distance = sqrt(dot(diff, diff));
 189                                         if (distance < closest) {
 190                                                 closest = distance;
 191                                         }
 192                                 }
 193                         }
 194                 }
 195         }
 196
 197         // closest ranges (0, 1), so normalizing to (-1,1) is trivial
 198         // though heavily biased towards lower numbers
 199         return 2.0f * closest - 1.0f;
 200 }
 201
 202 }