KSeExpr 6.0.0.0
Noise.cpp
Go to the documentation of this file.
1// SPDX-FileCopyrightText: 2011-2019 Disney Enterprises, Inc.
2// SPDX-License-Identifier: LicenseRef-Apache-2.0
3// SPDX-FileCopyrightText: 2020 L. E. Segovia <amy@amyspark.me>
4// SPDX-License-Identifier: GPL-3.0-or-later
5
6#include "ExprBuiltins.h"
7#include <iostream>
8
9namespace
10{
11#include "NoiseTables.h"
12}
13#include "Noise.h"
14#include "Utils.h"
15namespace KSeExpr
16{
18double s_curve(double t)
19{
20 return t * t * t * (t * (6 * t - 15) + 10);
21}
22
24template<int d> unsigned char hashReduceChar(int index[d])
25{
26 uint32_t seed = 0;
27 // blend with seed (constants from Numerical Recipes, attrib. from Knuth)
28 for (int k = 0; k < d; k++) {
29 static const uint32_t M = 1664525, C = 1013904223;
30 seed = seed * M + index[k] + C;
31 }
32 // tempering (from Matsumoto)
33 seed ^= (seed >> 11);
34 seed ^= (seed << 7) & 0x9d2c5680UL;
35 seed ^= (seed << 15) & 0xefc60000UL;
36 seed ^= (seed >> 18);
37 // compute one byte by mixing third and first bytes
38 return (((seed & 0xff0000) >> 4) + (seed & 0xff)) & 0xff;
39}
40
42template<int d> uint32_t hashReduce(uint32_t index[d])
43{
44 union {
45 uint32_t i;
46 unsigned char c[4];
47 } u1, u2;
48 // blend with seed (constants from Numerical Recipes, attrib. from Knuth)
49 u1.i = 0;
50 for (int k = 0; k < d; k++) {
51 static const uint32_t M = 1664525, C = 1013904223;
52 u1.i = u1.i * M + index[k] + C;
53 }
54 // tempering (from Matsumoto)
55 u1.i ^= (u1.i >> 11);
56 u1.i ^= (u1.i << 7) & 0x9d2c5680U;
57 u1.i ^= (u1.i << 15) & 0xefc60000U;
58 u1.i ^= (u1.i >> 18);
59 // permute bytes (shares perlin noise permutation table)
60 u2.c[3] = p[u1.c[0]];
61 u2.c[2] = p[u1.c[1] + u2.c[3]];
62 u2.c[1] = p[u1.c[2] + u2.c[2]];
63 u2.c[0] = p[u1.c[3] + u2.c[1]];
64
65 return u2.i;
66}
67
69template<int d, class T, bool periodic> T noiseHelper(const T *X, const int *period = 0)
70{
71 // find lattice index
72 T weights[2][d]; // lower and upper weights
73 int index[d];
74 for (int k = 0; k < d; k++) {
76 index[k] = (int)f;
77 if (periodic) {
78 index[k] %= period[k];
79 if (index[k] < 0)
80 index[k] += period[k];
81 }
82 weights[0][k] = X[k] - f;
83 weights[1][k] = weights[0][k] - 1; // dist to cell with index one above
84 }
85 // compute function values propagated from zero from each node
86 const int num = 1 << d;
87 T vals[num];
88 for (int dummy = 0; dummy < num; dummy++) {
89 int latticeIndex[d];
90 int offset[d];
91 for (int k = 0; k < d; k++) {
92 offset[k] = ((dummy & (1 << k)) != 0);
93 latticeIndex[k] = index[k] + offset[k];
94 }
95 // hash to get representative gradient vector
96 int lookup = hashReduceChar<d>(latticeIndex);
97 T val = 0;
98 for (int k = 0; k < d; k++) {
99 double grad = NOISE_TABLES<d>::g[lookup][k];
100 double weight = weights[offset[k]][k];
101 val += grad * weight;
102 }
103 vals[dummy] = val;
104 }
105 // compute linear interpolation coefficients
106 T alphas[d];
107 for (int k = 0; k < d; k++)
108 alphas[k] = s_curve(weights[0][k]);
109 // perform multilinear interpolation (i.e. linear, bilinear, trilinear, quadralinear)
110 for (int newd = d - 1; newd >= 0; newd--) {
111 int newnum = 1 << newd;
112 int k = (d - newd - 1);
113 T alpha = alphas[k];
114 T beta = T(1) - alphas[k];
115 for (int dummy = 0; dummy < newnum; dummy++) {
116 int index = dummy * (1 << (d - newd));
117 int otherIndex = index + (1 << k);
118 // T alpha=s_curve(weights[0][k]);
119 vals[index] = beta * vals[index] + alpha * vals[otherIndex];
120 }
121 }
122 // return reduced version
123 return vals[0];
124}
125}
126
127namespace KSeExpr
128{
130template<int d_in, int d_out, class T> void CellNoise(const T *in, T *out)
131{
132 uint32_t index[d_in];
133 int dim = 0;
134 for (int k = 0; k < d_in; k++)
135 index[k] = uint32_t(KSeExpr::Utils::floor(in[k]));
136 while (1) {
137 out[dim] = hashReduce<d_in>(index) * (1.0 / 0xffffffffu);
138 if (++dim >= d_out)
139 break;
140 for (int k = 0; k < d_in; k++)
141 index[k] += 1000;
142 }
143}
144
146template<int d_in, int d_out, class T> void Noise(const T *in, T *out)
147{
148 T P[d_in];
149 for (int i = 0; i < d_in; i++)
150 P[i] = in[i];
151
152 int i = 0;
153 while (1) {
155 if (++i >= d_out)
156 break;
157 for (int k = 0; k < d_out; k++)
158 P[k] += (T)1000;
159 }
160}
161
163template<int d_in, int d_out, class T> void PNoise(const T *in, const int *period, T *out)
164{
165 T P[d_in];
166 for (int i = 0; i < d_in; i++)
167 P[i] = in[i];
168
169 int i = 0;
170 while (1) {
172 if (++i >= d_out)
173 break;
174 for (int k = 0; k < d_out; k++)
175 P[k] += (T)1000;
176 }
177}
178
181template<int d_in, int d_out, bool turbulence, class T> void FBM(const T *in, T *out, int octaves, T lacunarity, T gain)
182{
183 T P[d_in];
184 for (int i = 0; i < d_in; i++)
185 P[i] = in[i];
186
187 T scale = 1;
188 for (int k = 0; k < d_out; k++)
189 out[k] = 0;
190 int octave = 0;
191 while (1) {
194 if (turbulence)
195 for (int k = 0; k < d_out; k++)
196 out[k] += fabs(localResult[k]) * scale;
197 else
198 for (int k = 0; k < d_out; k++)
199 out[k] += localResult[k] * scale;
200 if (++octave >= octaves)
201 break;
202 scale *= gain;
203 for (int k = 0; k < d_in; k++) {
204 P[k] *= lacunarity;
205 P[k] += (T)1234;
206 }
207 }
208}
209
210// Explicit instantiations
211template void CellNoise<3, 1, double>(const double *, double *);
212template void CellNoise<3, 3, double>(const double *, double *);
213template void Noise<1, 1, double>(const double *, double *);
214template void Noise<2, 1, double>(const double *, double *);
215template void Noise<3, 1, double>(const double *, double *);
216template void PNoise<3, 1, double>(const double *, const int *, double *);
217template void Noise<4, 1, double>(const double *, double *);
218template void Noise<3, 3, double>(const double *, double *);
219template void Noise<4, 3, double>(const double *, double *);
220template void FBM<3, 1, false, double>(const double *, double *, int, double, double);
221template void FBM<3, 1, true, double>(const double *, double *, int, double, double);
222template void FBM<3, 3, false, double>(const double *, double *, int, double, double);
223template void FBM<3, 3, true, double>(const double *, double *, int, double, double);
224template void FBM<4, 1, false, double>(const double *, double *, int, double, double);
225template void FBM<4, 3, false, double>(const double *, double *, int, double, double);
226}
227
228#ifdef MAINTEST
229int main(int argc, char *argv[])
230{
231 typedef double T;
232 T sum = 0;
233 for (int i = 0; i < 10000000; i++) {
234 T foo[3] = {.3, .3, .3};
235 // for(int k=0;k<3;k++) foo[k]=(double)rand()/double(RAND_MAX)*100.;
237 }
238}
239#endif
static constexpr std::array< int, 514 > p
Definition NoiseTables.h:10
int main(int, char *[])
KSeExpr_DEFAULT double_t floor(double_t val)
Definition Utils.cpp:168
double s_curve(double t)
This is the Quintic interpolant from Perlin's Improved Noise Paper.
Definition Noise.cpp:18
template void Noise< 2, 1, double >(const double *, double *)
template void FBM< 4, 1, false, double >(const double *, double *, int, double, double)
template void PNoise< 3, 1, double >(const double *, const int *, double *)
void Noise(const T *in, T *out)
Noise with d_in dimensional domain, d_out dimensional abcissa.
Definition Noise.cpp:146
template void FBM< 3, 3, false, double >(const double *, double *, int, double, double)
template void FBM< 3, 1, true, double >(const double *, double *, int, double, double)
template void Noise< 1, 1, double >(const double *, double *)
template void FBM< 3, 1, false, double >(const double *, double *, int, double, double)
template void CellNoise< 3, 1, double >(const double *, double *)
void CellNoise(const T *in, T *out)
Computes cellular noise (non-interpolated piecewise constant cell random values)
Definition Noise.cpp:130
template void CellNoise< 3, 3, double >(const double *, double *)
T noiseHelper(const T *X, const int *period=0)
Noise with d_in dimensional domain, 1 dimensional abcissa.
Definition Noise.cpp:69
unsigned char hashReduceChar(int index[d])
Does a hash reduce to a character.
Definition Noise.cpp:24
void PNoise(const T *in, const int *period, T *out)
Periodic Noise with d_in dimensional domain, d_out dimensional abcissa.
Definition Noise.cpp:163
template void FBM< 4, 3, false, double >(const double *, double *, int, double, double)
double turbulence(int n, const Vec3d *args)
void FBM(const T *in, T *out, int octaves, T lacunarity, T gain)
Fractional Brownian Motion. If turbulence is true then turbulence computed.
Definition Noise.cpp:181
uint32_t hashReduce(uint32_t index[d])
Does a hash reduce to an integer.
Definition Noise.cpp:42
template void Noise< 3, 3, double >(const double *, double *)
template void FBM< 3, 3, true, double >(const double *, double *, int, double, double)
template void Noise< 3, 1, double >(const double *, double *)
template void Noise< 4, 1, double >(const double *, double *)
template void Noise< 4, 3, double >(const double *, double *)