Vuo  2.0.2
VuoGradientNoiseCommon.c
Go to the documentation of this file.
1 
10 #include "VuoGradientNoiseCommon.h"
11 #include <math.h>
12 
13 
14 #include "module.h"
15 
16 #ifdef VUO_COMPILER
18  "title" : "VuoGradientNoiseCommon",
19  "dependencies" : [
20  "VuoReal",
21  "VuoPoint2d",
22  "VuoPoint3d",
23  "VuoPoint4d"
24  ]
25  });
26 #endif
27 
28 
32 static unsigned char perm[] =
33 {
34  151, 160, 137, 91,
35  90, 15, 131, 13, 201, 95, 96, 53, 194, 233, 7, 225, 140, 36, 103,
36  30, 69, 142, 8, 99, 37, 240, 21, 10, 23, 190, 6, 148, 247, 120,
37  234, 75, 0, 26, 197, 62, 94, 252, 219, 203, 117, 35, 11, 32, 57,
38  177, 33, 88, 237, 149, 56, 87, 174, 20, 125, 136, 171, 168, 68,
39  175, 74, 165, 71, 134, 139, 48, 27, 166, 77, 146, 158, 231, 83,
40  111, 229, 122, 60, 211, 133, 230, 220, 105, 92, 41, 55, 46, 245,
41  40, 244, 102, 143, 54, 65, 25, 63, 161, 1, 216, 80, 73, 209, 76,
42  132, 187, 208, 89, 18, 169, 200, 196, 135, 130, 116, 188, 159, 86,
43  164, 100, 109, 198, 173, 186, 3, 64, 52, 217, 226, 250, 124, 123,
44  5, 202, 38, 147, 118, 126, 255, 82, 85, 212, 207, 206, 59, 227, 47,
45  16, 58, 17, 182, 189, 28, 42, 223, 183, 170, 213, 119, 248, 152, 2,
46  44, 154, 163, 70, 221, 153, 101, 155, 167, 43, 172, 9, 129, 22, 39,
47  253, 19, 98, 108, 110, 79, 113, 224, 232, 178, 185, 112, 104, 218,
48  246, 97, 228, 251, 34, 242, 193, 238, 210, 144, 12, 191, 179, 162,
49  241, 81, 51, 145, 235, 249, 14, 239, 107, 49, 192, 214, 31, 181,
50  199, 106, 157, 184, 84, 204, 176, 115, 121, 50, 45, 127, 4, 150,
51  254, 138, 236, 205, 93, 222, 114, 67, 29, 24, 72, 243, 141, 128,
52  195, 78, 66, 215, 61, 156, 180,
53  151, 160, 137, 91,
54  90, 15, 131, 13, 201, 95, 96, 53, 194, 233, 7, 225, 140, 36, 103,
55  30, 69, 142, 8, 99, 37, 240, 21, 10, 23, 190, 6, 148, 247, 120,
56  234, 75, 0, 26, 197, 62, 94, 252, 219, 203, 117, 35, 11, 32, 57,
57  177, 33, 88, 237, 149, 56, 87, 174, 20, 125, 136, 171, 168, 68,
58  175, 74, 165, 71, 134, 139, 48, 27, 166, 77, 146, 158, 231, 83,
59  111, 229, 122, 60, 211, 133, 230, 220, 105, 92, 41, 55, 46, 245,
60  40, 244, 102, 143, 54, 65, 25, 63, 161, 1, 216, 80, 73, 209, 76,
61  132, 187, 208, 89, 18, 169, 200, 196, 135, 130, 116, 188, 159, 86,
62  164, 100, 109, 198, 173, 186, 3, 64, 52, 217, 226, 250, 124, 123,
63  5, 202, 38, 147, 118, 126, 255, 82, 85, 212, 207, 206, 59, 227, 47,
64  16, 58, 17, 182, 189, 28, 42, 223, 183, 170, 213, 119, 248, 152, 2,
65  44, 154, 163, 70, 221, 153, 101, 155, 167, 43, 172, 9, 129, 22, 39,
66  253, 19, 98, 108, 110, 79, 113, 224, 232, 178, 185, 112, 104, 218,
67  246, 97, 228, 251, 34, 242, 193, 238, 210, 144, 12, 191, 179, 162,
68  241, 81, 51, 145, 235, 249, 14, 239, 107, 49, 192, 214, 31, 181,
69  199, 106, 157, 184, 84, 204, 176, 115, 121, 50, 45, 127, 4, 150,
70  254, 138, 236, 205, 93, 222, 114, 67, 29, 24, 72, 243, 141, 128,
71  195, 78, 66, 215, 61, 156, 180
72 };
73 
77 static unsigned char simplex[64][4] =
78 {
79  {0,1,2,3},{0,1,3,2},{0,0,0,0},{0,2,3,1},{0,0,0,0},{0,0,0,0},{0,0,0,0},{1,2,3,0},
80  {0,2,1,3},{0,0,0,0},{0,3,1,2},{0,3,2,1},{0,0,0,0},{0,0,0,0},{0,0,0,0},{1,3,2,0},
81  {0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},
82  {1,2,0,3},{0,0,0,0},{1,3,0,2},{0,0,0,0},{0,0,0,0},{0,0,0,0},{2,3,0,1},{2,3,1,0},
83  {1,0,2,3},{1,0,3,2},{0,0,0,0},{0,0,0,0},{0,0,0,0},{2,0,3,1},{0,0,0,0},{2,1,3,0},
84  {0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},{0,0,0,0},
85  {2,0,1,3},{0,0,0,0},{0,0,0,0},{0,0,0,0},{3,0,1,2},{3,0,2,1},{0,0,0,0},{3,1,2,0},
86  {2,1,0,3},{0,0,0,0},{0,0,0,0},{0,0,0,0},{3,1,0,2},{0,0,0,0},{3,2,0,1},{3,2,1,0}
87 };
88 
89 
90 static inline VuoReal grad1d(int hash, VuoReal x);
91 static inline VuoReal grad2dPerlin(int hash, VuoReal x, VuoReal y);
92 static inline VuoReal grad2dSimplex(int hash, VuoReal x, VuoReal y);
93 static inline VuoReal grad3d(int hash, VuoReal x, VuoReal y, VuoReal z);
94 static inline VuoReal grad4dPerlin(int hash, VuoReal x, VuoReal y, VuoReal z, VuoReal w);
95 static inline float grad4dSimplex(int hash, float x, float y, float z, float t);
96 
97 static inline VuoReal fade(VuoReal t);
98 static inline VuoReal lerp(VuoReal t, VuoReal a, VuoReal b);
99 
100 
105 {
106  VuoInteger X;
107  VuoReal u;
108  VuoInteger A, B;
109 
110  X = (int)floor(x) & 255; // Integer part of x
111  x = x - floor(x); // Fractional part of x
112  u = fade(x);
113 
114  A = perm[X ];
115  B = perm[X+1];
116 
117  return 0.25 * lerp(u,
118  grad1d(perm[A], x),
119  grad1d(perm[B], x-1));
120 }
121 
126 {
129 }
130 
135 {
139 }
140 
145 {
150 }
151 
156 {
157  VuoReal x = point.x;
158  VuoReal y = point.y;
159 
160  VuoInteger X, Y;
161  VuoReal u, v;
162  VuoInteger A, AA, AB, B, BA, BB;
163 
164  X = (int)floor(x) & 255;
165  Y = (int)floor(y) & 255;
166 
167  x = x - floor(x);
168  y = y - floor(y);
169 
170  u = fade(x);
171  v = fade(y);
172 
173  A = perm[X ] + Y;
174  AA = perm[A ] ;
175  AB = perm[A+1] ;
176  B = perm[X+1] + Y;
177  BA = perm[B ] ;
178  BB = perm[B+1] ;
179 
180  return lerp(v,
181  lerp(u,
182  grad2dPerlin(perm[AA], x , y ),
183  grad2dPerlin(perm[BA], x-1, y )),
184  lerp(u,
185  grad2dPerlin(perm[AB], x , y-1),
186  grad2dPerlin(perm[BB], x-1, y-1)));
187 }
188 
193 {
196 }
197 
202 {
206 }
207 
212 {
217 }
218 
223 {
224  VuoReal x = point.x;
225  VuoReal y = point.y;
226  VuoReal z = point.z;
227 
228  VuoInteger X, Y, Z;
229  VuoReal u, v, w;
230  VuoInteger A, AA, AB, B, BA, BB;
231 
232  X = (int)floor(x) & 255; // Integer part of x
233  Y = (int)floor(y) & 255; // Integer part of y
234  Z = (int)floor(z) & 255; // Integer part of z
235 
236  x = x - floor(x); // Fractional part of x
237  y = y - floor(y); // Fractional part of y
238  z = z - floor(z); // Fractional part of z
239  u = fade(x);
240  v = fade(y);
241  w = fade(z);
242 
243  A = perm[X ] + Y;
244  AA = perm[A ] + Z;
245  AB = perm[A+1] + Z;
246  B = perm[X+1] + Y;
247  BA = perm[B ] + Z;
248  BB = perm[B+1] + Z;
249 
250  return lerp(w,
251  lerp(v,
252  lerp(u,
253  grad3d(perm[AA ], x , y , z),
254  grad3d(perm[BA ], x-1, y , z)),
255  lerp(u,
256  grad3d(perm[AB ], x , y-1, z),
257  grad3d(perm[BB ], x-1, y-1, z))),
258  lerp(v,
259  lerp(u,
260  grad3d(perm[AA+1], x , y , z-1),
261  grad3d(perm[BA+1], x-1, y , z-1)),
262  lerp(u,
263  grad3d(perm[AB+1], x , y-1, z-1),
264  grad3d(perm[BB+1], x-1, y-1, z-1))));
265 }
266 
271 {
274 }
275 
280 {
284 }
285 
290 {
295 }
296 
301 {
302  VuoReal x = point.x;
303  VuoReal y = point.y;
304  VuoReal z = point.z;
305  VuoReal w = point.w;
306 
307  VuoInteger X, Y, Z, W;
308  VuoReal a, b, c, d;
309  VuoInteger A, B,
310  AA, AB, BA, BB,
311  AAA, AAB, ABA, ABB,
312  BAA, BAB, BBA, BBB;
313 
314  X = (int)floor(x) & 255;
315  Y = (int)floor(y) & 255;
316  Z = (int)floor(z) & 255;
317  W = (int)floor(w) & 255;
318 
319  x = x - floor(x);
320  y = y - floor(y);
321  z = z - floor(z);
322  w = w - floor(w);
323 
324  a = fade(x);
325  b = fade(y);
326  c = fade(z);
327  d = fade(w);
328 
329  A = perm[X ]+Y;
330  AA = perm[A ]+Z;
331  AB = perm[A +1]+Z;
332  B = perm[X +1]+Y;
333  BA = perm[B ]+Z;
334  BB = perm[B +1]+Z;
335  AAA = perm[AA ]+W;
336  AAB = perm[AA+1]+W;
337  ABA = perm[AB ]+W;
338  ABB = perm[AB+1]+W;
339  BAA = perm[BA ]+W;
340  BAB = perm[BA+1]+W;
341  BBA = perm[BB ]+W;
342  BBB = perm[BB+1]+W;
343 
344  return lerp(d,
345  lerp(c,
346  lerp(b,
347  lerp(a,
348  grad4dPerlin(perm[AAA ], x , y , z , w ),
349  grad4dPerlin(perm[BAA ], x-1, y , z , w )),
350  lerp(b,
351  grad4dPerlin(perm[ABA ], x , y , z , w ),
352  grad4dPerlin(perm[BBA ], x-1, y-1, z , w ))),
353  lerp(b,
354  lerp(a,
355  grad4dPerlin(perm[AAB ], x , y , z-1, w ),
356  grad4dPerlin(perm[BAB ], x-1, y , z-1, w )),
357  lerp(a,
358  grad4dPerlin(perm[ABB ], x , y-1, z-1, w ),
359  grad4dPerlin(perm[BBB ], x-1, y-1, z-1, w )))),
360  lerp(c,
361  lerp(b,
362  lerp(a,
363  grad4dPerlin(perm[AAA+1], x , y , z , w-1),
364  grad4dPerlin(perm[BAA+1], x-1, y , z , w-1)),
365  lerp(b,
366  grad4dPerlin(perm[ABA+1], x , y , z , w-1),
367  grad4dPerlin(perm[BBA+1], x-1, y-1, z , w-1))),
368  lerp(b,
369  lerp(a,
370  grad4dPerlin(perm[AAB+1], x , y , z-1, w-1),
371  grad4dPerlin(perm[BAB+1], x-1, y , z-1, w-1)),
372  lerp(a,
373  grad4dPerlin(perm[ABB+1], x , y-1, z-1, w-1),
374  grad4dPerlin(perm[BBB+1], x-1, y-1, z-1, w-1)))));
375 }
376 
381 {
384 }
385 
390 {
394 }
395 
400 {
405 }
406 
411 {
412  VuoInteger i0 = (int)floor(x);
413  VuoInteger i1 = i0 + 1;
414  VuoReal x0 = x - i0;
415  VuoReal x1 = x0 - 1.0f;
416 
417  VuoReal n0, n1;
418 
419  float t0 = 1.0f - x0 * x0;
420  t0 *= t0;
421  n0 = t0 * t0 * grad1d(perm[i0 & 0xff], x0);
422 
423  float t1 = 1.0f - x1 * x1;
424  t1 *= t1;
425  n1 = t1 * t1 * grad1d(perm[i1 & 0xff], x1);
426  return 0.395f * (n0 + n1);
427 }
428 
433 {
436 }
437 
442 {
446 }
447 
452 {
457 }
458 
463 {
464  VuoReal x = point.x;
465  VuoReal y = point.y;
466 
467  VuoReal f2 = 0.366025403; // 0.5 * (sqrt(3.0) - 1.0)
468  VuoReal g2 = 0.211324865; // (3.0 - sqrt(3.0)) / 6.0
469 
470  // Noise contributions from three corners
471  VuoReal n0, n1, n2;
472 
473  // Skew the input space to determine which simplex cell we're in
474  VuoReal s = (x + y) * f2; // Hairy factor for 2D
475  VuoReal xs = x + s;
476  VuoReal ys = y + s;
477 
478  VuoInteger i = (int)floor(xs);
479  VuoInteger j = (int)floor(ys);
480 
481  VuoReal t = (VuoReal) (i + j) * g2;
482  // Unskew the cell origin back to (x, y) space
483  VuoReal X0 = i - t;
484  VuoReal Y0 = j - t;
485  // The x, y distances from the cell origin
486  VuoReal x0 = x - X0;
487  VuoReal y0 = y - Y0;
488 
489  // For the 2D case, the simplex shape is an equilateral triangle.
490  // Determine which simplex we are in.
491 
492  // Offsets for second (middle) corner of simplex in (i, j) coords
493  VuoInteger i1, j1;
494  // Lower triangle, XY order: (0,0)->(1,0)->(1,1)
495  if (x0 > y0)
496  {
497  i1 = 1;
498  j1 = 0;
499  }
500  // Upper triangle, YX order: (0,0)->(0,1)->(1,1)
501  else
502  {
503  i1 = 0;
504  j1 = 1;
505  }
506 
507  // Offsets for middle corner in (x, y) unskewed coords
508  VuoReal x1 = x0 - i1 + g2;
509  VuoReal y1 = y0 - j1 + g2;
510  // Ofsets for last corner in (x, y) unskewed coords
511  VuoReal x2 = x0 - 1.0f + 2.0f * g2;
512  VuoReal y2 = y0 - 1.0f + 2.0f * g2;
513 
514  // Wrap the integer indices at 256, to avoid indexing perm[] out of bounds
515  VuoInteger ii = i & 0xff;
516  VuoInteger jj = j & 0xff;
517 
518  // Calculate the contribution from the three corners
519  VuoReal t0 = 0.5f - x0 * x0 - y0 * y0;
520  if (t0 < 0.0f)
521  n0 = 0.0f;
522  else
523  {
524  t0 *= t0;
525  n0 = t0 * t0 * grad2dSimplex(perm[ii + perm[jj]], x0, y0);
526  }
527 
528  VuoReal t1 = 0.5f - x1 * x1 - y1 * y1;
529  if (t1 < 0.0f)
530  n1 = 0.0f;
531  else
532  {
533  t1 *= t1;
534  n1 = t1 * t1 * grad2dSimplex(perm[ii + i1 + perm[jj + j1]], x1, y1);
535  }
536 
537  VuoReal t2 = 0.5f - x2 * x2 - y2 * y2;
538  if (t2 < 0.0f)
539  n2 = 0.0f;
540  else
541  {
542  t2 *= t2;
543  n2 = t2 * t2 * grad2dSimplex(perm[ii + 1 + perm[jj + 1]], x2, y2);
544  }
545 
546  // Add contributions from each corner to the final noise value.
547  // The result should be scaled to return values in the interval [-1, 1].
548  return 46.0f * (n0 + n1 + n2);
549 }
550 
555 {
558 }
559 
564 {
568 }
569 
574 {
579 }
580 
585 {
586  VuoReal x = point.x;
587  VuoReal y = point.y;
588  VuoReal z = point.z;
589 
590  // Simple skewing factors for the 3D case
591  VuoReal F3 = 0.333333333;
592  VuoReal G3 = 0.166666667;
593 
594  VuoReal n0, n1, n2, n3; // Noise contributions from the four corners
595 
596  // Skew the input space to determine which simplex cell we're in
597  VuoReal s = (x + y + z) * F3; // Very nice and simple skew factor for 3D
598  VuoReal xs = x + s;
599  VuoReal ys = y + s;
600  VuoReal zs = z + s;
601  VuoInteger i = (int)floor(xs);
602  VuoInteger j = (int)floor(ys);
603  VuoInteger k = (int)floor(zs);
604 
605  VuoReal t = (VuoReal)(i + j + k) * G3;
606  VuoReal X0 = i - t; // Unskew the cell origin back to (x,y,z) space
607  VuoReal Y0 = j - t;
608  VuoReal Z0 = k - t;
609  VuoReal x0 = x - X0; // The x,y,z distances from the cell origin
610  VuoReal y0 = y - Y0;
611  VuoReal z0 = z - Z0;
612 
613  // For the 3D case, the simplex shape is a slightly irregular tetrahedron.
614  // Determine which simplex we are in.
615  VuoInteger i1, j1, k1; // Offsets for second corner of simplex in (i,j,k) coords
616  VuoInteger i2, j2, k2; // Offsets for third corner of simplex in (i,j,k) coords
617 
618  if (x0 >= y0)
619  {
620  if (y0 >= z0)
621  {
622  i1 = 1;
623  j1 = 0;
624  k1 = 0;
625  i2 = 1;
626  j2 = 1;
627  k2 = 0;
628  } // X Y Z order
629  else if (x0 >= z0)
630  {
631  i1 = 1;
632  j1 = 0;
633  k1 = 0;
634  i2 = 1;
635  j2 = 0;
636  k2 = 1;
637  } // X Z Y order
638  else
639  {
640  i1 = 0;
641  j1 = 0;
642  k1 = 1;
643  i2 = 1;
644  j2 = 0;
645  k2 = 1;
646  } // Z X Y order
647  }
648  else
649  { // x0 < y0
650  if (y0 < z0)
651  {
652  i1 = 0;
653  j1 = 0;
654  k1 = 1;
655  i2 = 0;
656  j2 = 1;
657  k2 = 1;
658  } // Z Y X order
659  else if (x0 < z0)
660  {
661  i1 = 0;
662  j1 = 1;
663  k1 = 0;
664  i2 = 0;
665  j2 = 1;
666  k2 = 1;
667  } // Y Z X order
668  else {
669  i1 = 0;
670  j1 = 1;
671  k1 = 0;
672  i2 = 1;
673  j2 = 1;
674  k2 = 0;
675  } // Y X Z order
676  }
677 
678  // A step of (1,0,0) in (i,j,k) means a step of (1-c,-c,-c) in (x,y,z),
679  // a step of (0,1,0) in (i,j,k) means a step of (-c,1-c,-c) in (x,y,z), and
680  // a step of (0,0,1) in (i,j,k) means a step of (-c,-c,1-c) in (x,y,z), where
681  // c = 1/6.
682 
683  VuoReal x1 = x0 - i1 + G3; // Offsets for second corner in (x,y,z) coords
684  VuoReal y1 = y0 - j1 + G3;
685  VuoReal z1 = z0 - k1 + G3;
686  VuoReal x2 = x0 - i2 + 2.0f * G3; // Offsets for third corner in (x,y,z) coords
687  VuoReal y2 = y0 - j2 + 2.0f * G3;
688  VuoReal z2 = z0 - k2 + 2.0f * G3;
689  VuoReal x3 = x0 - 1.0f + 3.0f * G3; // Offsets for last corner in (x,y,z) coords
690  VuoReal y3 = y0 - 1.0f + 3.0f * G3;
691  VuoReal z3 = z0 - 1.0f + 3.0f * G3;
692 
693  // Wrap the VuoIntegereger indices at 256, to avoid indexing perm[] out of bounds
694  VuoInteger ii = i & 0xff;
695  VuoInteger jj = j & 0xff;
696  VuoInteger kk = k & 0xff;
697 
698  // Calculate the contribution from the four corners
699  VuoReal t0 = 0.6f - x0 * x0 - y0 * y0 - z0 * z0;
700  if (t0 < 0.0f)
701  n0 = 0.0f;
702  else
703  {
704  t0 *= t0;
705  n0 = t0 * t0 * grad3d(perm[ii+perm[jj+perm[kk]]], x0, y0, z0);
706  }
707 
708  VuoReal t1 = 0.6f - x1 * x1 - y1 * y1 - z1 * z1;
709  if (t1 < 0.0f)
710  n1 = 0.0f;
711  else
712  {
713  t1 *= t1;
714  n1 = t1 * t1 * grad3d(perm[ii+i1+perm[jj+j1+perm[kk+k1]]], x1, y1, z1);
715  }
716 
717  VuoReal t2 = 0.6f - x2 * x2 - y2 * y2 - z2 * z2;
718  if (t2 < 0.0f)
719  n2 = 0.0f;
720  else
721  {
722  t2 *= t2;
723  n2 = t2 * t2 * grad3d(perm[ii+i2+perm[jj+j2+perm[kk+k2]]], x2, y2, z2);
724  }
725 
726  VuoReal t3 = 0.6f - x3 * x3 - y3 * y3 - z3 * z3;
727  if (t3 < 0.0f)
728  n3 = 0.0f;
729  else
730  {
731  t3 *= t3;
732  n3 = t3 * t3 * grad3d(perm[ii+1+perm[jj+1+perm[kk+1]]], x3, y3, z3);
733  }
734 
735  // Add contributions from each corner to get the final noise value.
736  // The result is scaled to stay just inside [-1,1]
737  return 32.0f * (n0 + n1 + n2 + n3); // TODO: The scale factor is preliminary!
738 }
739 
744 {
747 }
748 
753 {
757 }
758 
763 {
768 }
769 
774 {
775  float x = point.x;
776  float y = point.y;
777  float z = point.z;
778  float w = point.w;
779 
780  // The skewing and unskewing factors are hairy again for the 4D case
781  float F4 = 0.309016994f; // F4 = (Math.sqrt(5.0)-1.0)/4.0
782  float G4 = 0.138196601f; // G4 = (5.0-Math.sqrt(5.0))/20.0
783 
784  float n0, n1, n2, n3, n4; // Noise contributions from the five corners
785 
786  // Skew the (x,y,z,w) space to determine which cell of 24 simplices we're in
787  float s = (x + y + z + w) * F4; // Factor for 4D skewing
788  float xs = x + s;
789  float ys = y + s;
790  float zs = z + s;
791  float ws = w + s;
792  int i = (int)floorf(xs);
793  int j = (int)floorf(ys);
794  int k = (int)floorf(zs);
795  int l = (int)floorf(ws);
796 
797  float t = (i + j + k + l) * G4; // Factor for 4D unskewing
798  float X0 = i - t; // Unskew the cell origin back to (x,y,z,w) space
799  float Y0 = j - t;
800  float Z0 = k - t;
801  float W0 = l - t;
802 
803  float x0 = x - X0; // The x,y,z,w distances from the cell origin
804  float y0 = y - Y0;
805  float z0 = z - Z0;
806  float w0 = w - W0;
807 
808  // For the 4D case, the simplex is a 4D shape I won't even try to describe.
809  // To find out which of the 24 possible simplices we're in, we need to
810  // determine the magnitude ordering of x0, y0, z0 and w0.
811  // The method below is a good way of finding the ordering of x,y,z,w and
812  // then find the correct traversal order for the simplex we’re in.
813  // First, six pair-wise comparisons are performed between each possible pair
814  // of the four coordinates, and the results are used to add up binary bits
815  // for an integer index.
816  int c1 = (x0 > y0) ? 32 : 0;
817  int c2 = (x0 > z0) ? 16 : 0;
818  int c3 = (y0 > z0) ? 8 : 0;
819  int c4 = (x0 > w0) ? 4 : 0;
820  int c5 = (y0 > w0) ? 2 : 0;
821  int c6 = (z0 > w0) ? 1 : 0;
822  int c = c1 + c2 + c3 + c4 + c5 + c6;
823 
824  int i1, j1, k1, l1; // The integer offsets for the second simplex corner
825  int i2, j2, k2, l2; // The integer offsets for the third simplex corner
826  int i3, j3, k3, l3; // The integer offsets for the fourth simplex corner
827 
828  // simplex[c] is a 4-vector with the numbers 0, 1, 2 and 3 in some order.
829  // Many values of c will never occur, since e.g. x>y>z>w makes x<z, y<w and x<w
830  // impossible. Only the 24 indices which have non-zero entries make any sense.
831  // We use a thresholding to set the coordinates in turn from the largest magnitude.
832  // The number 3 in the "simplex" array is at the position of the largest coordinate.
833  i1 = simplex[c][0] >= 3 ? 1 : 0;
834  j1 = simplex[c][1] >= 3 ? 1 : 0;
835  k1 = simplex[c][2] >= 3 ? 1 : 0;
836  l1 = simplex[c][3] >= 3 ? 1 : 0;
837  // The number 2 in the "simplex" array is at the second largest coordinate.
838  i2 = simplex[c][0] >= 2 ? 1 : 0;
839  j2 = simplex[c][1] >= 2 ? 1 : 0;
840  k2 = simplex[c][2] >= 2 ? 1 : 0;
841  l2 = simplex[c][3] >= 2 ? 1 : 0;
842  // The number 1 in the "simplex" array is at the second smallest coordinate.
843  i3 = simplex[c][0]>=1 ? 1 : 0;
844  j3 = simplex[c][1]>=1 ? 1 : 0;
845  k3 = simplex[c][2]>=1 ? 1 : 0;
846  l3 = simplex[c][3]>=1 ? 1 : 0;
847  // The fifth corner has all coordinate offsets = 1, so no need to look that up.
848 
849  float x1 = x0 - i1 + G4; // Offsets for second corner in (x,y,z,w) coords
850  float y1 = y0 - j1 + G4;
851  float z1 = z0 - k1 + G4;
852  float w1 = w0 - l1 + G4;
853  float x2 = x0 - i2 + 2.0f*G4; // Offsets for third corner in (x,y,z,w) coords
854  float y2 = y0 - j2 + 2.0f*G4;
855  float z2 = z0 - k2 + 2.0f*G4;
856  float w2 = w0 - l2 + 2.0f*G4;
857  float x3 = x0 - i3 + 3.0f*G4; // Offsets for fourth corner in (x,y,z,w) coords
858  float y3 = y0 - j3 + 3.0f*G4;
859  float z3 = z0 - k3 + 3.0f*G4;
860  float w3 = w0 - l3 + 3.0f*G4;
861  float x4 = x0 - 1.0f + 4.0f*G4; // Offsets for last corner in (x,y,z,w) coords
862  float y4 = y0 - 1.0f + 4.0f*G4;
863  float z4 = z0 - 1.0f + 4.0f*G4;
864  float w4 = w0 - 1.0f + 4.0f*G4;
865 
866  // Wrap the integer indices at 256, to avoid indexing perm[] out of bounds
867  int ii = i & 0xff;
868  int jj = j & 0xff;
869  int kk = k & 0xff;
870  int ll = l & 0xff;
871 
872  // Calculate the contribution from the five corners
873  float t0 = 0.6f - x0*x0 - y0*y0 - z0*z0 - w0*w0;
874  if (t0 < 0.0f)
875  n0 = 0.0f;
876  else
877  {
878  t0 *= t0;
879  n0 = t0 * t0 * grad4dSimplex(perm[ii + perm[jj + perm[kk + perm[ll]]]], x0, y0, z0, w0);
880  }
881 
882  float t1 = 0.6f - x1*x1 - y1*y1 - z1*z1 - w1*w1;
883  if(t1 < 0.0f)
884  n1 = 0.0f;
885  else
886  {
887  t1 *= t1;
888  n1 = t1 * t1 * grad4dSimplex(perm[ii + i1 + perm[jj + j1 + perm[kk + k1 + perm[ll + l1]]]], x1, y1, z1, w1);
889  }
890 
891  float t2 = 0.6f - x2*x2 - y2*y2 - z2*z2 - w2*w2;
892  if (t2 < 0.0f)
893  n2 = 0.0f;
894  else
895  {
896  t2 *= t2;
897  n2 = t2 * t2 * grad4dSimplex(perm[ii + i2 + perm[jj + j2 + perm[kk + k2 + perm[ll + l2]]]], x2, y2, z2, w2);
898  }
899 
900  float t3 = 0.6f - x3*x3 - y3*y3 - z3*z3 - w3*w3;
901  if (t3 < 0.0f)
902  n3 = 0.0f;
903  else
904  {
905  t3 *= t3;
906  n3 = t3 * t3 * grad4dSimplex(perm[ii + i3 + perm[jj + j3 + perm[kk + k3 + perm[ll + l3]]]], x3, y3, z3, w3);
907  }
908 
909  float t4 = 0.6f - x4*x4 - y4*y4 - z4*z4 - w4*w4;
910  if (t4 < 0.0f)
911  n4 = 0.0f;
912  else
913  {
914  t4 *= t4;
915  n4 = t4 * t4 * grad4dSimplex(perm[ii + 1 + perm[jj + 1 + perm[kk + 1 + perm[ll + 1]]]], x4, y4, z4, w4);
916  }
917 
918  // Sum up and scale the result to cover the range [-1,1]
919  return 27.0f * (n0 + n1 + n2 + n3 + n4); // TODO: The scale factor is preliminary!
920 }
921 
926 {
929 }
930 
935 {
939 }
940 
945 {
950 }
951 
955 static inline VuoReal grad1d(int hash, VuoReal x)
956 {
957  int h = hash & 15;
958 
959  /*
960  * Gradient value 1.0, 2.0, ..., 8.0
961  * and a random sign for the gradient
962  */
963  float grad = 1.0 + (h & 7);
964  if (h & 8)
965  {
966  grad = -grad;
967  }
968 
969  // Multiply the gradient with the distance
970  return grad * x;
971 }
972 
976 static inline VuoReal grad2dPerlin(int hash, VuoReal x, VuoReal y)
977 {
978  VuoInteger h = hash & 15;
979  VuoReal u = h < 8 ? x : y;
980  VuoReal v = h < 4 ? y : h==12 || h==14 ? x : 0;
981  return ((h&1) == 0 ? u : -u) + ((h&2) == 0 ? v : -v);
982 }
983 
987 static inline VuoReal grad2dSimplex(int hash, VuoReal x, VuoReal y)
988 {
989  VuoInteger h = hash & 7;
990  VuoReal u = h < 4 ? x : y;
991  VuoReal v = h < 4 ? y : x;
992  return ((h&1) ? -u : u) + ((h&2) ? -2.0f*v : 2.0f*v);
993 }
994 
998 static inline VuoReal grad3d(int hash, VuoReal x, VuoReal y, VuoReal z)
999 {
1000  VuoInteger h = hash & 15;
1001  VuoReal u = h<8 ? x : y;
1002  VuoReal v = h<4 ? y : h==12||h==14 ? x : z;
1003  return ((h&1) == 0 ? u : - u) + ((h&2) == 0 ? v : -v);
1004 }
1005 
1009 static inline VuoReal grad4dPerlin(int hash, VuoReal x, VuoReal y, VuoReal z, VuoReal w)
1010 {
1011  VuoInteger h = hash & 31; // CONVERT LO 5 BITS OF HASH TO 32 GRAD DIRECTIONS.
1012  VuoReal a=y,b=z,c=w; // X,Y,Z
1013  switch (h >> 3) { // OR, DEPENDING ON HIGH ORDER 2 BITS:
1014  case 1:
1015  a=w;
1016  b=x;
1017  c=y;
1018  break; // W,X,Y
1019  case 2:
1020  a=z;
1021  b=w;
1022  c=x;
1023  break; // Z,W,X
1024  case 3:
1025  a=y;
1026  b=z;
1027  c=w;
1028  break; // Y,Z,W
1029  }
1030 
1031  return ((h&4)==0 ? -a:a) + ((h&2)==0 ? -b:b) + ((h&1)==0 ? -c:c);
1032 }
1033 
1037 static inline float grad4dSimplex(int hash, float x, float y, float z, float t)
1038 {
1039  int h = hash & 31;
1040  float u = h < 24 ? x : y;
1041  float v = h < 16 ? y : z;
1042  float w = h < 8 ? z : t;
1043  return ((h&1) ? -u : u) + ((h&2) ? -v : v) + ((h&4) ? -w : w);
1044 }
1045 
1049 static inline VuoReal fade(VuoReal t)
1050 {
1051  return (t * t * t * (t * (t * 6 - 15) + 10));
1052 }
1053 
1057 static inline VuoReal lerp(VuoReal t, VuoReal a, VuoReal b)
1058 {
1059  return (a + t * (b - a));
1060 }
1061 
1062