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