From 4d3be19898dccae8b216b888b4ef1418f1f2483c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Johannes=20Dahlstr=C3=B6m?= Date: Tue, 27 Aug 2024 14:52:25 +0300 Subject: [PATCH] WIP Perlin noise support --- core/src/math.rs | 1 + core/src/math/noise.rs | 256 +++++++++++++++++++++++++++++++++++++++++ 2 files changed, 257 insertions(+) create mode 100644 core/src/math/noise.rs diff --git a/core/src/math.rs b/core/src/math.rs index c3e41915..962b011e 100644 --- a/core/src/math.rs +++ b/core/src/math.rs @@ -30,6 +30,7 @@ pub mod approx; pub mod color; pub mod float; pub mod mat; +pub mod noise; pub mod rand; pub mod space; pub mod spline; diff --git a/core/src/math/noise.rs b/core/src/math/noise.rs new file mode 100644 index 00000000..9becd7de --- /dev/null +++ b/core/src/math/noise.rs @@ -0,0 +1,256 @@ +//! Procedural noise generation. + +use super::{ + spline::smoothstep, + vary::lerp, + vec::{splat, vec2, vec3, Vec2, Vec3, Vector}, +}; + +/// Returns the Perlin noise value corresponding to the 2D point. +pub fn perlin2(pt: Vec2) -> f32 { + use super::float::f32; + let f = pt - pt.map(f32::floor); + let [fx, fy] = f.0; + + let [g00, g01, g10, g11] = grads2(pt); + + let d00 = g00.dot(&vec2(fx, fy)); + let d01 = g01.dot(&vec2(fx, fy - 1.0)); + let d10 = g10.dot(&vec2(fx - 1.0, fy)); + let d11 = g11.dot(&vec2(fx - 1.0, fy - 1.0)); + + let t = f.map(smoothstep); + let (a, b) = lerp(t.x(), (d00, d01), (d10, d11)); + lerp(t.y(), a, b) +} + +/// Returns the Perlin gradient vector corresponding to the 2D point. +/// +/// The vector is computed by taking the gradient vectors of all four +/// surrounding grid points and smoothly interpolating between them. +pub fn perlin2v(pt: Vec2) -> Vec2 { + use super::float::f32; + let fract = pt - pt.map(f32::floor); + + let [g00, g01, g10, g11] = grads2(pt); + + let t = fract.map(smoothstep); + let (a, b) = lerp(t.x(), (g00, g01), (g10, g11)); + lerp(t.y(), a, b) +} + +/// Returns the Perlin noise value corresponding to the 3D point. +pub fn perlin3(pt: Vec3) -> f32 { + use super::float::f32; + let pt0 = pt.map(f32::floor); + let f = pt - pt0; + let i_f: Vec3 = f - splat(1.0); + + let [fx, fy, fz] = f.0; + let [ifx, ify, ifz] = i_f.0; + + let [g000, g001, g010, g011, g100, g101, g110, g111] = grads3(pt0); + + let d000 = g000.dot(&vec3(fx, fy, fz)); + let d001 = g001.dot(&vec3(fx, fy, ifz)); + let d010 = g010.dot(&vec3(fx, ify, fz)); + let d011 = g011.dot(&vec3(fx, ify, ifz)); + let d100 = g100.dot(&vec3(ifx, fy, fz)); + let d101 = g101.dot(&vec3(ifx, fy, ifz)); + let d110 = g110.dot(&vec3(ifx, ify, fz)); + let d111 = g111.dot(&vec3(ifx, ify, ifz)); + + let t: Vec3 = f.map(smoothstep); + + let x0 = Vector::<_>::new([d000, d010, d001, d011]); + let x1 = Vector::<_>::new([d100, d110, d101, d111]); + let [y0, y1, y2, y3] = lerp(t.x(), x0, x1).0; + let (z0, z1) = lerp(t.y(), (y0, y2), (y1, y3)); + lerp(t.z(), z0, z1) +} + +/// Returns the Perlin gradient vector corresponding to the 3D point. +pub fn perlin3v(pt: Vec3) -> Vec3 { + use super::float::f32; + let fract = pt - pt.map(f32::floor); + + let [g000, g001, g010, g011, g100, g101, g110, g111] = grads3(pt); + + let t: Vec3 = fract.map(smoothstep); + + let x0 = ((g000, g010), (g001, g011)); + let x1 = ((g100, g110), (g101, g111)); + let ((y0, y1), (y2, y3)) = lerp(t.x(), x0, x1); + let (z0, z1) = lerp(t.y(), (y0, y2), (y1, y3)); + lerp(t.z(), z0, z1) +} + +fn perm(x: i32) -> i32 { + PERM[(x & 0xFF) as usize] as i32 +} + +fn grad2(x: f32, y: f32) -> Vec2 { + let perm = perm(x as i32 + perm(y as i32)); + GRADS_2[perm as usize & 0x7] +} + +fn grads2(pt: Vec2) -> [Vec2; 4] { + use super::float::f32; + let [x0, y0] = pt.map(f32::floor).0; + + let g00 = grad2(x0, y0); + let g01 = grad2(x0, y0 + 1.0); + let g10 = grad2(x0 + 1.0, y0); + let g11 = grad2(x0 + 1.0, y0 + 1.0); + + [g00, g01, g10, g11] +} + +fn grad3(x: f32, y: f32, z: f32) -> Vec3 { + let perm = perm(x as i32 + perm(y as i32 + perm(z as i32))); + GRADS_3[perm as usize & 0xF] +} + +fn grads3(pt0: Vec3) -> [Vec3; 8] { + let pt1 = pt0 + splat(1.0); + let [x0, y0, z0] = pt0.0; + let [x1, y1, z1] = pt1.0; + + // + // 011 +--------------+ 111 + // / | / | + // / | / | + // 010 +--------------+ 110 | + // | | | | + // | +--------|-----+ 101 + // | / | / + // | / | / + // 000 +--------------+ 100 + // + + let g000 = grad3(x0, y0, z0); + let g001 = grad3(x0, y0, z1); + let g010 = grad3(x0, y1, z0); + let g011 = grad3(x0, y1, z1); + let g100 = grad3(x1, y0, z0); + let g101 = grad3(x1, y0, z1); + let g110 = grad3(x1, y1, z0); + let g111 = grad3(x1, y1, z1); + + [g000, g001, g010, g011, g100, g101, g110, g111] +} + +const A: f32 = 1.237; +const B: f32 = 0.513; +/// Gradient vectors for 2D noise. +static GRADS_2: [Vec2; 8] = [ + vec2(A, B), + vec2(B, A), + vec2(-B, A), + vec2(-A, B), + vec2(-A, -B), + vec2(-B, -A), + vec2(B, -A), + vec2(A, -B), +]; + +/// Gradient vectors for 3D noise. +static GRADS_3: [Vec3; 16] = [ + // + vec3(0.0, -1.0, -1.0), + vec3(0.0, -1.0, 1.0), + vec3(0.0, 1.0, -1.0), + vec3(0.0, 1.0, 1.0), + // + vec3(-1.0, 0.0, -1.0), + vec3(-1.0, 0.0, 1.0), + vec3(1.0, 0.0, -1.0), + vec3(1.0, 0.0, 1.0), + // + vec3(-1.0, -1.0, 0.0), + vec3(-1.0, 1.0, 0.0), + vec3(1.0, -1.0, 0.0), + vec3(1.0, 1.0, 0.0), + // + vec3(1.0, 1.0, 0.0), + vec3(-1.0, 1.0, 0.0), + vec3(0.0, -1.0, 1.0), + vec3(0.0, -1.0, -1.0), +]; + +/// Permutation table for calculating a pseudo-random index for each grid point. +#[rustfmt::skip] +static PERM: [u8; 256] = [ + 156, 2, 157, 90, 75, 199, 55, 167, 62, 92, 101, 253, 66, 134, 113, 83, + 1, 136, 78, 106, 254, 105, 248, 176, 234, 5, 195, 226, 49, 71, 87, 44, + 122, 94, 219, 140, 72, 159, 237, 212, 8, 162, 200, 124, 125, 69, 165, 74, + 245, 42, 89, 216, 158, 108, 238, 184, 217, 73, 126, 210, 14, 111, 19, 188, + 186, 45, 38, 223, 35, 112, 214, 26, 145, 95, 99, 193, 250, 189, 152, 182, + 166, 247, 148, 213, 168, 70, 96, 249, 127, 132, 4, 137, 41, 60, 102, 28, + 27, 240, 227, 155, 211, 230, 9, 80, 178, 3, 68, 153, 143, 84, 179, 181, + 12, 97, 103, 16, 225, 146, 63, 82, 203, 175, 163, 147, 11, 116, 185, 215, + 57, 120, 208, 129, 115, 198, 37, 201, 39, 98, 20, 183, 56, 118, 109, 142, + 138, 65, 117, 114, 160, 25, 43, 191, 204, 161, 22, 251, 139, 79, 131, 231, + 76, 0, 205, 206, 244, 51, 174, 13, 110, 85, 209, 77, 64, 53, 48, 221, + 133, 93, 224, 24, 33, 164, 23, 47, 171, 128, 243, 18, 52, 119, 149, 100, + 246, 233, 31, 192, 252, 190, 15, 172, 91, 229, 144, 54, 61, 58, 220, 36, + 222, 29, 50, 88, 121, 173, 232, 194, 239, 197, 32, 180, 107, 46, 7, 130, + 169, 81, 218, 67, 21, 170, 187, 59, 86, 235, 154, 123, 150, 177, 135, 228, + 104, 242, 6, 151, 255, 34, 30, 141, 202, 196, 236, 207, 241, 40, 17, 10 +]; + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn perlin2_statistics() { + let count = 1000u32; + let div = 10.0; + let [mut avg, mut min, mut max, mut var] = [0.0f32, 1.0, -1.0, 0.0]; + for i in 0..count { + for j in 0..count { + let v = perlin2(vec2(i as f32, j as f32) / div); + avg += v; + var += v * v; + min = min.min(v); + max = max.max(v); + } + } + let c = count.pow(2) as f32; + avg /= c; + var /= c; + + assert_eq!(avg, -0.0005586566); + assert_eq!(min, -0.8853002); + assert_eq!(max, 0.8853002); + assert_eq!(var, 0.07994927); + } + #[test] + fn perlin3_statistics() { + let count = 100u32; + let div = 10.0; + let [mut avg, mut min, mut max, mut var] = [0.0f32, 1.0, -1.0, 0.0]; + for i in 0..count { + for j in 0..count { + for k in 0..count { + let pt = vec3(i, j, k).map(|c| c as f32) / div; + let v = perlin3(pt); + + avg += v; + var += v * v; + min = min.min(v); + max = max.max(v); + } + } + } + let c = count.pow(3) as f32; + avg /= c; + var /= c; + + assert_eq!(avg, 0.00039595732); + assert_eq!(min, -0.8853568); + assert_eq!(max, 0.8781377); + assert_eq!(var, 0.06190137); + } +}