Blender V4.3
noise.cc
Go to the documentation of this file.
1/* SPDX-FileCopyrightText: 2009-2010 Sony Pictures Imageworks Inc., et al.
2 * All Rights Reserved. (BSD-3-Clause).
3 * SPDX-FileCopyrightText: 2011 Blender Authors (GPL-2.0-or-later).
4 *
5 * SPDX-License-Identifier: GPL-2.0-or-later AND BSD-3-Clause */
6
7#include <algorithm>
8#include <cmath>
9#include <cstdint>
10
11#include "BLI_math_base.hh"
12#include "BLI_math_base_safe.h"
14#include "BLI_math_numbers.hh"
15#include "BLI_math_vector.hh"
16#include "BLI_noise.hh"
17#include "BLI_utildefines.h"
18
19namespace blender::noise {
20
21/* -------------------------------------------------------------------- */
28{
29 return (x << k) | (x >> (32 - k));
30}
31
33{
34 a -= c;
35 a ^= hash_bit_rotate(c, 4);
36 c += b;
37 b -= a;
38 b ^= hash_bit_rotate(a, 6);
39 a += c;
40 c -= b;
41 c ^= hash_bit_rotate(b, 8);
42 b += a;
43 a -= c;
44 a ^= hash_bit_rotate(c, 16);
45 c += b;
46 b -= a;
47 b ^= hash_bit_rotate(a, 19);
48 a += c;
49 c -= b;
50 c ^= hash_bit_rotate(b, 4);
51 b += a;
52}
53
55{
56 c ^= b;
57 c -= hash_bit_rotate(b, 14);
58 a ^= c;
59 a -= hash_bit_rotate(c, 11);
60 b ^= a;
61 b -= hash_bit_rotate(a, 25);
62 c ^= b;
63 c -= hash_bit_rotate(b, 16);
64 a ^= c;
65 a -= hash_bit_rotate(c, 4);
66 b ^= a;
67 b -= hash_bit_rotate(a, 14);
68 c ^= b;
69 c -= hash_bit_rotate(b, 24);
70}
71
73{
74 uint32_t a, b, c;
75 a = b = c = 0xdeadbeef + (1 << 2) + 13;
76
77 a += kx;
78 hash_bit_final(a, b, c);
79
80 return c;
81}
82
84{
85 uint32_t a, b, c;
86 a = b = c = 0xdeadbeef + (2 << 2) + 13;
87
88 b += ky;
89 a += kx;
90 hash_bit_final(a, b, c);
91
92 return c;
93}
94
96{
97 uint32_t a, b, c;
98 a = b = c = 0xdeadbeef + (3 << 2) + 13;
99
100 c += kz;
101 b += ky;
102 a += kx;
103 hash_bit_final(a, b, c);
104
105 return c;
106}
107
109{
110 uint32_t a, b, c;
111 a = b = c = 0xdeadbeef + (4 << 2) + 13;
112
113 a += kx;
114 b += ky;
115 c += kz;
116 hash_bit_mix(a, b, c);
117
118 a += kw;
119 hash_bit_final(a, b, c);
120
121 return c;
122}
123
125{
126 union {
127 uint32_t i;
128 float f;
129 } u;
130 u.f = f;
131 return u.i;
132}
133
135{
136 return hash(float_as_uint(kx));
137}
138
140{
141 return hash(float_as_uint(k.x), float_as_uint(k.y));
142}
143
145{
146 return hash(float_as_uint(k.x), float_as_uint(k.y), float_as_uint(k.z));
147}
148
150{
151 return hash(float_as_uint(k.x), float_as_uint(k.y), float_as_uint(k.z), float_as_uint(k.w));
152}
153
155{
156 return hash(hash_float(k.x), hash_float(k.y), hash_float(k.z), hash_float(k.w));
157}
158
159/* Hashing a number of uint32_t into a float in the range [0, 1]. */
160
162{
163 return float(k) / float(0xFFFFFFFFu);
164}
165
167{
168 return uint_to_float_01(hash(kx));
169}
170
172{
173 return uint_to_float_01(hash(kx, ky));
174}
175
177{
178 return uint_to_float_01(hash(kx, ky, kz));
179}
180
182{
183 return uint_to_float_01(hash(kx, ky, kz, kw));
184}
185
186/* Hashing a number of floats into a float in the range [0, 1]. */
187
189{
190 return uint_to_float_01(hash_float(k));
191}
192
194{
195 return uint_to_float_01(hash_float(k));
196}
197
199{
200 return uint_to_float_01(hash_float(k));
201}
202
204{
205 return uint_to_float_01(hash_float(k));
206}
207
212
214{
215 return float2(hash_float_to_float(float3(k.x, k.y, k.z)),
216 hash_float_to_float(float3(k.z, k.x, k.y)));
217}
218
220{
221 return float2(hash_float_to_float(float4(k.x, k.y, k.z, k.w)),
222 hash_float_to_float(float4(k.z, k.x, k.w, k.y)));
223}
224
226{
227 return float3(hash_float_to_float(k),
229 hash_float_to_float(float2(k, 2.0)));
230}
231
233{
234 return float3(hash_float_to_float(k),
235 hash_float_to_float(float3(k.x, k.y, 1.0)),
236 hash_float_to_float(float3(k.x, k.y, 2.0)));
237}
238
240{
241 return float3(hash_float_to_float(k),
242 hash_float_to_float(float4(k.x, k.y, k.z, 1.0)),
243 hash_float_to_float(float4(k.x, k.y, k.z, 2.0)));
244}
245
247{
248 return float3(hash_float_to_float(k),
249 hash_float_to_float(float4(k.z, k.x, k.w, k.y)),
250 hash_float_to_float(float4(k.w, k.z, k.y, k.x)));
251}
252
254{
255 return float4(hash_float_to_float(k),
256 hash_float_to_float(float4(k.w, k.x, k.y, k.z)),
257 hash_float_to_float(float4(k.z, k.w, k.x, k.y)),
258 hash_float_to_float(float4(k.y, k.z, k.w, k.x)));
259}
260
263/* -------------------------------------------------------------------- */
273/* Linear Interpolation. */
274template<typename T> T static mix(T v0, T v1, float x)
275{
276 return (1 - x) * v0 + x * v1;
277}
278
279/* Bilinear Interpolation:
280 *
281 * v2 v3
282 * @ + + + + @ y
283 * + + ^
284 * + + |
285 * + + |
286 * @ + + + + @ @------> x
287 * v0 v1
288 */
289BLI_INLINE float mix(float v0, float v1, float v2, float v3, float x, float y)
290{
291 float x1 = 1.0 - x;
292 return (1.0 - y) * (v0 * x1 + v1 * x) + y * (v2 * x1 + v3 * x);
293}
294
295/* Trilinear Interpolation:
296 *
297 * v6 v7
298 * @ + + + + + + @
299 * +\ +\
300 * + \ + \
301 * + \ + \
302 * + \ v4 + \ v5
303 * + @ + + + +++ + @ z
304 * + + + + y ^
305 * v2 @ + +++ + + + @ v3 + \ |
306 * \ + \ + \ |
307 * \ + \ + \|
308 * \ + \ + +---------> x
309 * \+ \+
310 * @ + + + + + + @
311 * v0 v1
312 */
313BLI_INLINE float mix(float v0,
314 float v1,
315 float v2,
316 float v3,
317 float v4,
318 float v5,
319 float v6,
320 float v7,
321 float x,
322 float y,
323 float z)
324{
325 float x1 = 1.0 - x;
326 float y1 = 1.0 - y;
327 float z1 = 1.0 - z;
328 return z1 * (y1 * (v0 * x1 + v1 * x) + y * (v2 * x1 + v3 * x)) +
329 z * (y1 * (v4 * x1 + v5 * x) + y * (v6 * x1 + v7 * x));
330}
331
332/* Quadrilinear Interpolation. */
333BLI_INLINE float mix(float v0,
334 float v1,
335 float v2,
336 float v3,
337 float v4,
338 float v5,
339 float v6,
340 float v7,
341 float v8,
342 float v9,
343 float v10,
344 float v11,
345 float v12,
346 float v13,
347 float v14,
348 float v15,
349 float x,
350 float y,
351 float z,
352 float w)
353{
354 return mix(mix(v0, v1, v2, v3, v4, v5, v6, v7, x, y, z),
355 mix(v8, v9, v10, v11, v12, v13, v14, v15, x, y, z),
356 w);
357}
358
359BLI_INLINE float fade(float t)
360{
361 return t * t * t * (t * (t * 6.0 - 15.0) + 10.0);
362}
363
364BLI_INLINE float negate_if(float value, uint32_t condition)
365{
366 return (condition != 0u) ? -value : value;
367}
368
370{
371 uint32_t h = hash & 15u;
372 float g = 1u + (h & 7u);
373 return negate_if(g, h & 8u) * x;
374}
375
376BLI_INLINE float noise_grad(uint32_t hash, float x, float y)
377{
378 uint32_t h = hash & 7u;
379 float u = h < 4u ? x : y;
380 float v = 2.0 * (h < 4u ? y : x);
381 return negate_if(u, h & 1u) + negate_if(v, h & 2u);
382}
383
384BLI_INLINE float noise_grad(uint32_t hash, float x, float y, float z)
385{
386 uint32_t h = hash & 15u;
387 float u = h < 8u ? x : y;
388 float vt = ELEM(h, 12u, 14u) ? x : z;
389 float v = h < 4u ? y : vt;
390 return negate_if(u, h & 1u) + negate_if(v, h & 2u);
391}
392
393BLI_INLINE float noise_grad(uint32_t hash, float x, float y, float z, float w)
394{
395 uint32_t h = hash & 31u;
396 float u = h < 24u ? x : y;
397 float v = h < 16u ? y : z;
398 float s = h < 8u ? z : w;
399 return negate_if(u, h & 1u) + negate_if(v, h & 2u) + negate_if(s, h & 4u);
400}
401
402BLI_INLINE float floor_fraction(float x, int &i)
403{
404 float x_floor = math::floor(x);
405 i = int(x_floor);
406 return x - x_floor;
407}
408
409BLI_INLINE float perlin_noise(float position)
410{
411 int X;
412
413 float fx = floor_fraction(position, X);
414
415 float u = fade(fx);
416
417 float r = mix(noise_grad(hash(X), fx), noise_grad(hash(X + 1), fx - 1.0), u);
418
419 return r;
420}
421
423{
424 int X, Y;
425
426 float fx = floor_fraction(position.x, X);
427 float fy = floor_fraction(position.y, Y);
428
429 float u = fade(fx);
430 float v = fade(fy);
431
432 float r = mix(noise_grad(hash(X, Y), fx, fy),
433 noise_grad(hash(X + 1, Y), fx - 1.0, fy),
434 noise_grad(hash(X, Y + 1), fx, fy - 1.0),
435 noise_grad(hash(X + 1, Y + 1), fx - 1.0, fy - 1.0),
436 u,
437 v);
438
439 return r;
440}
441
443{
444 int X, Y, Z;
445
446 float fx = floor_fraction(position.x, X);
447 float fy = floor_fraction(position.y, Y);
448 float fz = floor_fraction(position.z, Z);
449
450 float u = fade(fx);
451 float v = fade(fy);
452 float w = fade(fz);
453
454 float r = mix(noise_grad(hash(X, Y, Z), fx, fy, fz),
455 noise_grad(hash(X + 1, Y, Z), fx - 1, fy, fz),
456 noise_grad(hash(X, Y + 1, Z), fx, fy - 1, fz),
457 noise_grad(hash(X + 1, Y + 1, Z), fx - 1, fy - 1, fz),
458 noise_grad(hash(X, Y, Z + 1), fx, fy, fz - 1),
459 noise_grad(hash(X + 1, Y, Z + 1), fx - 1, fy, fz - 1),
460 noise_grad(hash(X, Y + 1, Z + 1), fx, fy - 1, fz - 1),
461 noise_grad(hash(X + 1, Y + 1, Z + 1), fx - 1, fy - 1, fz - 1),
462 u,
463 v,
464 w);
465
466 return r;
467}
468
470{
471 int X, Y, Z, W;
472
473 float fx = floor_fraction(position.x, X);
474 float fy = floor_fraction(position.y, Y);
475 float fz = floor_fraction(position.z, Z);
476 float fw = floor_fraction(position.w, W);
477
478 float u = fade(fx);
479 float v = fade(fy);
480 float t = fade(fz);
481 float s = fade(fw);
482
483 float r = mix(
484 noise_grad(hash(X, Y, Z, W), fx, fy, fz, fw),
485 noise_grad(hash(X + 1, Y, Z, W), fx - 1.0, fy, fz, fw),
486 noise_grad(hash(X, Y + 1, Z, W), fx, fy - 1.0, fz, fw),
487 noise_grad(hash(X + 1, Y + 1, Z, W), fx - 1.0, fy - 1.0, fz, fw),
488 noise_grad(hash(X, Y, Z + 1, W), fx, fy, fz - 1.0, fw),
489 noise_grad(hash(X + 1, Y, Z + 1, W), fx - 1.0, fy, fz - 1.0, fw),
490 noise_grad(hash(X, Y + 1, Z + 1, W), fx, fy - 1.0, fz - 1.0, fw),
491 noise_grad(hash(X + 1, Y + 1, Z + 1, W), fx - 1.0, fy - 1.0, fz - 1.0, fw),
492 noise_grad(hash(X, Y, Z, W + 1), fx, fy, fz, fw - 1.0),
493 noise_grad(hash(X + 1, Y, Z, W + 1), fx - 1.0, fy, fz, fw - 1.0),
494 noise_grad(hash(X, Y + 1, Z, W + 1), fx, fy - 1.0, fz, fw - 1.0),
495 noise_grad(hash(X + 1, Y + 1, Z, W + 1), fx - 1.0, fy - 1.0, fz, fw - 1.0),
496 noise_grad(hash(X, Y, Z + 1, W + 1), fx, fy, fz - 1.0, fw - 1.0),
497 noise_grad(hash(X + 1, Y, Z + 1, W + 1), fx - 1.0, fy, fz - 1.0, fw - 1.0),
498 noise_grad(hash(X, Y + 1, Z + 1, W + 1), fx, fy - 1.0, fz - 1.0, fw - 1.0),
499 noise_grad(hash(X + 1, Y + 1, Z + 1, W + 1), fx - 1.0, fy - 1.0, fz - 1.0, fw - 1.0),
500 u,
501 v,
502 t,
503 s);
504
505 return r;
506}
507
508/* Signed versions of perlin noise in the range [-1, 1]. The scale values were computed
509 * experimentally by the OSL developers to remap the noise output to the correct range. */
510
511float perlin_signed(float position)
512{
513 float precision_correction = 0.5f * float(math::abs(position) >= 1000000.0f);
514 /* Repeat Perlin noise texture every 100000.0 on each axis to prevent floating point
515 * representation issues. */
516 position = math::mod(position, 100000.0f) + precision_correction;
517
518 return perlin_noise(position) * 0.2500f;
519}
520
521float perlin_signed(float2 position)
522{
523 float2 precision_correction = 0.5f * float2(float(math::abs(position.x) >= 1000000.0f),
524 float(math::abs(position.y) >= 1000000.0f));
525 /* Repeat Perlin noise texture every 100000.0f on each axis to prevent floating point
526 * representation issues. This causes discontinuities every 100000.0f, however at such scales
527 * this usually shouldn't be noticeable. */
528 position = math::mod(position, 100000.0f) + precision_correction;
529
530 return perlin_noise(position) * 0.6616f;
531}
532
533float perlin_signed(float3 position)
534{
535 float3 precision_correction = 0.5f * float3(float(math::abs(position.x) >= 1000000.0f),
536 float(math::abs(position.y) >= 1000000.0f),
537 float(math::abs(position.z) >= 1000000.0f));
538 /* Repeat Perlin noise texture every 100000.0f on each axis to prevent floating point
539 * representation issues. This causes discontinuities every 100000.0f, however at such scales
540 * this usually shouldn't be noticeable. */
541 position = math::mod(position, 100000.0f) + precision_correction;
542
543 return perlin_noise(position) * 0.9820f;
544}
545
546float perlin_signed(float4 position)
547{
548 float4 precision_correction = 0.5f * float4(float(math::abs(position.x) >= 1000000.0f),
549 float(math::abs(position.y) >= 1000000.0f),
550 float(math::abs(position.z) >= 1000000.0f),
551 float(math::abs(position.w) >= 1000000.0f));
552 /* Repeat Perlin noise texture every 100000.0f on each axis to prevent floating point
553 * representation issues. This causes discontinuities every 100000.0f, however at such scales
554 * this usually shouldn't be noticeable. */
555 position = math::mod(position, 100000.0f) + precision_correction;
556
557 return perlin_noise(position) * 0.8344f;
558}
559
560/* Positive versions of perlin noise in the range [0, 1]. */
561
562float perlin(float position)
563{
564 return perlin_signed(position) / 2.0f + 0.5f;
565}
566
567float perlin(float2 position)
568{
569 return perlin_signed(position) / 2.0f + 0.5f;
570}
571
572float perlin(float3 position)
573{
574 return perlin_signed(position) / 2.0f + 0.5f;
575}
576
577float perlin(float4 position)
578{
579 return perlin_signed(position) / 2.0f + 0.5f;
580}
581
582/* Fractal perlin noise. */
583
584/* fBM = Fractal Brownian Motion */
585template<typename T>
587 T p, const float detail, const float roughness, const float lacunarity, const bool normalize)
588{
589 float fscale = 1.0f;
590 float amp = 1.0f;
591 float maxamp = 0.0f;
592 float sum = 0.0f;
593
594 for (int i = 0; i <= int(detail); i++) {
595 float t = perlin_signed(fscale * p);
596 sum += t * amp;
597 maxamp += amp;
598 amp *= roughness;
599 fscale *= lacunarity;
600 }
601 float rmd = detail - std::floor(detail);
602 if (rmd != 0.0f) {
603 float t = perlin_signed(fscale * p);
604 float sum2 = sum + t * amp;
605 return normalize ? mix(0.5f * sum / maxamp + 0.5f, 0.5f * sum2 / (maxamp + amp) + 0.5f, rmd) :
606 mix(sum, sum2, rmd);
607 }
608 return normalize ? 0.5f * sum / maxamp + 0.5f : sum;
609}
610
611/* Explicit instantiation for Wave Texture. */
612template float perlin_fbm<float3>(float3 p,
613 const float detail,
614 const float roughness,
615 const float lacunarity,
616 const bool normalize);
617
618template<typename T>
619float perlin_multi_fractal(T p, const float detail, const float roughness, const float lacunarity)
620{
621 float value = 1.0f;
622 float pwr = 1.0f;
623
624 for (int i = 0; i <= int(detail); i++) {
625 value *= (pwr * perlin_signed(p) + 1.0f);
626 pwr *= roughness;
627 p *= lacunarity;
628 }
629
630 const float rmd = detail - floorf(detail);
631 if (rmd != 0.0f) {
632 value *= (rmd * pwr * perlin_signed(p) + 1.0f); /* correct? */
633 }
634
635 return value;
636}
637
638template<typename T>
640 T p, const float detail, const float roughness, const float lacunarity, const float offset)
641{
642 float pwr = roughness;
643
644 /* First unscaled octave of function; later octaves are scaled. */
645 float value = offset + perlin_signed(p);
646 p *= lacunarity;
647
648 for (int i = 1; i <= int(detail); i++) {
649 float increment = (perlin_signed(p) + offset) * pwr * value;
650 value += increment;
651 pwr *= roughness;
652 p *= lacunarity;
653 }
654
655 const float rmd = detail - floorf(detail);
656 if (rmd != 0.0f) {
657 float increment = (perlin_signed(p) + offset) * pwr * value;
658 value += rmd * increment;
659 }
660
661 return value;
662}
663
664template<typename T>
666 const float detail,
667 const float roughness,
668 const float lacunarity,
669 const float offset,
670 const float gain)
671{
672 float pwr = 1.0f;
673 float value = 0.0f;
674 float weight = 1.0f;
675
676 for (int i = 0; (weight > 0.001f) && (i <= int(detail)); i++) {
677 if (weight > 1.0f) {
678 weight = 1.0f;
679 }
680
681 float signal = (perlin_signed(p) + offset) * pwr;
682 pwr *= roughness;
683 value += weight * signal;
684 weight *= gain * signal;
685 p *= lacunarity;
686 }
687
688 const float rmd = detail - floorf(detail);
689 if ((rmd != 0.0f) && (weight > 0.001f)) {
690 if (weight > 1.0f) {
691 weight = 1.0f;
692 }
693 float signal = (perlin_signed(p) + offset) * pwr;
694 value += rmd * weight * signal;
695 }
696
697 return value;
698}
699
700template<typename T>
702 const float detail,
703 const float roughness,
704 const float lacunarity,
705 const float offset,
706 const float gain)
707{
708 float pwr = roughness;
709
710 float signal = offset - std::abs(perlin_signed(p));
711 signal *= signal;
712 float value = signal;
713 float weight = 1.0f;
714
715 for (int i = 1; i <= int(detail); i++) {
716 p *= lacunarity;
717 weight = std::clamp(signal * gain, 0.0f, 1.0f);
718 signal = offset - std::abs(perlin_signed(p));
719 signal *= signal;
720 signal *= weight;
721 value += signal * pwr;
722 pwr *= roughness;
723 }
724
725 return value;
726}
727
728enum {
734};
735
736template<typename T>
738 float detail,
739 float roughness,
740 float lacunarity,
741 float offset,
742 float gain,
743 int type,
744 bool normalize)
745{
746 switch (type) {
748 return perlin_multi_fractal<T>(p, detail, roughness, lacunarity);
749 }
751 return perlin_fbm<T>(p, detail, roughness, lacunarity, normalize);
752 }
754 return perlin_hybrid_multi_fractal<T>(p, detail, roughness, lacunarity, offset, gain);
755 }
757 return perlin_ridged_multi_fractal<T>(p, detail, roughness, lacunarity, offset, gain);
758 }
760 return perlin_hetero_terrain<T>(p, detail, roughness, lacunarity, offset);
761 }
762 default: {
763 return 0.0;
764 }
765 }
766}
767
768/* The following offset functions generate random offsets to be added to
769 * positions to act as a seed since the noise functions don't have seed values.
770 * The offset's components are in the range [100, 200], not too high to cause
771 * bad precision and not too small to be noticeable. We use float seed because
772 * OSL only supports float hashes and we need to maintain compatibility with it.
773 */
774
776{
777 return 100.0f + hash_float_to_float(seed) * 100.0f;
778}
779
781{
782 return float2(100.0f + hash_float_to_float(float2(seed, 0.0f)) * 100.0f,
783 100.0f + hash_float_to_float(float2(seed, 1.0f)) * 100.0f);
784}
785
787{
788 return float3(100.0f + hash_float_to_float(float2(seed, 0.0f)) * 100.0f,
789 100.0f + hash_float_to_float(float2(seed, 1.0f)) * 100.0f,
790 100.0f + hash_float_to_float(float2(seed, 2.0f)) * 100.0f);
791}
792
794{
795 return float4(100.0f + hash_float_to_float(float2(seed, 0.0f)) * 100.0f,
796 100.0f + hash_float_to_float(float2(seed, 1.0f)) * 100.0f,
797 100.0f + hash_float_to_float(float2(seed, 2.0f)) * 100.0f,
798 100.0f + hash_float_to_float(float2(seed, 3.0f)) * 100.0f);
799}
800
801/* Perlin noises to be added to the position to distort other noises. */
802
803BLI_INLINE float perlin_distortion(float position, float strength)
804{
805 return perlin_signed(position + random_float_offset(0.0)) * strength;
806}
807
808BLI_INLINE float2 perlin_distortion(float2 position, float strength)
809{
810 return float2(perlin_signed(position + random_float2_offset(0.0f)) * strength,
811 perlin_signed(position + random_float2_offset(1.0f)) * strength);
812}
813
814BLI_INLINE float3 perlin_distortion(float3 position, float strength)
815{
816 return float3(perlin_signed(position + random_float3_offset(0.0f)) * strength,
817 perlin_signed(position + random_float3_offset(1.0f)) * strength,
818 perlin_signed(position + random_float3_offset(2.0f)) * strength);
819}
820
821BLI_INLINE float4 perlin_distortion(float4 position, float strength)
822{
823 return float4(perlin_signed(position + random_float4_offset(0.0f)) * strength,
824 perlin_signed(position + random_float4_offset(1.0f)) * strength,
825 perlin_signed(position + random_float4_offset(2.0f)) * strength,
826 perlin_signed(position + random_float4_offset(3.0f)) * strength);
827}
828
829/* Distorted fractal perlin noise. */
830
831template<typename T>
833 float detail,
834 float roughness,
835 float lacunarity,
836 float offset,
837 float gain,
838 float distortion,
839 int type,
840 bool normalize)
841{
842 position += perlin_distortion(position, distortion);
843 return perlin_select<T>(position, detail, roughness, lacunarity, offset, gain, type, normalize);
844}
845
846template float perlin_fractal_distorted<float>(float position,
847 float detail,
848 float roughness,
849 float lacunarity,
850 float offset,
851 float gain,
852 float distortion,
853 int type,
854 bool normalize);
856 float detail,
857 float roughness,
858 float lacunarity,
859 float offset,
860 float gain,
861 float distortion,
862 int type,
863 bool normalize);
865 float detail,
866 float roughness,
867 float lacunarity,
868 float offset,
869 float gain,
870 float distortion,
871 int type,
872 bool normalize);
874 float detail,
875 float roughness,
876 float lacunarity,
877 float offset,
878 float gain,
879 float distortion,
880 int type,
881 bool normalize);
882
883/* Distorted fractal perlin noise that outputs a float3. The arbitrary seeds are for
884 * compatibility with shading functions. */
885
887 float detail,
888 float roughness,
889 float lacunarity,
890 float offset,
891 float gain,
892 float distortion,
893 int type,
894 bool normalize)
895{
896 position += perlin_distortion(position, distortion);
897 return float3(
898 perlin_select<float>(position, detail, roughness, lacunarity, offset, gain, type, normalize),
900 detail,
901 roughness,
902 lacunarity,
903 offset,
904 gain,
905 type,
906 normalize),
908 detail,
909 roughness,
910 lacunarity,
911 offset,
912 gain,
913 type,
914 normalize));
915}
916
918 float detail,
919 float roughness,
920 float lacunarity,
921 float offset,
922 float gain,
923 float distortion,
924 int type,
925 bool normalize)
926{
927 position += perlin_distortion(position, distortion);
929 position, detail, roughness, lacunarity, offset, gain, type, normalize),
931 detail,
932 roughness,
933 lacunarity,
934 offset,
935 gain,
936 type,
937 normalize),
939 detail,
940 roughness,
941 lacunarity,
942 offset,
943 gain,
944 type,
945 normalize));
946}
947
949 float detail,
950 float roughness,
951 float lacunarity,
952 float offset,
953 float gain,
954 float distortion,
955 int type,
956 bool normalize)
957{
958 position += perlin_distortion(position, distortion);
960 position, detail, roughness, lacunarity, offset, gain, type, normalize),
962 detail,
963 roughness,
964 lacunarity,
965 offset,
966 gain,
967 type,
968 normalize),
970 detail,
971 roughness,
972 lacunarity,
973 offset,
974 gain,
975 type,
976 normalize));
977}
978
980 float detail,
981 float roughness,
982 float lacunarity,
983 float offset,
984 float gain,
985 float distortion,
986 int type,
987 bool normalize)
988{
989 position += perlin_distortion(position, distortion);
991 position, detail, roughness, lacunarity, offset, gain, type, normalize),
993 detail,
994 roughness,
995 lacunarity,
996 offset,
997 gain,
998 type,
999 normalize),
1001 detail,
1002 roughness,
1003 lacunarity,
1004 offset,
1005 gain,
1006 type,
1007 normalize));
1008}
1009
1012/* -------------------------------------------------------------------- */
1032/* Ensure to align with DNA. */
1033
1034enum {
1039};
1040
1041enum {
1047};
1048
1049/* ***** Distances ***** */
1050
1051float voronoi_distance(const float a, const float b)
1052{
1053 return std::abs(b - a);
1054}
1055
1057{
1058 switch (params.metric) {
1060 return math::distance(a, b);
1062 return std::abs(a.x - b.x) + std::abs(a.y - b.y);
1064 return std::max(std::abs(a.x - b.x), std::abs(a.y - b.y));
1066 return std::pow(std::pow(std::abs(a.x - b.x), params.exponent) +
1067 std::pow(std::abs(a.y - b.y), params.exponent),
1068 1.0f / params.exponent);
1069 default:
1071 break;
1072 }
1073 return 0.0f;
1074}
1075
1077{
1078 switch (params.metric) {
1080 return math::distance(a, b);
1082 return std::abs(a.x - b.x) + std::abs(a.y - b.y) + std::abs(a.z - b.z);
1084 return std::max(std::abs(a.x - b.x), std::max(std::abs(a.y - b.y), std::abs(a.z - b.z)));
1086 return std::pow(std::pow(std::abs(a.x - b.x), params.exponent) +
1087 std::pow(std::abs(a.y - b.y), params.exponent) +
1088 std::pow(std::abs(a.z - b.z), params.exponent),
1089 1.0f / params.exponent);
1090 default:
1092 break;
1093 }
1094 return 0.0f;
1095}
1096
1098{
1099 switch (params.metric) {
1101 return math::distance(a, b);
1103 return std::abs(a.x - b.x) + std::abs(a.y - b.y) + std::abs(a.z - b.z) + std::abs(a.w - b.w);
1105 return std::max(
1106 std::abs(a.x - b.x),
1107 std::max(std::abs(a.y - b.y), std::max(std::abs(a.z - b.z), std::abs(a.w - b.w))));
1109 return std::pow(std::pow(std::abs(a.x - b.x), params.exponent) +
1110 std::pow(std::abs(a.y - b.y), params.exponent) +
1111 std::pow(std::abs(a.z - b.z), params.exponent) +
1112 std::pow(std::abs(a.w - b.w), params.exponent),
1113 1.0f / params.exponent);
1114 default:
1116 break;
1117 }
1118 return 0.0f;
1119}
1120
1121/* **** 1D Voronoi **** */
1122
1123float4 voronoi_position(const float coord)
1124{
1125 return {0.0f, 0.0f, 0.0f, coord};
1126}
1127
1129{
1130 float cellPosition = floorf(coord);
1131 float localPosition = coord - cellPosition;
1132
1133 float minDistance = FLT_MAX;
1134 float targetOffset = 0.0f;
1135 float targetPosition = 0.0f;
1136 for (int i = -1; i <= 1; i++) {
1137 float cellOffset = i;
1138 float pointPosition = cellOffset +
1139 hash_float_to_float(cellPosition + cellOffset) * params.randomness;
1140 float distanceToPoint = voronoi_distance(pointPosition, localPosition);
1141 if (distanceToPoint < minDistance) {
1142 targetOffset = cellOffset;
1143 minDistance = distanceToPoint;
1144 targetPosition = pointPosition;
1145 }
1146 }
1147
1148 VoronoiOutput octave;
1149 octave.distance = minDistance;
1150 octave.color = hash_float_to_float3(cellPosition + targetOffset);
1151 octave.position = voronoi_position(targetPosition + cellPosition);
1152 return octave;
1153}
1154
1156 const float coord,
1157 const bool calc_color)
1158{
1159 float cellPosition = floorf(coord);
1160 float localPosition = coord - cellPosition;
1161
1162 float smoothDistance = 0.0f;
1163 float smoothPosition = 0.0f;
1164 float3 smoothColor = {0.0f, 0.0f, 0.0f};
1165 float h = -1.0f;
1166 for (int i = -2; i <= 2; i++) {
1167 float cellOffset = i;
1168 float pointPosition = cellOffset +
1169 hash_float_to_float(cellPosition + cellOffset) * params.randomness;
1170 float distanceToPoint = voronoi_distance(pointPosition, localPosition);
1171 h = h == -1.0f ?
1172 1.0f :
1173 smoothstep(
1174 0.0f, 1.0f, 0.5f + 0.5f * (smoothDistance - distanceToPoint) / params.smoothness);
1175 float correctionFactor = params.smoothness * h * (1.0f - h);
1176 smoothDistance = mix(smoothDistance, distanceToPoint, h) - correctionFactor;
1177 correctionFactor /= 1.0f + 3.0f * params.smoothness;
1178 float3 cellColor = hash_float_to_float3(cellPosition + cellOffset);
1179 if (calc_color) {
1180 /* Only compute Color output if necessary, as it is very expensive. */
1181 smoothColor = mix(smoothColor, cellColor, h) - correctionFactor;
1182 }
1183 smoothPosition = mix(smoothPosition, pointPosition, h) - correctionFactor;
1184 }
1185
1186 VoronoiOutput octave;
1187 octave.distance = smoothDistance;
1188 octave.color = smoothColor;
1189 octave.position = voronoi_position(cellPosition + smoothPosition);
1190 return octave;
1191}
1192
1194{
1195 float cellPosition = floorf(coord);
1196 float localPosition = coord - cellPosition;
1197
1198 float distanceF1 = FLT_MAX;
1199 float distanceF2 = FLT_MAX;
1200 float offsetF1 = 0.0f;
1201 float positionF1 = 0.0f;
1202 float offsetF2 = 0.0f;
1203 float positionF2 = 0.0f;
1204 for (int i = -1; i <= 1; i++) {
1205 float cellOffset = i;
1206 float pointPosition = cellOffset +
1207 hash_float_to_float(cellPosition + cellOffset) * params.randomness;
1208 float distanceToPoint = voronoi_distance(pointPosition, localPosition);
1209 if (distanceToPoint < distanceF1) {
1210 distanceF2 = distanceF1;
1211 distanceF1 = distanceToPoint;
1212 offsetF2 = offsetF1;
1213 offsetF1 = cellOffset;
1214 positionF2 = positionF1;
1215 positionF1 = pointPosition;
1216 }
1217 else if (distanceToPoint < distanceF2) {
1218 distanceF2 = distanceToPoint;
1219 offsetF2 = cellOffset;
1220 positionF2 = pointPosition;
1221 }
1222 }
1223
1224 VoronoiOutput octave;
1225 octave.distance = distanceF2;
1226 octave.color = hash_float_to_float3(cellPosition + offsetF2);
1227 octave.position = voronoi_position(positionF2 + cellPosition);
1228 return octave;
1229}
1230
1231float voronoi_distance_to_edge(const VoronoiParams &params, const float coord)
1232{
1233 float cellPosition = floorf(coord);
1234 float localPosition = coord - cellPosition;
1235
1236 float midPointPosition = hash_float_to_float(cellPosition) * params.randomness;
1237 float leftPointPosition = -1.0f + hash_float_to_float(cellPosition - 1.0f) * params.randomness;
1238 float rightPointPosition = 1.0f + hash_float_to_float(cellPosition + 1.0f) * params.randomness;
1239 float distanceToMidLeft = fabsf((midPointPosition + leftPointPosition) / 2.0f - localPosition);
1240 float distanceToMidRight = fabsf((midPointPosition + rightPointPosition) / 2.0f - localPosition);
1241
1242 return math::min(distanceToMidLeft, distanceToMidRight);
1243}
1244
1245float voronoi_n_sphere_radius(const VoronoiParams &params, const float coord)
1246{
1247 float cellPosition = floorf(coord);
1248 float localPosition = coord - cellPosition;
1249
1250 float closestPoint = 0.0f;
1251 float closestPointOffset = 0.0f;
1252 float minDistance = FLT_MAX;
1253 for (int i = -1; i <= 1; i++) {
1254 float cellOffset = i;
1255 float pointPosition = cellOffset +
1256 hash_float_to_float(cellPosition + cellOffset) * params.randomness;
1257 float distanceToPoint = fabsf(pointPosition - localPosition);
1258 if (distanceToPoint < minDistance) {
1259 minDistance = distanceToPoint;
1260 closestPoint = pointPosition;
1261 closestPointOffset = cellOffset;
1262 }
1263 }
1264
1265 minDistance = FLT_MAX;
1266 float closestPointToClosestPoint = 0.0f;
1267 for (int i = -1; i <= 1; i++) {
1268 if (i == 0) {
1269 continue;
1270 }
1271 float cellOffset = i + closestPointOffset;
1272 float pointPosition = cellOffset +
1273 hash_float_to_float(cellPosition + cellOffset) * params.randomness;
1274 float distanceToPoint = fabsf(closestPoint - pointPosition);
1275 if (distanceToPoint < minDistance) {
1276 minDistance = distanceToPoint;
1277 closestPointToClosestPoint = pointPosition;
1278 }
1279 }
1280
1281 return fabsf(closestPointToClosestPoint - closestPoint) / 2.0f;
1282}
1283
1284/* **** 2D Voronoi **** */
1285
1287{
1288 return {coord.x, coord.y, 0.0f, 0.0f};
1289}
1290
1292{
1293 float2 cellPosition = math::floor(coord);
1294 float2 localPosition = coord - cellPosition;
1295
1296 float minDistance = FLT_MAX;
1297 float2 targetOffset = {0.0f, 0.0f};
1298 float2 targetPosition = {0.0f, 0.0f};
1299 for (int j = -1; j <= 1; j++) {
1300 for (int i = -1; i <= 1; i++) {
1301 float2 cellOffset(i, j);
1302 float2 pointPosition = cellOffset +
1303 hash_float_to_float2(cellPosition + cellOffset) * params.randomness;
1304 float distanceToPoint = voronoi_distance(pointPosition, localPosition, params);
1305 if (distanceToPoint < minDistance) {
1306 targetOffset = cellOffset;
1307 minDistance = distanceToPoint;
1308 targetPosition = pointPosition;
1309 }
1310 }
1311 }
1312
1313 VoronoiOutput octave;
1314 octave.distance = minDistance;
1315 octave.color = hash_float_to_float3(cellPosition + targetOffset);
1316 octave.position = voronoi_position(targetPosition + cellPosition);
1317 return octave;
1318}
1319
1321 const float2 coord,
1322 const bool calc_color)
1323{
1324 float2 cellPosition = math::floor(coord);
1325 float2 localPosition = coord - cellPosition;
1326
1327 float smoothDistance = 0.0f;
1328 float3 smoothColor = {0.0f, 0.0f, 0.0f};
1329 float2 smoothPosition = {0.0f, 0.0f};
1330 float h = -1.0f;
1331 for (int j = -2; j <= 2; j++) {
1332 for (int i = -2; i <= 2; i++) {
1333 float2 cellOffset(i, j);
1334 float2 pointPosition = cellOffset +
1335 hash_float_to_float2(cellPosition + cellOffset) * params.randomness;
1336 float distanceToPoint = voronoi_distance(pointPosition, localPosition, params);
1337 h = h == -1.0f ?
1338 1.0f :
1339 smoothstep(0.0f,
1340 1.0f,
1341 0.5f + 0.5f * (smoothDistance - distanceToPoint) / params.smoothness);
1342 float correctionFactor = params.smoothness * h * (1.0f - h);
1343 smoothDistance = mix(smoothDistance, distanceToPoint, h) - correctionFactor;
1344 correctionFactor /= 1.0f + 3.0f * params.smoothness;
1345 if (calc_color) {
1346 /* Only compute Color output if necessary, as it is very expensive. */
1347 float3 cellColor = hash_float_to_float3(cellPosition + cellOffset);
1348 smoothColor = mix(smoothColor, cellColor, h) - correctionFactor;
1349 }
1350 smoothPosition = mix(smoothPosition, pointPosition, h) - correctionFactor;
1351 }
1352 }
1353
1354 VoronoiOutput octave;
1355 octave.distance = smoothDistance;
1356 octave.color = smoothColor;
1357 octave.position = voronoi_position(cellPosition + smoothPosition);
1358 return octave;
1359}
1360
1362{
1363 float2 cellPosition = math::floor(coord);
1364 float2 localPosition = coord - cellPosition;
1365
1366 float distanceF1 = FLT_MAX;
1367 float distanceF2 = FLT_MAX;
1368 float2 offsetF1 = {0.0f, 0.0f};
1369 float2 positionF1 = {0.0f, 0.0f};
1370 float2 offsetF2 = {0.0f, 0.0f};
1371 float2 positionF2 = {0.0f, 0.0f};
1372 for (int j = -1; j <= 1; j++) {
1373 for (int i = -1; i <= 1; i++) {
1374 float2 cellOffset(i, j);
1375 float2 pointPosition = cellOffset +
1376 hash_float_to_float2(cellPosition + cellOffset) * params.randomness;
1377 float distanceToPoint = voronoi_distance(pointPosition, localPosition, params);
1378 if (distanceToPoint < distanceF1) {
1379 distanceF2 = distanceF1;
1380 distanceF1 = distanceToPoint;
1381 offsetF2 = offsetF1;
1382 offsetF1 = cellOffset;
1383 positionF2 = positionF1;
1384 positionF1 = pointPosition;
1385 }
1386 else if (distanceToPoint < distanceF2) {
1387 distanceF2 = distanceToPoint;
1388 offsetF2 = cellOffset;
1389 positionF2 = pointPosition;
1390 }
1391 }
1392 }
1393
1394 VoronoiOutput octave;
1395 octave.distance = distanceF2;
1396 octave.color = hash_float_to_float3(cellPosition + offsetF2);
1397 octave.position = voronoi_position(positionF2 + cellPosition);
1398 return octave;
1399}
1400
1402{
1403 float2 cellPosition = math::floor(coord);
1404 float2 localPosition = coord - cellPosition;
1405
1406 float2 vectorToClosest = {0.0f, 0.0f};
1407 float minDistance = FLT_MAX;
1408 for (int j = -1; j <= 1; j++) {
1409 for (int i = -1; i <= 1; i++) {
1410 float2 cellOffset(i, j);
1411 float2 vectorToPoint = cellOffset +
1412 hash_float_to_float2(cellPosition + cellOffset) * params.randomness -
1413 localPosition;
1414 float distanceToPoint = math::dot(vectorToPoint, vectorToPoint);
1415 if (distanceToPoint < minDistance) {
1416 minDistance = distanceToPoint;
1417 vectorToClosest = vectorToPoint;
1418 }
1419 }
1420 }
1421
1422 minDistance = FLT_MAX;
1423 for (int j = -1; j <= 1; j++) {
1424 for (int i = -1; i <= 1; i++) {
1425 float2 cellOffset(i, j);
1426 float2 vectorToPoint = cellOffset +
1427 hash_float_to_float2(cellPosition + cellOffset) * params.randomness -
1428 localPosition;
1429 float2 perpendicularToEdge = vectorToPoint - vectorToClosest;
1430 if (math::dot(perpendicularToEdge, perpendicularToEdge) > 0.0001f) {
1431 float distanceToEdge = math::dot((vectorToClosest + vectorToPoint) / 2.0f,
1432 math::normalize(perpendicularToEdge));
1433 minDistance = math::min(minDistance, distanceToEdge);
1434 }
1435 }
1436 }
1437
1438 return minDistance;
1439}
1440
1442{
1443 float2 cellPosition = math::floor(coord);
1444 float2 localPosition = coord - cellPosition;
1445
1446 float2 closestPoint = {0.0f, 0.0f};
1447 float2 closestPointOffset = {0.0f, 0.0f};
1448 float minDistance = FLT_MAX;
1449 for (int j = -1; j <= 1; j++) {
1450 for (int i = -1; i <= 1; i++) {
1451 float2 cellOffset(i, j);
1452 float2 pointPosition = cellOffset +
1453 hash_float_to_float2(cellPosition + cellOffset) * params.randomness;
1454 float distanceToPoint = math::distance(pointPosition, localPosition);
1455 if (distanceToPoint < minDistance) {
1456 minDistance = distanceToPoint;
1457 closestPoint = pointPosition;
1458 closestPointOffset = cellOffset;
1459 }
1460 }
1461 }
1462
1463 minDistance = FLT_MAX;
1464 float2 closestPointToClosestPoint = {0.0f, 0.0f};
1465 for (int j = -1; j <= 1; j++) {
1466 for (int i = -1; i <= 1; i++) {
1467 if (i == 0 && j == 0) {
1468 continue;
1469 }
1470 float2 cellOffset = float2(i, j) + closestPointOffset;
1471 float2 pointPosition = cellOffset +
1472 hash_float_to_float2(cellPosition + cellOffset) * params.randomness;
1473 float distanceToPoint = math::distance(closestPoint, pointPosition);
1474 if (distanceToPoint < minDistance) {
1475 minDistance = distanceToPoint;
1476 closestPointToClosestPoint = pointPosition;
1477 }
1478 }
1479 }
1480
1481 return math::distance(closestPointToClosestPoint, closestPoint) / 2.0f;
1482}
1483
1484/* **** 3D Voronoi **** */
1485
1487{
1488 return {coord.x, coord.y, coord.z, 0.0f};
1489}
1490
1492{
1493 float3 cellPosition = math::floor(coord);
1494 float3 localPosition = coord - cellPosition;
1495
1496 float minDistance = FLT_MAX;
1497 float3 targetOffset = {0.0f, 0.0f, 0.0f};
1498 float3 targetPosition = {0.0f, 0.0f, 0.0f};
1499 for (int k = -1; k <= 1; k++) {
1500 for (int j = -1; j <= 1; j++) {
1501 for (int i = -1; i <= 1; i++) {
1502 float3 cellOffset(i, j, k);
1503 float3 pointPosition = cellOffset +
1504 hash_float_to_float3(cellPosition + cellOffset) * params.randomness;
1505 float distanceToPoint = voronoi_distance(pointPosition, localPosition, params);
1506 if (distanceToPoint < minDistance) {
1507 targetOffset = cellOffset;
1508 minDistance = distanceToPoint;
1509 targetPosition = pointPosition;
1510 }
1511 }
1512 }
1513 }
1514
1515 VoronoiOutput octave;
1516 octave.distance = minDistance;
1517 octave.color = hash_float_to_float3(cellPosition + targetOffset);
1518 octave.position = voronoi_position(targetPosition + cellPosition);
1519 return octave;
1520}
1521
1523 const float3 coord,
1524 const bool calc_color)
1525{
1526 float3 cellPosition = math::floor(coord);
1527 float3 localPosition = coord - cellPosition;
1528
1529 float smoothDistance = 0.0f;
1530 float3 smoothColor = {0.0f, 0.0f, 0.0f};
1531 float3 smoothPosition = {0.0f, 0.0f, 0.0f};
1532 float h = -1.0f;
1533 for (int k = -2; k <= 2; k++) {
1534 for (int j = -2; j <= 2; j++) {
1535 for (int i = -2; i <= 2; i++) {
1536 float3 cellOffset(i, j, k);
1537 float3 pointPosition = cellOffset +
1538 hash_float_to_float3(cellPosition + cellOffset) * params.randomness;
1539 float distanceToPoint = voronoi_distance(pointPosition, localPosition, params);
1540 h = h == -1.0f ?
1541 1.0f :
1542 smoothstep(0.0f,
1543 1.0f,
1544 0.5f + 0.5f * (smoothDistance - distanceToPoint) / params.smoothness);
1545 float correctionFactor = params.smoothness * h * (1.0f - h);
1546 smoothDistance = mix(smoothDistance, distanceToPoint, h) - correctionFactor;
1547 correctionFactor /= 1.0f + 3.0f * params.smoothness;
1548 if (calc_color) {
1549 /* Only compute Color output if necessary, as it is very expensive. */
1550 float3 cellColor = hash_float_to_float3(cellPosition + cellOffset);
1551 smoothColor = mix(smoothColor, cellColor, h) - correctionFactor;
1552 }
1553 smoothPosition = mix(smoothPosition, pointPosition, h) - correctionFactor;
1554 }
1555 }
1556 }
1557
1558 VoronoiOutput octave;
1559 octave.distance = smoothDistance;
1560 octave.color = smoothColor;
1561 octave.position = voronoi_position(cellPosition + smoothPosition);
1562 return octave;
1563}
1564
1566{
1567 float3 cellPosition = math::floor(coord);
1568 float3 localPosition = coord - cellPosition;
1569
1570 float distanceF1 = FLT_MAX;
1571 float distanceF2 = FLT_MAX;
1572 float3 offsetF1 = {0.0f, 0.0f, 0.0f};
1573 float3 positionF1 = {0.0f, 0.0f, 0.0f};
1574 float3 offsetF2 = {0.0f, 0.0f, 0.0f};
1575 float3 positionF2 = {0.0f, 0.0f, 0.0f};
1576 for (int k = -1; k <= 1; k++) {
1577 for (int j = -1; j <= 1; j++) {
1578 for (int i = -1; i <= 1; i++) {
1579 float3 cellOffset(i, j, k);
1580 float3 pointPosition = cellOffset +
1581 hash_float_to_float3(cellPosition + cellOffset) * params.randomness;
1582 float distanceToPoint = voronoi_distance(pointPosition, localPosition, params);
1583 if (distanceToPoint < distanceF1) {
1584 distanceF2 = distanceF1;
1585 distanceF1 = distanceToPoint;
1586 offsetF2 = offsetF1;
1587 offsetF1 = cellOffset;
1588 positionF2 = positionF1;
1589 positionF1 = pointPosition;
1590 }
1591 else if (distanceToPoint < distanceF2) {
1592 distanceF2 = distanceToPoint;
1593 offsetF2 = cellOffset;
1594 positionF2 = pointPosition;
1595 }
1596 }
1597 }
1598 }
1599
1600 VoronoiOutput octave;
1601 octave.distance = distanceF2;
1602 octave.color = hash_float_to_float3(cellPosition + offsetF2);
1603 octave.position = voronoi_position(positionF2 + cellPosition);
1604 return octave;
1605}
1606
1608{
1609 float3 cellPosition = math::floor(coord);
1610 float3 localPosition = coord - cellPosition;
1611
1612 float3 vectorToClosest = {0.0f, 0.0f, 0.0f};
1613 float minDistance = FLT_MAX;
1614 for (int k = -1; k <= 1; k++) {
1615 for (int j = -1; j <= 1; j++) {
1616 for (int i = -1; i <= 1; i++) {
1617 float3 cellOffset(i, j, k);
1618 float3 vectorToPoint = cellOffset +
1619 hash_float_to_float3(cellPosition + cellOffset) *
1620 params.randomness -
1621 localPosition;
1622 float distanceToPoint = math::dot(vectorToPoint, vectorToPoint);
1623 if (distanceToPoint < minDistance) {
1624 minDistance = distanceToPoint;
1625 vectorToClosest = vectorToPoint;
1626 }
1627 }
1628 }
1629 }
1630
1631 minDistance = FLT_MAX;
1632 for (int k = -1; k <= 1; k++) {
1633 for (int j = -1; j <= 1; j++) {
1634 for (int i = -1; i <= 1; i++) {
1635 float3 cellOffset(i, j, k);
1636 float3 vectorToPoint = cellOffset +
1637 hash_float_to_float3(cellPosition + cellOffset) *
1638 params.randomness -
1639 localPosition;
1640 float3 perpendicularToEdge = vectorToPoint - vectorToClosest;
1641 if (math::dot(perpendicularToEdge, perpendicularToEdge) > 0.0001f) {
1642 float distanceToEdge = math::dot((vectorToClosest + vectorToPoint) / 2.0f,
1643 math::normalize(perpendicularToEdge));
1644 minDistance = math::min(minDistance, distanceToEdge);
1645 }
1646 }
1647 }
1648 }
1649
1650 return minDistance;
1651}
1652
1654{
1655 float3 cellPosition = math::floor(coord);
1656 float3 localPosition = coord - cellPosition;
1657
1658 float3 closestPoint = {0.0f, 0.0f, 0.0f};
1659 float3 closestPointOffset = {0.0f, 0.0f, 0.0f};
1660 float minDistance = FLT_MAX;
1661 for (int k = -1; k <= 1; k++) {
1662 for (int j = -1; j <= 1; j++) {
1663 for (int i = -1; i <= 1; i++) {
1664 float3 cellOffset(i, j, k);
1665 float3 pointPosition = cellOffset +
1666 hash_float_to_float3(cellPosition + cellOffset) * params.randomness;
1667 float distanceToPoint = math::distance(pointPosition, localPosition);
1668 if (distanceToPoint < minDistance) {
1669 minDistance = distanceToPoint;
1670 closestPoint = pointPosition;
1671 closestPointOffset = cellOffset;
1672 }
1673 }
1674 }
1675 }
1676
1677 minDistance = FLT_MAX;
1678 float3 closestPointToClosestPoint = {0.0f, 0.0f, 0.0f};
1679 for (int k = -1; k <= 1; k++) {
1680 for (int j = -1; j <= 1; j++) {
1681 for (int i = -1; i <= 1; i++) {
1682 if (i == 0 && j == 0 && k == 0) {
1683 continue;
1684 }
1685 float3 cellOffset = float3(i, j, k) + closestPointOffset;
1686 float3 pointPosition = cellOffset +
1687 hash_float_to_float3(cellPosition + cellOffset) * params.randomness;
1688 float distanceToPoint = math::distance(closestPoint, pointPosition);
1689 if (distanceToPoint < minDistance) {
1690 minDistance = distanceToPoint;
1691 closestPointToClosestPoint = pointPosition;
1692 }
1693 }
1694 }
1695 }
1696
1697 return math::distance(closestPointToClosestPoint, closestPoint) / 2.0f;
1698}
1699
1700/* **** 4D Voronoi **** */
1701
1703{
1704 return coord;
1705}
1706
1708{
1709 float4 cellPosition = math::floor(coord);
1710 float4 localPosition = coord - cellPosition;
1711
1712 float minDistance = FLT_MAX;
1713 float4 targetOffset = {0.0f, 0.0f, 0.0f, 0.0f};
1714 float4 targetPosition = {0.0f, 0.0f, 0.0f, 0.0f};
1715 for (int u = -1; u <= 1; u++) {
1716 for (int k = -1; k <= 1; k++) {
1717 for (int j = -1; j <= 1; j++) {
1718 for (int i = -1; i <= 1; i++) {
1719 float4 cellOffset(i, j, k, u);
1720 float4 pointPosition = cellOffset + hash_float_to_float4(cellPosition + cellOffset) *
1721 params.randomness;
1722 float distanceToPoint = voronoi_distance(pointPosition, localPosition, params);
1723 if (distanceToPoint < minDistance) {
1724 targetOffset = cellOffset;
1725 minDistance = distanceToPoint;
1726 targetPosition = pointPosition;
1727 }
1728 }
1729 }
1730 }
1731 }
1732
1733 VoronoiOutput octave;
1734 octave.distance = minDistance;
1735 octave.color = hash_float_to_float3(cellPosition + targetOffset);
1736 octave.position = voronoi_position(targetPosition + cellPosition);
1737 return octave;
1738}
1739
1741 const float4 coord,
1742 const bool calc_color)
1743{
1744 float4 cellPosition = math::floor(coord);
1745 float4 localPosition = coord - cellPosition;
1746
1747 float smoothDistance = 0.0f;
1748 float3 smoothColor = {0.0f, 0.0f, 0.0f};
1749 float4 smoothPosition = {0.0f, 0.0f, 0.0f, 0.0f};
1750 float h = -1.0f;
1751 for (int u = -2; u <= 2; u++) {
1752 for (int k = -2; k <= 2; k++) {
1753 for (int j = -2; j <= 2; j++) {
1754 for (int i = -2; i <= 2; i++) {
1755 float4 cellOffset(i, j, k, u);
1756 float4 pointPosition = cellOffset + hash_float_to_float4(cellPosition + cellOffset) *
1757 params.randomness;
1758 float distanceToPoint = voronoi_distance(pointPosition, localPosition, params);
1759 h = h == -1.0f ?
1760 1.0f :
1761 smoothstep(0.0f,
1762 1.0f,
1763 0.5f + 0.5f * (smoothDistance - distanceToPoint) / params.smoothness);
1764 float correctionFactor = params.smoothness * h * (1.0f - h);
1765 smoothDistance = mix(smoothDistance, distanceToPoint, h) - correctionFactor;
1766 correctionFactor /= 1.0f + 3.0f * params.smoothness;
1767 if (calc_color) {
1768 /* Only compute Color output if necessary, as it is very expensive. */
1769 float3 cellColor = hash_float_to_float3(cellPosition + cellOffset);
1770 smoothColor = mix(smoothColor, cellColor, h) - correctionFactor;
1771 }
1772 smoothPosition = mix(smoothPosition, pointPosition, h) - correctionFactor;
1773 }
1774 }
1775 }
1776 }
1777
1778 VoronoiOutput octave;
1779 octave.distance = smoothDistance;
1780 octave.color = smoothColor;
1781 octave.position = voronoi_position(cellPosition + smoothPosition);
1782 return octave;
1783}
1784
1786{
1787 float4 cellPosition = math::floor(coord);
1788 float4 localPosition = coord - cellPosition;
1789
1790 float distanceF1 = FLT_MAX;
1791 float distanceF2 = FLT_MAX;
1792 float4 offsetF1 = {0.0f, 0.0f, 0.0f, 0.0f};
1793 float4 positionF1 = {0.0f, 0.0f, 0.0f, 0.0f};
1794 float4 offsetF2 = {0.0f, 0.0f, 0.0f, 0.0f};
1795 float4 positionF2 = {0.0f, 0.0f, 0.0f, 0.0f};
1796 for (int u = -1; u <= 1; u++) {
1797 for (int k = -1; k <= 1; k++) {
1798 for (int j = -1; j <= 1; j++) {
1799 for (int i = -1; i <= 1; i++) {
1800 float4 cellOffset(i, j, k, u);
1801 float4 pointPosition = cellOffset + hash_float_to_float4(cellPosition + cellOffset) *
1802 params.randomness;
1803 float distanceToPoint = voronoi_distance(pointPosition, localPosition, params);
1804 if (distanceToPoint < distanceF1) {
1805 distanceF2 = distanceF1;
1806 distanceF1 = distanceToPoint;
1807 offsetF2 = offsetF1;
1808 offsetF1 = cellOffset;
1809 positionF2 = positionF1;
1810 positionF1 = pointPosition;
1811 }
1812 else if (distanceToPoint < distanceF2) {
1813 distanceF2 = distanceToPoint;
1814 offsetF2 = cellOffset;
1815 positionF2 = pointPosition;
1816 }
1817 }
1818 }
1819 }
1820 }
1821
1822 VoronoiOutput octave;
1823 octave.distance = distanceF2;
1824 octave.color = hash_float_to_float3(cellPosition + offsetF2);
1825 octave.position = voronoi_position(positionF2 + cellPosition);
1826 return octave;
1827}
1828
1830{
1831 float4 cellPosition = math::floor(coord);
1832 float4 localPosition = coord - cellPosition;
1833
1834 float4 vectorToClosest = {0.0f, 0.0f, 0.0f, 0.0f};
1835 float minDistance = FLT_MAX;
1836 for (int u = -1; u <= 1; u++) {
1837 for (int k = -1; k <= 1; k++) {
1838 for (int j = -1; j <= 1; j++) {
1839 for (int i = -1; i <= 1; i++) {
1840 float4 cellOffset(i, j, k, u);
1841 float4 vectorToPoint = cellOffset +
1842 hash_float_to_float4(cellPosition + cellOffset) *
1843 params.randomness -
1844 localPosition;
1845 float distanceToPoint = math::dot(vectorToPoint, vectorToPoint);
1846 if (distanceToPoint < minDistance) {
1847 minDistance = distanceToPoint;
1848 vectorToClosest = vectorToPoint;
1849 }
1850 }
1851 }
1852 }
1853 }
1854
1855 minDistance = FLT_MAX;
1856 for (int u = -1; u <= 1; u++) {
1857 for (int k = -1; k <= 1; k++) {
1858 for (int j = -1; j <= 1; j++) {
1859 for (int i = -1; i <= 1; i++) {
1860 float4 cellOffset(i, j, k, u);
1861 float4 vectorToPoint = cellOffset +
1862 hash_float_to_float4(cellPosition + cellOffset) *
1863 params.randomness -
1864 localPosition;
1865 float4 perpendicularToEdge = vectorToPoint - vectorToClosest;
1866 if (math::dot(perpendicularToEdge, perpendicularToEdge) > 0.0001f) {
1867 float distanceToEdge = math::dot((vectorToClosest + vectorToPoint) / 2.0f,
1868 math::normalize(perpendicularToEdge));
1869 minDistance = math::min(minDistance, distanceToEdge);
1870 }
1871 }
1872 }
1873 }
1874 }
1875
1876 return minDistance;
1877}
1878
1880{
1881 float4 cellPosition = math::floor(coord);
1882 float4 localPosition = coord - cellPosition;
1883
1884 float4 closestPoint = {0.0f, 0.0f, 0.0f, 0.0f};
1885 float4 closestPointOffset = {0.0f, 0.0f, 0.0f, 0.0f};
1886 float minDistance = FLT_MAX;
1887 for (int u = -1; u <= 1; u++) {
1888 for (int k = -1; k <= 1; k++) {
1889 for (int j = -1; j <= 1; j++) {
1890 for (int i = -1; i <= 1; i++) {
1891 float4 cellOffset(i, j, k, u);
1892 float4 pointPosition = cellOffset + hash_float_to_float4(cellPosition + cellOffset) *
1893 params.randomness;
1894 float distanceToPoint = math::distance(pointPosition, localPosition);
1895 if (distanceToPoint < minDistance) {
1896 minDistance = distanceToPoint;
1897 closestPoint = pointPosition;
1898 closestPointOffset = cellOffset;
1899 }
1900 }
1901 }
1902 }
1903 }
1904
1905 minDistance = FLT_MAX;
1906 float4 closestPointToClosestPoint = {0.0f, 0.0f, 0.0f, 0.0f};
1907 for (int u = -1; u <= 1; u++) {
1908 for (int k = -1; k <= 1; k++) {
1909 for (int j = -1; j <= 1; j++) {
1910 for (int i = -1; i <= 1; i++) {
1911 if (i == 0 && j == 0 && k == 0 && u == 0) {
1912 continue;
1913 }
1914 float4 cellOffset = float4(i, j, k, u) + closestPointOffset;
1915 float4 pointPosition = cellOffset + hash_float_to_float4(cellPosition + cellOffset) *
1916 params.randomness;
1917 float distanceToPoint = math::distance(closestPoint, pointPosition);
1918 if (distanceToPoint < minDistance) {
1919 minDistance = distanceToPoint;
1920 closestPointToClosestPoint = pointPosition;
1921 }
1922 }
1923 }
1924 }
1925 }
1926
1927 return math::distance(closestPointToClosestPoint, closestPoint) / 2.0f;
1928}
1929
1930/* **** Fractal Voronoi **** */
1931
1932/* The fractalization logic is the same as for fBM Noise, except that some additions are replaced
1933 * by lerps. */
1934template<typename T>
1936 const T coord,
1937 const bool calc_color /* Only used to optimize Smooth F1 */)
1938{
1939 float amplitude = 1.0f;
1940 float max_amplitude = 0.0f;
1941 float scale = 1.0f;
1942
1944 const bool zero_input = params.detail == 0.0f || params.roughness == 0.0f;
1945
1946 for (int i = 0; i <= ceilf(params.detail); ++i) {
1947 VoronoiOutput octave = (params.feature == NOISE_SHD_VORONOI_F2) ?
1948 voronoi_f2(params, coord * scale) :
1949 (params.feature == NOISE_SHD_VORONOI_SMOOTH_F1 &&
1950 params.smoothness != 0.0f) ?
1951 voronoi_smooth_f1(params, coord * scale, calc_color) :
1952 voronoi_f1(params, coord * scale);
1953
1954 if (zero_input) {
1955 max_amplitude = 1.0f;
1956 output = octave;
1957 break;
1958 }
1959 if (i <= params.detail) {
1960 max_amplitude += amplitude;
1961 output.distance += octave.distance * amplitude;
1962 output.color += octave.color * amplitude;
1963 output.position = mix(output.position, octave.position / scale, amplitude);
1964 scale *= params.lacunarity;
1965 amplitude *= params.roughness;
1966 }
1967 else {
1968 float remainder = params.detail - floorf(params.detail);
1969 if (remainder != 0.0f) {
1970 max_amplitude = mix(max_amplitude, max_amplitude + amplitude, remainder);
1971 output.distance = mix(
1972 output.distance, output.distance + octave.distance * amplitude, remainder);
1973 output.color = mix(output.color, output.color + octave.color * amplitude, remainder);
1974 output.position = mix(
1975 output.position, mix(output.position, octave.position / scale, amplitude), remainder);
1976 }
1977 }
1978 }
1979
1980 if (params.normalize) {
1981 output.distance /= max_amplitude * params.max_distance;
1982 output.color /= max_amplitude;
1983 }
1984
1985 output.position = (params.scale != 0.0f) ? output.position / params.scale :
1986 float4{0.0f, 0.0f, 0.0f, 0.0f};
1987
1988 return output;
1989}
1990
1991/* The fractalization logic is the same as for fBM Noise, except that some additions are replaced
1992 * by lerps. */
1993template<typename T>
1995{
1996 float amplitude = 1.0f;
1997 float max_amplitude = params.max_distance;
1998 float scale = 1.0f;
1999 float distance = 8.0f;
2000
2001 const bool zero_input = params.detail == 0.0f || params.roughness == 0.0f;
2002
2003 for (int i = 0; i <= ceilf(params.detail); ++i) {
2004 const float octave_distance = voronoi_distance_to_edge(params, coord * scale);
2005
2006 if (zero_input) {
2007 distance = octave_distance;
2008 break;
2009 }
2010 if (i <= params.detail) {
2011 max_amplitude = mix(max_amplitude, params.max_distance / scale, amplitude);
2012 distance = mix(distance, math::min(distance, octave_distance / scale), amplitude);
2013 scale *= params.lacunarity;
2014 amplitude *= params.roughness;
2015 }
2016 else {
2017 float remainder = params.detail - floorf(params.detail);
2018 if (remainder != 0.0f) {
2019 float lerp_amplitude = mix(max_amplitude, params.max_distance / scale, amplitude);
2020 max_amplitude = mix(max_amplitude, lerp_amplitude, remainder);
2021 float lerp_distance = mix(
2022 distance, math::min(distance, octave_distance / scale), amplitude);
2023 distance = mix(distance, math::min(distance, lerp_distance), remainder);
2024 }
2025 }
2026 }
2027
2028 if (params.normalize) {
2029 distance /= max_amplitude;
2030 }
2031
2032 return distance;
2033}
2034
2035/* Explicit function template instantiation */
2036
2038 const float coord,
2039 const bool calc_color);
2041 const float2 coord,
2042 const bool calc_color);
2044 const float3 coord,
2045 const bool calc_color);
2047 const float4 coord,
2048 const bool calc_color);
2049
2051 const float coord);
2053 const float2 coord);
2055 const float3 coord);
2057 const float4 coord);
2060/* -------------------------------------------------------------------- */
2080/* The original Gabor noise paper specifies that the impulses count for each cell should be
2081 * computed by sampling a Poisson distribution whose mean is the impulse density. However,
2082 * Tavernier's paper showed that stratified Poisson point sampling is better assuming the weights
2083 * are sampled using a Bernoulli distribution, as shown in Figure (3). By stratified sampling, they
2084 * mean a constant number of impulses per cell, so the stratification is the grid itself in that
2085 * sense, as described in the supplementary material of the paper. */
2086static constexpr int gabor_impulses_count = 8;
2087
2088/* Computes a 2D Gabor kernel based on Equation (6) in the original Gabor noise paper. Where the
2089 * frequency argument is the F_0 parameter and the orientation argument is the w_0 parameter. We
2090 * assume the Gaussian envelope has a unit magnitude, that is, K = 1. That is because we will
2091 * eventually normalize the final noise value to the unit range, so the multiplication by the
2092 * magnitude will be canceled by the normalization. Further, we also assume a unit Gaussian width,
2093 * that is, a = 1. That is because it does not provide much artistic control. It follows that the
2094 * Gaussian will be truncated at pi.
2095 *
2096 * To avoid the discontinuities caused by the aforementioned truncation, the Gaussian is windowed
2097 * using a Hann window, that is because contrary to the claim made in the original Gabor paper,
2098 * truncating the Gaussian produces significant artifacts especially when differentiated for bump
2099 * mapping. The Hann window is C1 continuous and has limited effect on the shape of the Gaussian,
2100 * so it felt like an appropriate choice.
2101 *
2102 * Finally, instead of computing the Gabor value directly, we instead use the complex phasor
2103 * formulation described in section 3.1.1 in Tricard's paper. That's done to be able to compute the
2104 * phase and intensity of the Gabor noise after summation based on equations (8) and (9). The
2105 * return value of the Gabor kernel function is then a complex number whose real value is the
2106 * value computed in the original Gabor noise paper, and whose imaginary part is the sine
2107 * counterpart of the real part, which is the only extra computation in the new formulation.
2108 *
2109 * Note that while the original Gabor noise paper uses the cosine part of the phasor, that is, the
2110 * real part of the phasor, we use the sine part instead, that is, the imaginary part of the
2111 * phasor, as suggested by Tavernier's paper in "Section 3.3. Instance stationarity and
2112 * normalization", to ensure a zero mean, which should help with normalization. */
2114 const float frequency,
2115 const float orientation)
2116{
2117 const float distance_squared = math::length_squared(position);
2118 const float hann_window = 0.5f + 0.5f * math::cos(math::numbers::pi * distance_squared);
2119 const float gaussian_envelop = math::exp(-math::numbers::pi * distance_squared);
2120 const float windowed_gaussian_envelope = gaussian_envelop * hann_window;
2121
2122 const float2 frequency_vector = frequency * float2(cos(orientation), sin(orientation));
2123 const float angle = 2.0f * math::numbers::pi * math::dot(position, frequency_vector);
2124 const float2 phasor = float2(math::cos(angle), math::sin(angle));
2125
2126 return windowed_gaussian_envelope * phasor;
2127}
2128
2129/* Computes the approximate standard deviation of the zero mean normal distribution representing
2130 * the amplitude distribution of the noise based on Equation (9) in the original Gabor noise paper.
2131 * For simplicity, the Hann window is ignored and the orientation is fixed since the variance is
2132 * orientation invariant. We start integrating the squared Gabor kernel with respect to x:
2133 *
2134 * \int_{-\infty}^{-\infty} (e^{- \pi (x^2 + y^2)} cos(2 \pi f_0 x))^2 dx
2135 *
2136 * Which gives:
2137 *
2138 * \frac{(e^{2 \pi f_0^2}-1) e^{-2 \pi y^2 - 2 pi f_0^2}}{2^\frac{3}{2}}
2139 *
2140 * Then we similarly integrate with respect to y to get:
2141 *
2142 * \frac{1 - e^{-2 \pi f_0^2}}{4}
2143 *
2144 * Secondly, we note that the second moment of the weights distribution is 0.5 since it is a
2145 * fair Bernoulli distribution. So the final standard deviation expression is square root the
2146 * integral multiplied by the impulse density multiplied by the second moment.
2147 *
2148 * Note however that the integral is almost constant for all frequencies larger than one, and
2149 * converges to an upper limit as the frequency approaches infinity, so we replace the expression
2150 * with the following limit:
2151 *
2152 * \lim_{x \to \infty} \frac{1 - e^{-2 \pi f_0^2}}{4}
2153 *
2154 * To get an approximation of 0.25. */
2156{
2157 const float integral_of_gabor_squared = 0.25f;
2158 const float second_moment = 0.5f;
2159 return math::sqrt(gabor_impulses_count * second_moment * integral_of_gabor_squared);
2160}
2161
2162/* Computes the Gabor noise value at the given position for the given cell. This is essentially the
2163 * sum in Equation (8) in the original Gabor noise paper, where we sum Gabor kernels sampled at a
2164 * random position with a random weight. The orientation of the kernel is constant for anisotropic
2165 * noise while it is random for isotropic noise. The original Gabor noise paper mentions that the
2166 * weights should be uniformly distributed in the [-1, 1] range, however, Tavernier's paper showed
2167 * that using a Bernoulli distribution yields better results, so that is what we do. */
2169 const float2 position,
2170 const float frequency,
2171 const float isotropy,
2172 const float base_orientation)
2173
2174{
2175 float2 noise(0.0f);
2176 for (const int i : IndexRange(gabor_impulses_count)) {
2177 /* Compute unique seeds for each of the needed random variables. */
2178 const float3 seed_for_orientation(cell.x, cell.y, i * 3);
2179 const float3 seed_for_kernel_center(cell.x, cell.y, i * 3 + 1);
2180 const float3 seed_for_weight(cell.x, cell.y, i * 3 + 2);
2181
2182 /* For isotropic noise, add a random orientation amount, while for anisotropic noise, use the
2183 * base orientation. Linearly interpolate between the two cases using the isotropy factor. Note
2184 * that the random orientation range spans pi as opposed to two pi, that's because the Gabor
2185 * kernel is symmetric around pi. */
2186 const float random_orientation = (noise::hash_float_to_float(seed_for_orientation) - 0.5f) *
2188 const float orientation = base_orientation + random_orientation * isotropy;
2189
2190 const float2 kernel_center = noise::hash_float_to_float2(seed_for_kernel_center);
2191 const float2 position_in_kernel_space = position - kernel_center;
2192
2193 /* The kernel is windowed beyond the unit distance, so early exit with a zero for points that
2194 * are further than a unit radius. */
2195 if (math::length_squared(position_in_kernel_space) >= 1.0f) {
2196 continue;
2197 }
2198
2199 /* We either add or subtract the Gabor kernel based on a Bernoulli distribution of equal
2200 * probability. */
2201 const float weight = noise::hash_float_to_float(seed_for_weight) < 0.5f ? -1.0f : 1.0f;
2202
2203 noise += weight * compute_2d_gabor_kernel(position_in_kernel_space, frequency, orientation);
2204 }
2205 return noise;
2206}
2207
2208/* Computes the Gabor noise value by dividing the space into a grid and evaluating the Gabor noise
2209 * in the space of each cell of the 3x3 cell neighborhood. */
2210static float2 compute_2d_gabor_noise(const float2 coordinates,
2211 const float frequency,
2212 const float isotropy,
2213 const float base_orientation)
2214{
2215 const float2 cell_position = math::floor(coordinates);
2216 const float2 local_position = coordinates - cell_position;
2217
2218 float2 sum(0.0f);
2219 for (int j = -1; j <= 1; j++) {
2220 for (int i = -1; i <= 1; i++) {
2221 const float2 cell_offset = float2(i, j);
2222 const float2 current_cell_position = cell_position + cell_offset;
2223 const float2 position_in_cell_space = local_position - cell_offset;
2225 current_cell_position, position_in_cell_space, frequency, isotropy, base_orientation);
2226 }
2227 }
2228
2229 return sum;
2230}
2231
2232/* Identical to compute_2d_gabor_kernel, except it is evaluated in 3D space. Notice that Equation
2233 * (6) in the original Gabor noise paper computes the frequency vector using (cos(w_0), sin(w_0)),
2234 * which we also do in the 2D variant, however, for 3D, the orientation is already a unit frequency
2235 * vector, so we just need to scale it by the frequency value. */
2237 const float frequency,
2238 const float3 orientation)
2239{
2240 const float distance_squared = math::length_squared(position);
2241 const float hann_window = 0.5f + 0.5f * math::cos(math::numbers::pi * distance_squared);
2242 const float gaussian_envelop = math::exp(-math::numbers::pi * distance_squared);
2243 const float windowed_gaussian_envelope = gaussian_envelop * hann_window;
2244
2245 const float3 frequency_vector = frequency * orientation;
2246 const float angle = 2.0f * math::numbers::pi * math::dot(position, frequency_vector);
2247 const float2 phasor = float2(math::cos(angle), math::sin(angle));
2248
2249 return windowed_gaussian_envelope * phasor;
2250}
2251
2252/* Identical to compute_2d_gabor_standard_deviation except we do triple integration in 3D. The only
2253 * difference is the denominator in the integral expression, which is 2^{5 / 2} for the 3D case
2254 * instead of 4 for the 2D case. Similarly, the limit evaluates to 1 / (4 * sqrt(2)). */
2256{
2257 const float integral_of_gabor_squared = 1.0f / (4.0f * math::numbers::sqrt2);
2258 const float second_moment = 0.5f;
2259 return math::sqrt(gabor_impulses_count * second_moment * integral_of_gabor_squared);
2260}
2261
2262/* Computes the orientation of the Gabor kernel such that it is constant for anisotropic
2263 * noise while it is random for isotropic noise. We randomize in spherical coordinates for a
2264 * uniform distribution. */
2265static float3 compute_3d_orientation(const float3 orientation,
2266 const float isotropy,
2267 const float4 seed)
2268{
2269 /* Return the base orientation in case we are completely anisotropic. */
2270 if (isotropy == 0.0) {
2271 return orientation;
2272 }
2273
2274 /* Compute the orientation in spherical coordinates. */
2275 float inclination = math::acos(orientation.z);
2276 float azimuth = math::sign(orientation.y) *
2277 math::acos(orientation.x / math::length(float2(orientation.x, orientation.y)));
2278
2279 /* For isotropic noise, add a random orientation amount, while for anisotropic noise, use the
2280 * base orientation. Linearly interpolate between the two cases using the isotropy factor. Note
2281 * that the random orientation range is to pi as opposed to two pi, that's because the Gabor
2282 * kernel is symmetric around pi. */
2284 inclination += random_angles.x * isotropy;
2285 azimuth += random_angles.y * isotropy;
2286
2287 /* Convert back to Cartesian coordinates, */
2288 return float3(math::sin(inclination) * math::cos(azimuth),
2289 math::sin(inclination) * math::sin(azimuth),
2290 math::cos(inclination));
2291}
2292
2294 const float3 position,
2295 const float frequency,
2296 const float isotropy,
2297 const float3 base_orientation)
2298
2299{
2300 float2 noise(0.0f);
2301 for (const int i : IndexRange(gabor_impulses_count)) {
2302 /* Compute unique seeds for each of the needed random variables. */
2303 const float4 seed_for_orientation(cell.x, cell.y, cell.z, i * 3);
2304 const float4 seed_for_kernel_center(cell.x, cell.y, cell.z, i * 3 + 1);
2305 const float4 seed_for_weight(cell.x, cell.y, cell.z, i * 3 + 2);
2306
2307 const float3 orientation = compute_3d_orientation(
2308 base_orientation, isotropy, seed_for_orientation);
2309
2310 const float3 kernel_center = noise::hash_float_to_float3(seed_for_kernel_center);
2311 const float3 position_in_kernel_space = position - kernel_center;
2312
2313 /* The kernel is windowed beyond the unit distance, so early exit with a zero for points that
2314 * are further than a unit radius. */
2315 if (math::length_squared(position_in_kernel_space) >= 1.0f) {
2316 continue;
2317 }
2318
2319 /* We either add or subtract the Gabor kernel based on a Bernoulli distribution of equal
2320 * probability. */
2321 const float weight = noise::hash_float_to_float(seed_for_weight) < 0.5f ? -1.0f : 1.0f;
2322
2323 noise += weight * compute_3d_gabor_kernel(position_in_kernel_space, frequency, orientation);
2324 }
2325 return noise;
2326}
2327
2328/* Identical to compute_2d_gabor_noise but works in the 3D neighborhood of the noise. */
2329static float2 compute_3d_gabor_noise(const float3 coordinates,
2330 const float frequency,
2331 const float isotropy,
2332 const float3 base_orientation)
2333{
2334 const float3 cell_position = math::floor(coordinates);
2335 const float3 local_position = coordinates - cell_position;
2336
2337 float2 sum(0.0f);
2338 for (int k = -1; k <= 1; k++) {
2339 for (int j = -1; j <= 1; j++) {
2340 for (int i = -1; i <= 1; i++) {
2341 const float3 cell_offset = float3(i, j, k);
2342 const float3 current_cell_position = cell_position + cell_offset;
2343 const float3 position_in_cell_space = local_position - cell_offset;
2345 current_cell_position, position_in_cell_space, frequency, isotropy, base_orientation);
2346 }
2347 }
2348 }
2349
2350 return sum;
2351}
2352
2353void gabor(const float2 coordinates,
2354 const float scale,
2355 const float frequency,
2356 const float anisotropy,
2357 const float orientation,
2358 float *r_value,
2359 float *r_phase,
2360 float *r_intensity)
2361{
2362 const float2 scaled_coordinates = coordinates * scale;
2363 const float isotropy = 1.0f - math::clamp(anisotropy, 0.0f, 1.0f);
2364 const float sanitized_frequency = math::max(0.001f, frequency);
2365
2366 const float2 phasor = compute_2d_gabor_noise(
2367 scaled_coordinates, sanitized_frequency, isotropy, orientation);
2368 const float standard_deviation = compute_2d_gabor_standard_deviation();
2369
2370 /* Normalize the noise by dividing by six times the standard deviation, which was determined
2371 * empirically. */
2372 const float normalization_factor = 6.0f * standard_deviation;
2373
2374 /* As discussed in compute_2d_gabor_kernel, we use the imaginary part of the phasor as the Gabor
2375 * value. But remap to [0, 1] from [-1, 1]. */
2376 if (r_value) {
2377 *r_value = (phasor.y / normalization_factor) * 0.5f + 0.5f;
2378 }
2379
2380 /* Compute the phase based on equation (9) in Tricard's paper. But remap the phase into the
2381 * [0, 1] range. */
2382 if (r_phase) {
2383 *r_phase = (math::atan2(phasor.y, phasor.x) + math::numbers::pi) / (2.0f * math::numbers::pi);
2384 }
2385
2386 /* Compute the intensity based on equation (8) in Tricard's paper. */
2387 if (r_intensity) {
2388 *r_intensity = math::length(phasor) / normalization_factor;
2389 }
2390}
2391
2392void gabor(const float3 coordinates,
2393 const float scale,
2394 const float frequency,
2395 const float anisotropy,
2396 const float3 orientation,
2397 float *r_value,
2398 float *r_phase,
2399 float *r_intensity)
2400{
2401 const float3 scaled_coordinates = coordinates * scale;
2402 const float isotropy = 1.0f - math::clamp(anisotropy, 0.0f, 1.0f);
2403 const float sanitized_frequency = math::max(0.001f, frequency);
2404
2405 const float3 normalized_orientation = math::normalize(orientation);
2406 const float2 phasor = compute_3d_gabor_noise(
2407 scaled_coordinates, sanitized_frequency, isotropy, normalized_orientation);
2408 const float standard_deviation = compute_3d_gabor_standard_deviation();
2409
2410 /* Normalize the noise by dividing by six times the standard deviation, which was determined
2411 * empirically. */
2412 const float normalization_factor = 6.0f * standard_deviation;
2413
2414 /* As discussed in compute_2d_gabor_kernel, we use the imaginary part of the phasor as the Gabor
2415 * value. But remap to [0, 1] from [-1, 1]. */
2416 if (r_value) {
2417 *r_value = (phasor.y / normalization_factor) * 0.5f + 0.5f;
2418 }
2419
2420 /* Compute the phase based on equation (9) in Tricard's paper. But remap the phase into the
2421 * [0, 1] range. */
2422 if (r_phase) {
2423 *r_phase = (math::atan2(phasor.y, phasor.x) + math::numbers::pi) / (2.0f * math::numbers::pi);
2424 }
2425
2426 /* Compute the intensity based on equation (8) in Tricard's paper. */
2427 if (r_intensity) {
2428 *r_intensity = math::length(phasor) / normalization_factor;
2429 }
2430}
2431
2434} // namespace blender::noise
#define BLI_assert_unreachable()
Definition BLI_assert.h:97
#define BLI_INLINE
#define ELEM(...)
#define X
#define Z
#define Y
ATTR_WARN_UNUSED_RESULT const BMVert * v2
ATTR_WARN_UNUSED_RESULT const BMVert * v
#define output
SIMD_FORCE_INLINE const btScalar & z() const
Return the z value.
Definition btQuadWord.h:117
SIMD_FORCE_INLINE const btScalar & w() const
Return the w value.
Definition btQuadWord.h:119
static T sum(const btAlignedObjectArray< T > &items)
static unsigned long seed
Definition btSoftBody.h:39
SIMD_FORCE_INLINE btVector3 & normalize()
Normalize this vector x^2 + y^2 + z^2 = 1.
Definition btVector3.h:303
local_group_size(16, 16) .push_constant(Type b
#define ceilf(x)
#define floorf(x)
#define fabsf(x)
draw_view in_light_buf[] float
draw_view push_constant(Type::INT, "radiance_src") .push_constant(Type capture_info_buf storage_buf(1, Qualifier::READ, "ObjectBounds", "bounds_buf[]") .push_constant(Type draw_view int
#define mix(a, b, c)
Definition hash.h:36
uiWidgetBaseParameters params[MAX_WIDGET_BASE_BATCH]
MINLINE float smoothstep(float edge0, float edge1, float x)
ccl_device_inline float3 cos(float3 v)
T length_squared(const VecBase< T, Size > &a)
T cos(const AngleRadianBase< T > &a)
T clamp(const T &a, const T &min, const T &max)
T acos(const T &a)
T sqrt(const T &a)
T sign(const T &a)
T floor(const T &a)
T distance(const T &a, const T &b)
T length(const VecBase< T, Size > &a)
T exp(const T &x)
T dot(const QuaternionBase< T > &a, const QuaternionBase< T > &b)
T min(const T &a, const T &b)
T atan2(const T &y, const T &x)
MatBase< T, NumCol, NumRow > normalize(const MatBase< T, NumCol, NumRow > &a)
T sin(const AngleRadianBase< T > &a)
T max(const T &a, const T &b)
T mod(const T &a, const T &b)
T abs(const T &a)
template VoronoiOutput fractal_voronoi_x_fx< float4 >(const VoronoiParams &params, const float4 coord, const bool calc_color)
template float fractal_voronoi_distance_to_edge< float >(const VoronoiParams &params, const float coord)
float3 perlin_float3_fractal_distorted(float position, float detail, float roughness, float lacunarity, float offset, float gain, float distortion, int type, bool normalize)
Definition noise.cc:886
BLI_INLINE float negate_if(float value, uint32_t condition)
Definition noise.cc:364
BLI_INLINE float uint_to_float_01(uint32_t k)
Definition noise.cc:161
float perlin_ridged_multi_fractal(T p, const float detail, const float roughness, const float lacunarity, const float offset, const float gain)
Definition noise.cc:701
static float2 compute_3d_gabor_noise(const float3 coordinates, const float frequency, const float isotropy, const float3 base_orientation)
Definition noise.cc:2329
uint32_t hash_float(float kx)
Definition noise.cc:134
static float compute_3d_gabor_standard_deviation()
Definition noise.cc:2255
template float perlin_fbm< float3 >(float3 p, const float detail, const float roughness, const float lacunarity, const bool normalize)
static float2 compute_2d_gabor_noise_cell(const float2 cell, const float2 position, const float frequency, const float isotropy, const float base_orientation)
Definition noise.cc:2168
BLI_INLINE float perlin_distortion(float position, float strength)
Definition noise.cc:803
float hash_float_to_float(float k)
Definition noise.cc:188
float perlin_signed(float position)
Definition noise.cc:511
void gabor(const float2 coordinates, const float scale, const float frequency, const float anisotropy, const float orientation, float *r_value, float *r_phase, float *r_intensity)
Definition noise.cc:2353
@ NOISE_SHD_VORONOI_N_SPHERE_RADIUS
Definition noise.cc:1046
@ NOISE_SHD_VORONOI_SMOOTH_F1
Definition noise.cc:1044
@ NOISE_SHD_VORONOI_F2
Definition noise.cc:1043
@ NOISE_SHD_VORONOI_DISTANCE_TO_EDGE
Definition noise.cc:1045
@ NOISE_SHD_VORONOI_F1
Definition noise.cc:1042
template float perlin_fractal_distorted< float >(float position, float detail, float roughness, float lacunarity, float offset, float gain, float distortion, int type, bool normalize)
@ NOISE_SHD_VORONOI_MANHATTAN
Definition noise.cc:1036
@ NOISE_SHD_VORONOI_MINKOWSKI
Definition noise.cc:1038
@ NOISE_SHD_VORONOI_EUCLIDEAN
Definition noise.cc:1035
@ NOISE_SHD_VORONOI_CHEBYCHEV
Definition noise.cc:1037
BLI_INLINE void hash_bit_mix(uint32_t &a, uint32_t &b, uint32_t &c)
Definition noise.cc:32
float4 voronoi_position(const float coord)
Definition noise.cc:1123
BLI_INLINE float fade(float t)
Definition noise.cc:359
VoronoiOutput voronoi_smooth_f1(const VoronoiParams &params, const float coord, const bool calc_color)
Definition noise.cc:1155
float voronoi_n_sphere_radius(const VoronoiParams &params, const float coord)
Definition noise.cc:1245
static float2 compute_2d_gabor_noise(const float2 coordinates, const float frequency, const float isotropy, const float base_orientation)
Definition noise.cc:2210
float perlin(float position)
Definition noise.cc:562
BLI_INLINE float random_float_offset(float seed)
Definition noise.cc:775
BLI_INLINE float floor_fraction(float x, int &i)
Definition noise.cc:402
float voronoi_distance_to_edge(const VoronoiParams &params, const float coord)
Definition noise.cc:1231
BLI_INLINE uint32_t float_as_uint(float f)
Definition noise.cc:124
template VoronoiOutput fractal_voronoi_x_fx< float3 >(const VoronoiParams &params, const float3 coord, const bool calc_color)
template float perlin_fractal_distorted< float2 >(float2 position, float detail, float roughness, float lacunarity, float offset, float gain, float distortion, int type, bool normalize)
BLI_INLINE float3 random_float3_offset(float seed)
Definition noise.cc:786
static float2 compute_3d_gabor_noise_cell(const float3 cell, const float3 position, const float frequency, const float isotropy, const float3 base_orientation)
Definition noise.cc:2293
static float2 compute_2d_gabor_kernel(const float2 position, const float frequency, const float orientation)
Definition noise.cc:2113
VoronoiOutput fractal_voronoi_x_fx(const VoronoiParams &params, const T coord, const bool calc_color)
Definition noise.cc:1935
float voronoi_distance(const float a, const float b)
Definition noise.cc:1051
static float2 compute_3d_gabor_kernel(const float3 position, const float frequency, const float3 orientation)
Definition noise.cc:2236
template float perlin_fractal_distorted< float4 >(float4 position, float detail, float roughness, float lacunarity, float offset, float gain, float distortion, int type, bool normalize)
BLI_INLINE float2 random_float2_offset(float seed)
Definition noise.cc:780
float perlin_fbm(T p, float detail, float roughness, float lacunarity, bool normalize)
Definition noise.cc:586
VoronoiOutput voronoi_f2(const VoronoiParams &params, const float coord)
Definition noise.cc:1193
float3 hash_float_to_float3(float k)
Definition noise.cc:225
template VoronoiOutput fractal_voronoi_x_fx< float >(const VoronoiParams &params, const float coord, const bool calc_color)
BLI_INLINE float4 random_float4_offset(float seed)
Definition noise.cc:793
BLI_INLINE float perlin_noise(float position)
Definition noise.cc:409
float fractal_voronoi_distance_to_edge(const VoronoiParams &params, const T coord)
Definition noise.cc:1994
BLI_INLINE void hash_bit_final(uint32_t &a, uint32_t &b, uint32_t &c)
Definition noise.cc:54
float4 hash_float_to_float4(float4 k)
Definition noise.cc:253
float perlin_hybrid_multi_fractal(T p, const float detail, const float roughness, const float lacunarity, const float offset, const float gain)
Definition noise.cc:665
static float compute_2d_gabor_standard_deviation()
Definition noise.cc:2155
@ NOISE_SHD_PERLIN_MULTIFRACTAL
Definition noise.cc:729
@ NOISE_SHD_PERLIN_FBM
Definition noise.cc:730
@ NOISE_SHD_PERLIN_RIDGED_MULTIFRACTAL
Definition noise.cc:732
@ NOISE_SHD_PERLIN_HYBRID_MULTIFRACTAL
Definition noise.cc:731
@ NOISE_SHD_PERLIN_HETERO_TERRAIN
Definition noise.cc:733
float perlin_select(T p, float detail, float roughness, float lacunarity, float offset, float gain, int type, bool normalize)
Definition noise.cc:737
BLI_INLINE float noise_grad(uint32_t hash, float x)
Definition noise.cc:369
VoronoiOutput voronoi_f1(const VoronoiParams &params, const float coord)
Definition noise.cc:1128
float perlin_multi_fractal(T p, const float detail, const float roughness, const float lacunarity)
Definition noise.cc:619
template float fractal_voronoi_distance_to_edge< float3 >(const VoronoiParams &params, const float3 coord)
static constexpr int gabor_impulses_count
Definition noise.cc:2086
float perlin_fractal_distorted(T position, float detail, float roughness, float lacunarity, float offset, float gain, float distortion, int type, bool normalize)
Definition noise.cc:832
template float perlin_fractal_distorted< float3 >(float3 position, float detail, float roughness, float lacunarity, float offset, float gain, float distortion, int type, bool normalize)
template VoronoiOutput fractal_voronoi_x_fx< float2 >(const VoronoiParams &params, const float2 coord, const bool calc_color)
template float fractal_voronoi_distance_to_edge< float4 >(const VoronoiParams &params, const float4 coord)
BLI_INLINE uint32_t hash_bit_rotate(uint32_t x, uint32_t k)
Definition noise.cc:27
float perlin_hetero_terrain(T p, const float detail, const float roughness, const float lacunarity, const float offset)
Definition noise.cc:639
float2 hash_float_to_float2(float2 k)
Definition noise.cc:208
static float3 compute_3d_orientation(const float3 orientation, const float isotropy, const float4 seed)
Definition noise.cc:2265
float hash_to_float(uint32_t kx)
Definition noise.cc:166
template float fractal_voronoi_distance_to_edge< float2 >(const VoronoiParams &params, const float2 coord)
VecBase< float, 4 > float4
VecBase< float, 2 > float2
VecBase< float, 3 > float3
float distance(float a, float b)
#define hash
Definition noise.c:154
#define FLT_MAX
Definition stdcycles.h:14
unsigned int uint32_t
Definition stdint.h:80