perlin.fut

import "lib/github.com/diku-dk/cpprandom/random"

-- We'll need a random number generator. We don't need any good guarantees, 
-- so we use the fastest one.

module dist = uniform_real_distribution f32 minstd_rand

-- Let's define some aliases for readability. 

let split = minstd_rand.split_rng
let from_seed x = minstd_rand.rng_from_seed [x]
let sample rng = dist.rand (0,1) rng

-- Now we can generate a random grid!

let randGrid (x: i64) (y: i64) (rng: minstd_rand.rng): [y][x]f32 = 
  let y_rngs = split y rng
  let grid_rngs = map (\rng -> split x rng) y_rngs
  in map (map (\rng -> (sample rng).1)) grid_rngs

-- > :img randGrid 100i64 100i64 (from_seed 10)

-- However, we probably want something nicer. First, some math utilities.

type vec = (f32, f32)
let v0: vec = (0.0, 0.0)

let frac (x: f32): f32 = x - f32.floor x
let floori (x: f32): i64 = i64.f32 (f32.floor x)
let dot (x: vec) (y: vec): f32 = x.0 * y.0 + x.1 * y.1
let sub (x: vec) (y: vec): vec = (x.0 - y.0, x.1 - y.1)
let len (x: vec): f32 = f32.sqrt (x.0 * x.0 + x.1 * x.1)
let norm (x: vec): vec = (x.0 / len x, x.1 / len x)

-- This perlin noise function returns a closure which takes a pair of values 
-- between 0 and 1, and gives the value of the perlin noise.

let interpolate (a0: f32) (a1: f32) (w: f32) = (a1 - a0) * (3.0 - w * 2.0) * w * w + a0

type perlin_noise [n] = {
  vecs: [n][n]vec,
  subdivs: i64
}

-- This function allows us to store an array with element type perlin_noise. 
-- To do that, we need to have them all be the same size. Expand pads the perlin
-- noise arrays with zero vectors.

-- The type system is not expressive enough to track the size changes to nested
-- arrays in this function. Therefore, we have to use dynamic size assertions
let expand [n] (to: i64) (from: perlin_noise [n]): perlin_noise [to] = 
  let extended_rows: [n][to]vec = map (\arr -> concat arr (replicate (n - to) v0) :> [to]vec) from.vecs
  let extra_rows = replicate to (replicate to v0)
  in {subdivs = n, vecs = concat extended_rows extra_rows :> [to][to]vec}

let gen_perlin_noise (subdivs: i64) (rng: minstd_rand.rng): perlin_noise [subdivs] = 
  let rngs = split 2 rng
  let xss = randGrid subdivs subdivs rngs[0]
  let yss = randGrid subdivs subdivs rngs[1]
  -- We need to zip the x components of the vectors with the y components, and while
  -- we do this, we should also shift and normalize the vectors.
  let vecs = map 
    (\(xs, ys) -> map (\(x,y) -> norm (x-0.5, y-0.5)) (zip xs ys))
    (zip xss yss)
  in { vecs, subdivs }

let sample_perlin_noise [n] (noise: perlin_noise [n]) (pt: (f32, f32)): f32 = 
    let x = pt.0 * f32.i64 (noise.subdivs - 1)
    let y = pt.1 * f32.i64 (noise.subdivs - 1)
    let pts = [
      (f32.floor x, f32.floor y),
      (f32.floor x, f32.floor y + 1),
      (f32.floor x + 1, f32.floor y),
      (f32.floor x + 1, f32.floor y + 1)
    ]
    let pt_vecs = [
      noise.vecs[floori x, floori y],
      noise.vecs[floori x, floori y + 1], 
      noise.vecs[floori x + 1, floori y],
      noise.vecs[floori x + 1, floori y + 1]
    ]
    let dxs = map (sub (x, y)) pts
    let dots = map2 dot pt_vecs dxs
    let x0 = interpolate dots[0] dots[1] (frac y)
    let x1 = interpolate dots[2] dots[3] (frac y)
    in (1.0 + interpolate x0 x1 (frac x)) / 2.0

let sample_grid (x: i64) (y: i64) (f: (f32, f32) -> f32): [y][x]f32 = 
  map (\yi -> 
    map (\xi -> 
      f (f32.i64 xi / f32.i64 x, f32.i64 yi / f32.i64 y)
    ) 
    (iota x)
  ) 
  (iota y)

let example = sample_grid 500i64 500i64 (sample_perlin_noise (gen_perlin_noise 20 (from_seed 12)))

-- > :img example

-- Now we can combine perlin noise at different scales to get the rough looking
-- perlin noise we're all used to!

type perlin_octaves [o] [n] = {
  octaves: [o](perlin_noise [n])
}

let gen_perlin_octaves (start_amplitude: i64) (num_amplitudes: i64) (rng: minstd_rand.rng) 
  : perlin_octaves [num_amplitudes] [] = 
  let maxsize = start_amplitude * 2 ** (num_amplitudes -1)
  let rngs = split num_amplitudes rng
  let vecs: [num_amplitudes](perlin_noise [maxsize]) = map2 
    (\freq rng -> 
      expand maxsize (gen_perlin_noise freq rng)
    )
    (map (\i -> start_amplitude * 2 ** i) (iota num_amplitudes))
    rngs
  in {octaves = vecs}

let sample_perlin_octaves [n] [m] (octaves: perlin_octaves [n] [m]) (persistence: f32) (pt: (f32, f32)): f32 = 
  (loop (r, m, p) = (0.0, 0.0, 1.0) for perlin_noise in octaves.octaves do
    (r + p * sample_perlin_noise perlin_noise pt, r + p, p * persistence)).0

let example2 = sample_grid 500i64 500i64 (sample_perlin_octaves (gen_perlin_octaves 2 5 (from_seed 12)) 1.0)

-- > :img example2