2
* Mesa 3-D graphics library
5
* Copyright (C) 2006 Brian Paul All Rights Reserved.
7
* Permission is hereby granted, free of charge, to any person obtaining a
8
* copy of this software and associated documentation files (the "Software"),
9
* to deal in the Software without restriction, including without limitation
10
* the rights to use, copy, modify, merge, publish, distribute, sublicense,
11
* and/or sell copies of the Software, and to permit persons to whom the
12
* Software is furnished to do so, subject to the following conditions:
14
* The above copyright notice and this permission notice shall be included
15
* in all copies or substantial portions of the Software.
17
* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
18
* OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
19
* FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
20
* BRIAN PAUL BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN
21
* AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
22
* CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
27
* Copyright (c) 2003-2005, Stefan Gustavson
29
* Contact: stegu@itn.liu.se
34
* \brief C implementation of Perlin Simplex Noise over 1, 2, 3 and 4 dims.
35
* \author Stefan Gustavson (stegu@itn.liu.se)
38
* This implementation is "Simplex Noise" as presented by
39
* Ken Perlin at a relatively obscure and not often cited course
40
* session "Real-Time Shading" at Siggraph 2001 (before real
41
* time shading actually took on), under the title "hardware noise".
42
* The 3D function is numerically equivalent to his Java reference
43
* code available in the PDF course notes, although I re-implemented
44
* it from scratch to get more readable code. The 1D, 2D and 4D cases
45
* were implemented from scratch by me from Ken Perlin's text.
47
* This file has no dependencies on any other file, not even its own
48
* header file. The header file is made for use by external code only.
52
#include "main/imports.h"
53
#include "prog_noise.h"
55
#define FASTFLOOR(x) ( ((x)>0) ? ((int)x) : (((int)x)-1) )
58
* ---------------------------------------------------------------------
63
* Permutation table. This is just a random jumble of all numbers 0-255,
64
* repeated twice to avoid wrapping the index at 255 for each lookup.
65
* This needs to be exactly the same for all instances on all platforms,
66
* so it's easiest to just keep it as static explicit data.
67
* This also removes the need for any initialisation of this class.
69
* Note that making this an int[] instead of a char[] might make the
70
* code run faster on platforms with a high penalty for unaligned single
71
* byte addressing. Intel x86 is generally single-byte-friendly, but
72
* some other CPUs are faster with 4-aligned reads.
73
* However, a char[] is smaller, which avoids cache trashing, and that
74
* is probably the most important aspect on most architectures.
75
* This array is accessed a *lot* by the noise functions.
76
* A vector-valued noise over 3D accesses it 96 times, and a
77
* float-valued 4D noise 64 times. We want this to fit in the cache!
79
unsigned char perm[512] = { 151, 160, 137, 91, 90, 15,
80
131, 13, 201, 95, 96, 53, 194, 233, 7, 225, 140, 36, 103, 30, 69, 142, 8,
81
99, 37, 240, 21, 10, 23,
82
190, 6, 148, 247, 120, 234, 75, 0, 26, 197, 62, 94, 252, 219, 203, 117, 35,
84
88, 237, 149, 56, 87, 174, 20, 125, 136, 171, 168, 68, 175, 74, 165, 71,
85
134, 139, 48, 27, 166,
86
77, 146, 158, 231, 83, 111, 229, 122, 60, 211, 133, 230, 220, 105, 92, 41,
88
102, 143, 54, 65, 25, 63, 161, 1, 216, 80, 73, 209, 76, 132, 187, 208, 89,
90
135, 130, 116, 188, 159, 86, 164, 100, 109, 198, 173, 186, 3, 64, 52, 217,
92
5, 202, 38, 147, 118, 126, 255, 82, 85, 212, 207, 206, 59, 227, 47, 16, 58,
94
223, 183, 170, 213, 119, 248, 152, 2, 44, 154, 163, 70, 221, 153, 101, 155,
96
129, 22, 39, 253, 19, 98, 108, 110, 79, 113, 224, 232, 178, 185, 112, 104,
98
251, 34, 242, 193, 238, 210, 144, 12, 191, 179, 162, 241, 81, 51, 145, 235,
100
49, 192, 214, 31, 181, 199, 106, 157, 184, 84, 204, 176, 115, 121, 50, 45,
102
138, 236, 205, 93, 222, 114, 67, 29, 24, 72, 243, 141, 128, 195, 78, 66,
104
151, 160, 137, 91, 90, 15,
105
131, 13, 201, 95, 96, 53, 194, 233, 7, 225, 140, 36, 103, 30, 69, 142, 8,
106
99, 37, 240, 21, 10, 23,
107
190, 6, 148, 247, 120, 234, 75, 0, 26, 197, 62, 94, 252, 219, 203, 117, 35,
109
88, 237, 149, 56, 87, 174, 20, 125, 136, 171, 168, 68, 175, 74, 165, 71,
110
134, 139, 48, 27, 166,
111
77, 146, 158, 231, 83, 111, 229, 122, 60, 211, 133, 230, 220, 105, 92, 41,
112
55, 46, 245, 40, 244,
113
102, 143, 54, 65, 25, 63, 161, 1, 216, 80, 73, 209, 76, 132, 187, 208, 89,
115
135, 130, 116, 188, 159, 86, 164, 100, 109, 198, 173, 186, 3, 64, 52, 217,
117
5, 202, 38, 147, 118, 126, 255, 82, 85, 212, 207, 206, 59, 227, 47, 16, 58,
118
17, 182, 189, 28, 42,
119
223, 183, 170, 213, 119, 248, 152, 2, 44, 154, 163, 70, 221, 153, 101, 155,
121
129, 22, 39, 253, 19, 98, 108, 110, 79, 113, 224, 232, 178, 185, 112, 104,
123
251, 34, 242, 193, 238, 210, 144, 12, 191, 179, 162, 241, 81, 51, 145, 235,
125
49, 192, 214, 31, 181, 199, 106, 157, 184, 84, 204, 176, 115, 121, 50, 45,
127
138, 236, 205, 93, 222, 114, 67, 29, 24, 72, 243, 141, 128, 195, 78, 66,
132
* ---------------------------------------------------------------------
136
* Helper functions to compute gradients-dot-residualvectors (1D to 4D)
137
* Note that these generate gradients of more than unit length. To make
138
* a close match with the value range of classic Perlin noise, the final
139
* noise values need to be rescaled to fit nicely within [-1,1].
140
* (The simplex noise functions as such also have different scaling.)
141
* Note also that these noise functions are the most practical and useful
142
* signed version of Perlin noise. To return values according to the
143
* RenderMan specification from the SL noise() and pnoise() functions,
144
* the noise values need to be scaled and offset to [0,1], like this:
145
* float SLnoise = (SimplexNoise1234::noise(x,y,z) + 1.0) * 0.5;
149
grad1(int hash, float x)
152
float grad = 1.0f + (h & 7); /* Gradient value 1.0, 2.0, ..., 8.0 */
154
grad = -grad; /* Set a random sign for the gradient */
155
return (grad * x); /* Multiply the gradient with the distance */
159
grad2(int hash, float x, float y)
161
int h = hash & 7; /* Convert low 3 bits of hash code */
162
float u = h < 4 ? x : y; /* into 8 simple gradient directions, */
163
float v = h < 4 ? y : x; /* and compute the dot product with (x,y). */
164
return ((h & 1) ? -u : u) + ((h & 2) ? -2.0f * v : 2.0f * v);
168
grad3(int hash, float x, float y, float z)
170
int h = hash & 15; /* Convert low 4 bits of hash code into 12 simple */
171
float u = h < 8 ? x : y; /* gradient directions, and compute dot product. */
172
float v = h < 4 ? y : h == 12 || h == 14 ? x : z; /* Fix repeats at h = 12 to 15 */
173
return ((h & 1) ? -u : u) + ((h & 2) ? -v : v);
177
grad4(int hash, float x, float y, float z, float t)
179
int h = hash & 31; /* Convert low 5 bits of hash code into 32 simple */
180
float u = h < 24 ? x : y; /* gradient directions, and compute dot product. */
181
float v = h < 16 ? y : z;
182
float w = h < 8 ? z : t;
183
return ((h & 1) ? -u : u) + ((h & 2) ? -v : v) + ((h & 4) ? -w : w);
187
* A lookup table to traverse the simplex around a given point in 4D.
188
* Details can be found where this table is used, in the 4D noise method.
189
* TODO: This should not be required, backport it from Bill's GLSL code!
191
static unsigned char simplex[64][4] = {
192
{0, 1, 2, 3}, {0, 1, 3, 2}, {0, 0, 0, 0}, {0, 2, 3, 1},
193
{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {1, 2, 3, 0},
194
{0, 2, 1, 3}, {0, 0, 0, 0}, {0, 3, 1, 2}, {0, 3, 2, 1},
195
{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {1, 3, 2, 0},
196
{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0},
197
{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0},
198
{1, 2, 0, 3}, {0, 0, 0, 0}, {1, 3, 0, 2}, {0, 0, 0, 0},
199
{0, 0, 0, 0}, {0, 0, 0, 0}, {2, 3, 0, 1}, {2, 3, 1, 0},
200
{1, 0, 2, 3}, {1, 0, 3, 2}, {0, 0, 0, 0}, {0, 0, 0, 0},
201
{0, 0, 0, 0}, {2, 0, 3, 1}, {0, 0, 0, 0}, {2, 1, 3, 0},
202
{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0},
203
{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0},
204
{2, 0, 1, 3}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0},
205
{3, 0, 1, 2}, {3, 0, 2, 1}, {0, 0, 0, 0}, {3, 1, 2, 0},
206
{2, 1, 0, 3}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0},
207
{3, 1, 0, 2}, {0, 0, 0, 0}, {3, 2, 0, 1}, {3, 2, 1, 0}
211
/** 1D simplex noise */
213
_mesa_noise1(GLfloat x)
215
int i0 = FASTFLOOR(x);
218
float x1 = x0 - 1.0f;
219
float t1 = 1.0f - x1 * x1;
222
float t0 = 1.0f - x0 * x0;
223
/* if(t0 < 0.0f) t0 = 0.0f; // this never happens for the 1D case */
225
n0 = t0 * t0 * grad1(perm[i0 & 0xff], x0);
227
/* if(t1 < 0.0f) t1 = 0.0f; // this never happens for the 1D case */
229
n1 = t1 * t1 * grad1(perm[i1 & 0xff], x1);
230
/* The maximum value of this noise is 8*(3/4)^4 = 2.53125 */
231
/* A factor of 0.395 would scale to fit exactly within [-1,1], but */
232
/* we want to match PRMan's 1D noise, so we scale it down some more. */
233
return 0.25f * (n0 + n1);
237
/** 2D simplex noise */
239
_mesa_noise2(GLfloat x, GLfloat y)
241
#define F2 0.366025403f /* F2 = 0.5*(sqrt(3.0)-1.0) */
242
#define G2 0.211324865f /* G2 = (3.0-Math.sqrt(3.0))/6.0 */
244
float n0, n1, n2; /* Noise contributions from the three corners */
246
/* Skew the input space to determine which simplex cell we're in */
247
float s = (x + y) * F2; /* Hairy factor for 2D */
250
int i = FASTFLOOR(xs);
251
int j = FASTFLOOR(ys);
253
float t = (float) (i + j) * G2;
254
float X0 = i - t; /* Unskew the cell origin back to (x,y) space */
256
float x0 = x - X0; /* The x,y distances from the cell origin */
259
float x1, y1, x2, y2;
263
/* For the 2D case, the simplex shape is an equilateral triangle. */
264
/* Determine which simplex we are in. */
265
int i1, j1; /* Offsets for second (middle) corner of simplex in (i,j) coords */
269
} /* lower triangle, XY order: (0,0)->(1,0)->(1,1) */
273
} /* upper triangle, YX order: (0,0)->(0,1)->(1,1) */
275
/* A step of (1,0) in (i,j) means a step of (1-c,-c) in (x,y), and */
276
/* a step of (0,1) in (i,j) means a step of (-c,1-c) in (x,y), where */
277
/* c = (3-sqrt(3))/6 */
279
x1 = x0 - i1 + G2; /* Offsets for middle corner in (x,y) unskewed coords */
281
x2 = x0 - 1.0f + 2.0f * G2; /* Offsets for last corner in (x,y) unskewed coords */
282
y2 = y0 - 1.0f + 2.0f * G2;
284
/* Wrap the integer indices at 256, to avoid indexing perm[] out of bounds */
288
/* Calculate the contribution from the three corners */
289
t0 = 0.5f - x0 * x0 - y0 * y0;
294
n0 = t0 * t0 * grad2(perm[ii + perm[jj]], x0, y0);
297
t1 = 0.5f - x1 * x1 - y1 * y1;
302
n1 = t1 * t1 * grad2(perm[ii + i1 + perm[jj + j1]], x1, y1);
305
t2 = 0.5f - x2 * x2 - y2 * y2;
310
n2 = t2 * t2 * grad2(perm[ii + 1 + perm[jj + 1]], x2, y2);
313
/* Add contributions from each corner to get the final noise value. */
314
/* The result is scaled to return values in the interval [-1,1]. */
315
return 40.0f * (n0 + n1 + n2); /* TODO: The scale factor is preliminary! */
319
/** 3D simplex noise */
321
_mesa_noise3(GLfloat x, GLfloat y, GLfloat z)
323
/* Simple skewing factors for the 3D case */
324
#define F3 0.333333333f
325
#define G3 0.166666667f
327
float n0, n1, n2, n3; /* Noise contributions from the four corners */
329
/* Skew the input space to determine which simplex cell we're in */
330
float s = (x + y + z) * F3; /* Very nice and simple skew factor for 3D */
334
int i = FASTFLOOR(xs);
335
int j = FASTFLOOR(ys);
336
int k = FASTFLOOR(zs);
338
float t = (float) (i + j + k) * G3;
339
float X0 = i - t; /* Unskew the cell origin back to (x,y,z) space */
342
float x0 = x - X0; /* The x,y,z distances from the cell origin */
346
float x1, y1, z1, x2, y2, z2, x3, y3, z3;
348
float t0, t1, t2, t3;
350
/* For the 3D case, the simplex shape is a slightly irregular tetrahedron. */
351
/* Determine which simplex we are in. */
352
int i1, j1, k1; /* Offsets for second corner of simplex in (i,j,k) coords */
353
int i2, j2, k2; /* Offsets for third corner of simplex in (i,j,k) coords */
355
/* This code would benefit from a backport from the GLSL version! */
409
/* A step of (1,0,0) in (i,j,k) means a step of (1-c,-c,-c) in
410
* (x,y,z), a step of (0,1,0) in (i,j,k) means a step of
411
* (-c,1-c,-c) in (x,y,z), and a step of (0,0,1) in (i,j,k) means a
412
* step of (-c,-c,1-c) in (x,y,z), where c = 1/6.
415
x1 = x0 - i1 + G3; /* Offsets for second corner in (x,y,z) coords */
418
x2 = x0 - i2 + 2.0f * G3; /* Offsets for third corner in (x,y,z) coords */
419
y2 = y0 - j2 + 2.0f * G3;
420
z2 = z0 - k2 + 2.0f * G3;
421
x3 = x0 - 1.0f + 3.0f * G3;/* Offsets for last corner in (x,y,z) coords */
422
y3 = y0 - 1.0f + 3.0f * G3;
423
z3 = z0 - 1.0f + 3.0f * G3;
425
/* Wrap the integer indices at 256 to avoid indexing perm[] out of bounds */
430
/* Calculate the contribution from the four corners */
431
t0 = 0.6f - x0 * x0 - y0 * y0 - z0 * z0;
436
n0 = t0 * t0 * grad3(perm[ii + perm[jj + perm[kk]]], x0, y0, z0);
439
t1 = 0.6f - x1 * x1 - y1 * y1 - z1 * z1;
445
t1 * t1 * grad3(perm[ii + i1 + perm[jj + j1 + perm[kk + k1]]], x1,
449
t2 = 0.6f - x2 * x2 - y2 * y2 - z2 * z2;
455
t2 * t2 * grad3(perm[ii + i2 + perm[jj + j2 + perm[kk + k2]]], x2,
459
t3 = 0.6f - x3 * x3 - y3 * y3 - z3 * z3;
465
t3 * t3 * grad3(perm[ii + 1 + perm[jj + 1 + perm[kk + 1]]], x3, y3,
469
/* Add contributions from each corner to get the final noise value.
470
* The result is scaled to stay just inside [-1,1]
472
return 32.0f * (n0 + n1 + n2 + n3); /* TODO: The scale factor is preliminary! */
476
/** 4D simplex noise */
478
_mesa_noise4(GLfloat x, GLfloat y, GLfloat z, GLfloat w)
480
/* The skewing and unskewing factors are hairy again for the 4D case */
481
#define F4 0.309016994f /* F4 = (Math.sqrt(5.0)-1.0)/4.0 */
482
#define G4 0.138196601f /* G4 = (5.0-Math.sqrt(5.0))/20.0 */
484
float n0, n1, n2, n3, n4; /* Noise contributions from the five corners */
486
/* Skew the (x,y,z,w) space to determine which cell of 24 simplices we're in */
487
float s = (x + y + z + w) * F4; /* Factor for 4D skewing */
492
int i = FASTFLOOR(xs);
493
int j = FASTFLOOR(ys);
494
int k = FASTFLOOR(zs);
495
int l = FASTFLOOR(ws);
497
float t = (i + j + k + l) * G4; /* Factor for 4D unskewing */
498
float X0 = i - t; /* Unskew the cell origin back to (x,y,z,w) space */
503
float x0 = x - X0; /* The x,y,z,w distances from the cell origin */
508
/* For the 4D case, the simplex is a 4D shape I won't even try to describe.
509
* To find out which of the 24 possible simplices we're in, we need to
510
* determine the magnitude ordering of x0, y0, z0 and w0.
511
* The method below is a good way of finding the ordering of x,y,z,w and
512
* then find the correct traversal order for the simplex we're in.
513
* First, six pair-wise comparisons are performed between each possible pair
514
* of the four coordinates, and the results are used to add up binary bits
515
* for an integer index.
517
int c1 = (x0 > y0) ? 32 : 0;
518
int c2 = (x0 > z0) ? 16 : 0;
519
int c3 = (y0 > z0) ? 8 : 0;
520
int c4 = (x0 > w0) ? 4 : 0;
521
int c5 = (y0 > w0) ? 2 : 0;
522
int c6 = (z0 > w0) ? 1 : 0;
523
int c = c1 + c2 + c3 + c4 + c5 + c6;
525
int i1, j1, k1, l1; /* The integer offsets for the second simplex corner */
526
int i2, j2, k2, l2; /* The integer offsets for the third simplex corner */
527
int i3, j3, k3, l3; /* The integer offsets for the fourth simplex corner */
529
float x1, y1, z1, w1, x2, y2, z2, w2, x3, y3, z3, w3, x4, y4, z4, w4;
531
float t0, t1, t2, t3, t4;
534
* simplex[c] is a 4-vector with the numbers 0, 1, 2 and 3 in some
535
* order. Many values of c will never occur, since e.g. x>y>z>w
536
* makes x<z, y<w and x<w impossible. Only the 24 indices which
537
* have non-zero entries make any sense. We use a thresholding to
538
* set the coordinates in turn from the largest magnitude. The
539
* number 3 in the "simplex" array is at the position of the
540
* largest coordinate.
542
i1 = simplex[c][0] >= 3 ? 1 : 0;
543
j1 = simplex[c][1] >= 3 ? 1 : 0;
544
k1 = simplex[c][2] >= 3 ? 1 : 0;
545
l1 = simplex[c][3] >= 3 ? 1 : 0;
546
/* The number 2 in the "simplex" array is at the second largest coordinate. */
547
i2 = simplex[c][0] >= 2 ? 1 : 0;
548
j2 = simplex[c][1] >= 2 ? 1 : 0;
549
k2 = simplex[c][2] >= 2 ? 1 : 0;
550
l2 = simplex[c][3] >= 2 ? 1 : 0;
551
/* The number 1 in the "simplex" array is at the second smallest coordinate. */
552
i3 = simplex[c][0] >= 1 ? 1 : 0;
553
j3 = simplex[c][1] >= 1 ? 1 : 0;
554
k3 = simplex[c][2] >= 1 ? 1 : 0;
555
l3 = simplex[c][3] >= 1 ? 1 : 0;
556
/* The fifth corner has all coordinate offsets = 1, so no need to look that up. */
558
x1 = x0 - i1 + G4; /* Offsets for second corner in (x,y,z,w) coords */
562
x2 = x0 - i2 + 2.0f * G4; /* Offsets for third corner in (x,y,z,w) coords */
563
y2 = y0 - j2 + 2.0f * G4;
564
z2 = z0 - k2 + 2.0f * G4;
565
w2 = w0 - l2 + 2.0f * G4;
566
x3 = x0 - i3 + 3.0f * G4; /* Offsets for fourth corner in (x,y,z,w) coords */
567
y3 = y0 - j3 + 3.0f * G4;
568
z3 = z0 - k3 + 3.0f * G4;
569
w3 = w0 - l3 + 3.0f * G4;
570
x4 = x0 - 1.0f + 4.0f * G4; /* Offsets for last corner in (x,y,z,w) coords */
571
y4 = y0 - 1.0f + 4.0f * G4;
572
z4 = z0 - 1.0f + 4.0f * G4;
573
w4 = w0 - 1.0f + 4.0f * G4;
575
/* Wrap the integer indices at 256, to avoid indexing perm[] out of bounds */
581
/* Calculate the contribution from the five corners */
582
t0 = 0.6f - x0 * x0 - y0 * y0 - z0 * z0 - w0 * w0;
588
t0 * t0 * grad4(perm[ii + perm[jj + perm[kk + perm[ll]]]], x0, y0,
592
t1 = 0.6f - x1 * x1 - y1 * y1 - z1 * z1 - w1 * w1;
599
grad4(perm[ii + i1 + perm[jj + j1 + perm[kk + k1 + perm[ll + l1]]]],
603
t2 = 0.6f - x2 * x2 - y2 * y2 - z2 * z2 - w2 * w2;
610
grad4(perm[ii + i2 + perm[jj + j2 + perm[kk + k2 + perm[ll + l2]]]],
614
t3 = 0.6f - x3 * x3 - y3 * y3 - z3 * z3 - w3 * w3;
621
grad4(perm[ii + i3 + perm[jj + j3 + perm[kk + k3 + perm[ll + l3]]]],
625
t4 = 0.6f - x4 * x4 - y4 * y4 - z4 * z4 - w4 * w4;
632
grad4(perm[ii + 1 + perm[jj + 1 + perm[kk + 1 + perm[ll + 1]]]], x4,
636
/* Sum up and scale the result to cover the range [-1,1] */
637
return 27.0f * (n0 + n1 + n2 + n3 + n4); /* TODO: The scale factor is preliminary! */