code / hill-noise

📦 Download .tar.gz
Lines837 Lua675 GLSL103 Markdown59
(541 lines)
Raw
   1 -- Hill Noise library

   2 -- Copyright 2017 Bruce Hill

   3 --

   4 -- This library provides functions for continuous pseudorandom 1D, 2D, and 3D noise functions

   5 -- with approximately uniform distribution on the interval (0,1) and customizable resolution,

   6 -- levels of detail, and shape characteristics.

   7 --

   8 -- Functions: make1d([opts]), make2d([opts]), make3d([opts])

   9 --   each of these returns a noise function, and a gradient function of the noise function

  10 -- Shader functions: make1dShader([opts]), make2dShader([opts]), make3dShader([opts])

  11 --   each of these returns a noise shader, and a gradient shader of the noise function

  12 --

  13 -- Usage:

  14 --     local Noise = require 'hill_noise'

  15 --     local noise,noise_gradient = Noise.make3d{resolution=9}

  16 --     local x,y,z = 1,2,3

  17 --     local n = noise(x,y,z)

  18 --     local dndx,dndy,dndz = noise_gradient(x,y,z)

  19 --

  20 -- The noise construction functions all take the following options:

  21 --     * random: a random number function (default: math.random)

  22 --     * seed: if "random" is not provided, and the love.math module is loaded, a new

  23 --       pseudorandom number generator will be used with the specified seed

  24 --     * amplitudes: a list of sine wave amplitudes to use

  25 --     * resolution: if "amplitudes" is not provided, the number of sine waves to use

  26 --       (default: 5,7,9 for 1d,2d,3d functions, and 10,15,20 for 1d,2d,3d shaders)

  27 --     * distribution: if "amplitudes" is not provided, a function to evenly sample

  28 --       amplitudes from on the interval (0,1) (default: 0.1^x)

  29 --

  30 -- The noise and noise gradient functions all run in O(resolution), and more resolution

  31 -- may be needed for higher dimensional noise and additional detail.

  32 --

  33 -- For the LOVE game engine, there are versions that produce 1D, 2D, and 3D noise shaders

  34 -- as well, which are much faster and run on the GPU.

  35 --     local noise_shader, noise_gradient_shader = Noise.make3dShader{resolution=11}

  36 --     local canvas = love.graphics.newCanvas()

  37 --     -- At each location on the texture, the noise function is sampled at a linear

  38 --     --   interpolation between range_min and range_max, using the texture coordinates.

  39 --     --   For 3D noise, the "z" value is passed as a parameter.

  40 --     noise_shader:send("range_min", {0,20})

  41 --     noise_shader:send("range_max", {0,20})

  42 --     noise_shader:send('z', love.timer.getTime())

  43 --     love.graphics.setShader(noise_shader)

  44 --     love.graphics.draw(canvas)

  45 --     love.graphics.setShader()

  46 -- The noise gradient shaders populate the RGB channels with:

  47 --   * for 1D: sigmoid(2*dn/dx), sigmoid(2*dn/dx), sigmoid(2*dn/dx)

  48 --   * for 2D: sigmoid(1.5*dn/dx), sigmoid(1.5*dn/dy), 0

  49 --   * for 3D: sigmoid(4*dn/dx), sigmoid(4*dn/dy), sigmoid(4*dn/dz)

  50 -- (coefficients were arbitrarily chosen to best squash typical gradients for the default

  51 --  parameters evenly into (0,1))

  52 -- Additionally, the 1D noise shader has an extra option "draw_outline" that, if set to

  53 -- true, makes the 1D shader render output as white or black if the y-coordinate is

  54 -- above/below the  noise value for a given x-coordinate.

  55 -- 

  56 -- Permission is hereby granted, free of charge, to any person obtaining a copy of

  57 -- this software and associated documentation files (the "Software"), to deal in

  58 -- the Software without restriction, including without limitation the rights to

  59 -- use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies

  60 -- of the Software, and to permit persons to whom the Software is furnished to do

  61 -- so, subject to the following conditions:

  62 -- 

  63 -- The above copyright notice and this permission notice shall be included in all

  64 -- copies or substantial portions of the Software.

  65 -- 

  66 -- THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR

  67 -- IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,

  68 -- FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE

  69 -- AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER

  70 -- LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,

  71 -- OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE

  72 -- SOFTWARE.

  74 local noise = {}

  75 local exp,sin,cos,floor,log,acos,sqrt,abs = math.exp,math.sin,math.cos,math.floor,math.log,math.acos,math.sqrt,math.abs

  76 local GR, PI, TAU, SQRT5, LOG_GR = (sqrt(5)+1)/2, math.pi, 2*math.pi, sqrt(5), log((sqrt(5)+1)/2)

  78 local function cdf(x)

  79     return .5 + .5*(x<0 and -1 or 1)*sqrt(1.-exp(-4./TAU * x*x))

  80 end

  82 local function cdf_prime(x, dx)

  83     return (0.31831 * exp(-2/PI * x*x) * abs(x)*dx)/sqrt(1-exp(-2/PI*x*x))

  84 end

  86 local default_resolutions = {["1d"]=5, ["2d"]=7, ["3d"]=9, ["1dShader"]=10, ["2dShader"]=15, ["3dShader"]=20}

  87 local function _defaultArgs(opts, kind)

  88     opts = opts or {}

  89     local amplitudes, random = opts.amplitudes, opts.random

  90     if random then

  91         assert(type(random) == 'function', "Random number function must be a function.")

  92     elseif not random then

  93         if opts.seed and love and love.math then

  94             local rng = love.math.newRandomGenerator(opts.seed)

  95             random = function(...) return rng:random(...) end

  96         else

  97             random = math.random

  98         end

  99     end

 100     if not amplitudes then

 101         local resolution = opts.resolution or default_resolutions[kind]

 102         local distribution = opts.distribution or (function(x) return .1^x end)

 103         amplitudes = {}

 104         for i=1,resolution do

 105             amplitudes[i] = distribution((i-.5)/resolution)

 106         end

 107     end

 108     return amplitudes, random

 109 end

 111 local function calculateSigma(amplitudes)

 112     local sigma = 0

 113     for _,a in ipairs(amplitudes) do

 114         sigma = sigma + a*a

 115     end

 116     return math.sqrt(sigma/2)

 117 end

 119 local function makeOffsets(count, random)

 120     local offsets = {}

 121     for i=1,count do

 122         offsets[i] = random()*TAU

 123     end

 124     return offsets

 125 end

 127 noise.make1d = function(opts)

 128     local amplitudes, random = _defaultArgs(opts, "1d")

 129     local sigma = calculateSigma(amplitudes)

 130     local offsets = makeOffsets(#amplitudes, random)

 131     local function noise(x)

 132         local n = 0

 133         for i,a in ipairs(amplitudes) do

 134             n = n + a*cos(x/a + offsets[i])

 135         end

 136         n = n/sigma

 137         return cdf(n)

 138     end

 140     local function gradient(x)

 141         local n, dndx = 0, 0

 142         for i,a in ipairs(amplitudes) do

 143             n = n + a*cos(x/a + offsets[i])

 144             dndx = dndx - sin(x/a + offsets[i])

 145         end

 146         dndx = dndx/sigma

 147         n = n/sigma

 148         return cdf_prime(n,dndx)

 149     end

 150     return noise, gradient

 151 end

 153 noise.make2d = function(opts)

 154     local amplitudes, random = _defaultArgs(opts, "2d")

 155     local sigma = calculateSigma(amplitudes)

 156     local offsets = makeOffsets(2*#amplitudes, random)

 157     sigma = sigma/sqrt(2)

 158     local function noise(x,y)

 159         local n = 0

 160         for i,a in ipairs(amplitudes) do

 161             local angle = ((i*GR) % 1)*TAU

 162             local sina, cosa = sin(angle), cos(angle)

 163             local u = x*cosa - y*sina

 164             local v = -x*sina - y*cosa

 165             local k, w = offsets[2*i], offsets[2*i-1]

 166             n = n + a*(cos(u/a + k) + cos(v/a + w))

 167         end

 168         n = n/(2*sigma)

 169         return cdf(n)

 170     end

 171     local function gradient(x,y)

 172         local n, dx, dy = 0,0,0

 173         for i,a in ipairs(amplitudes) do

 174             local angle = ((i*GR) % 1)*TAU

 175             local sina, cosa = sin(angle), cos(angle)

 176             local u = x*cosa - y*sina

 177             local v = -x*sina - y*cosa

 178             local k, w = offsets[2*i], offsets[2*i-1]

 179             n = n + a*(cos(u/a + k) + cos(v/a + w))

 180             dx = dx + (sina*sin(v/a + w) - cosa*sin(u/a + k))

 181             dy = dy + (sina*sin(u/a + k) + cosa*sin(v/a + w))

 182         end

 183         n = n/(2*sigma)

 184         dx = dx/(2*sigma)

 185         dy = dy/(2*sigma)

 186         return cdf_prime(n,dx), cdf_prime(n,dy)

 187     end

 188     return noise, gradient

 189 end

 191 noise.make3d = function(opts)

 192     local amplitudes, random = _defaultArgs(opts, "3d")

 193     local resolution = #amplitudes

 194     local sigma = calculateSigma(amplitudes)

 195     local offsets = makeOffsets(3*#amplitudes, random)

 196     sigma = sigma/sqrt(3)

 197     local function noise(x,y,z)

 198         -- Find the biggest fibonacci number F_n such that F_n < resolution

 199         local n = floor(log((resolution-1)*SQRT5 + .5)/LOG_GR)

 200         local dec = floor(.5 + (GR^n)/SQRT5) -- F_n, using closed form Fibonacci

 201         local inc = floor(.5 + dec/GR) -- F_(n-1)

 203         local n,i,j = 0,0,0

 204         for i=0,resolution-1 do

 205             if j >= dec then

 206                 j = j - dec

 207             else

 208                 j = j + inc

 209                 if j >= resolution then

 210                     j = j - dec

 211                 end

 212             end

 213             -- Convert golden ratio sequence into polar coordinate unit vector

 214             local phi = ((i*GR) % 1)*TAU

 215             local theta = acos(-1+2*((j*GR) % 1))

 216             -- Make an orthonormal basis, where n1 is from polar phi/theta,

 217             -- n2 is roated 90 degrees along phi, and n3 is the cross product of the two

 218             local n1x,n1y,n1z = sin(phi)*cos(theta), sin(phi)*sin(theta), cos(phi)

 219             local n2x,n2y,n2z = sin(phi+TAU/4.)*cos(theta), sin(phi+TAU/4.)*sin(theta), cos(phi+TAU/4.)

 220             -- Cross product

 221             local n3x,n3y,n3z = n1y*n2z - n1z*n2y,

 222                         n1z*n2x - n1x*n2z,

 223                         n1x*n2y - n1y*n2x

 224             -- Convert pos from x/y/z coordinates to n1/n2/n3 coordinates

 225             local u = n1x*x + n1y*y + n1z*z

 226             local v = n2x*x + n2y*y + n2z*z

 227             local w = n3x*x + n3y*y + n3z*z

 228             -- Pull the amplitude from the shuffled array index ("j"), not "i",

 229             -- otherwise neighboring unit vectors will have similar amplitudes!

 230             local a = amplitudes[j+1]

 231             -- Noise is the average of cosine of distance along each axis, shifted by offsets and scaled by amplitude.

 232             n = n + a*(cos(u/a + offsets[3*i+1]) + cos(v/a + offsets[3*i+2]) + cos(w/a + offsets[3*i+3]))

 233         end

 234         return cdf(n/(3*sigma))

 235     end

 236     local function gradient(x,y,z)

 237         -- Find the biggest fibonacci number F_n such that F_n < resolution

 238         local n = floor(log((resolution-1)*SQRT5 + .5)/LOG_GR)

 239         local dec = floor(.5 + (GR^n)/SQRT5) -- F_n, using closed form Fibonacci

 240         local inc = floor(.5 + dec/GR) -- F_(n-1)

 242         local n,i,j = 0,0,0

 243         local dndx,dndy,dndz = 0,0,0

 244         for i=0,resolution-1 do

 245             if j >= dec then

 246                 j = j - dec

 247             else

 248                 j = j + inc

 249                 if j >= resolution then

 250                     j = j - dec

 251                 end

 252             end

 253             -- Convert golden ratio sequence into polar coordinate unit vector

 254             local phi = ((i*GR) % 1)*TAU

 255             local theta = acos(-1+2*((j*GR) % 1))

 256             -- Make an orthonormal basis, where n1 is from polar phi/theta,

 257             -- n2 is roated 90 degrees along phi, and n3 is the cross product of the two

 258             local n1x,n1y,n1z = sin(phi)*cos(theta), sin(phi)*sin(theta), cos(phi)

 259             local n2x,n2y,n2z = sin(phi+TAU/4.)*cos(theta), sin(phi+TAU/4.)*sin(theta), cos(phi+TAU/4.)

 260             -- Cross product

 261             local n3x,n3y,n3z = n1y*n2z - n1z*n2y,

 262                         n1z*n2x - n1x*n2z,

 263                         n1x*n2y - n1y*n2x

 264             -- Convert pos from x/y/z coordinates to n1/n2/n3 coordinates

 265             local u = n1x*x + n1y*y + n1z*z

 266             local v = n2x*x + n2y*y + n2z*z

 267             local w = n3x*x + n3y*y + n3z*z

 268             -- Pull the amplitude from the shuffled array index ("j"), not "i",

 269             -- otherwise neighboring unit vectors will have similar amplitudes!

 270             local a = amplitudes[j+1]

 271             -- Noise is the average of cosine of distance along each axis, shifted by offsets and scaled by amplitude.

 272             n = n + a*(cos(u/a + offsets[3*i+1]) + cos(v/a + offsets[3*i+2]) + cos(w/a + offsets[3*i+3]))

 274             local kx,ky,kz = -sin(u/a+offsets[3*i+1]),-sin(v/a+offsets[3*i+2]),-sin(w/a + offsets[3*i+3])

 275             dndx = dndx + (n1x*kx + n2x*ky + n3x*kz)

 276             dndy = dndy + (n1y*kx + n2y*ky + n3y*kz)

 277             dndz = dndz + (n1z*kx + n2z*ky + n3z*kz)

 278         end

 279         n = n / (3*sigma)

 280         return cdf_prime(n, dndx/(3*sigma)), cdf_prime(n, dndy/(3*sigma)), cdf_prime(n, dndz/(3*sigma))

 281     end

 282     return noise, gradient

 283 end

 285 if love and love.graphics then

 286     local shader1d = [[

 287     #define TAU 6.283185307179586476925286766559005768394338798750211641949

 288     #define MAX_RESOLUTION 64

 289     #define SIGMOID(x) (1./(1.+exp(-(x))))

 290     #define CDF(x) (.5 + .5*sign(x)*sqrt(1.-exp(-4./TAU * (x)*(x))))

 291     #define CDF_PRIME(x, dx) ((0.31831 * exp(-4./TAU * (x)*(x)) * abs(x)*(dx))/sqrt(1.0-exp(-4./TAU * (x)*(x))))

 292     extern int resolution;

 293     extern float sigma, range_min, range_max;

 294     extern float amplitudes[MAX_RESOLUTION];

 295     extern float offsets[MAX_RESOLUTION];

 296     extern bool draw_outline = false;

 298     vec4 effect(vec4 color, Image texture, vec2 texture_coords, vec2 pixel_coords)

 300         float x = mix(range_min,range_max,texture_coords.x);

 301         float noise = 0.;

 302     #ifdef GRADIENT

 303         float dndx = 0.;

 304     #endif

 305         for (int i=0; i < resolution; i++) {

 306             float a = amplitudes[i];

 307             noise += a*cos(x/a + offsets[i]);

 308     #ifdef GRADIENT

 309             dndx -= sin(x/a + offsets[i]);

 310     #endif

 312         noise /= sigma;

 313     #ifdef GRADIENT

 314         dndx = CDF_PRIME(noise, dndx/sigma);

 315         dndx = SIGMOID(dndx*2.);

 316         if (draw_outline) {

 317             return texture_coords.y > (1.-dndx) ? vec4(1.,1.,1.,1.) : vec4(0.,0.,0.,1.);

 318         } else {

 319             return vec4(dndx,dndx,dndx,1.);

 321     #else

 322         noise = CDF(noise);

 323         if (draw_outline) {

 324             return texture_coords.y > (1.-noise) ? vec4(1.,1.,1.,1.) : vec4(0.,0.,0.,1.);

 325         } else {

 326             return vec4(noise,noise,noise,1.);

 328     #endif

 332     noise.make1dShader = function(opts)

 333         local amplitudes, random = _defaultArgs(opts, "1dShader")

 334         local resolution = #amplitudes

 335         local sigma = calculateSigma(amplitudes)

 336         local offsets = makeOffsets(#amplitudes, random)

 337         local shaderCode = "#define RESOLUTION "..tostring(resolution).."\n"..shader1d

 338         local noiseShader = love.graphics.newShader(shaderCode)

 339         local gradShader = love.graphics.newShader("#define GRADIENT\n"..shaderCode)

 340         do -- Dumb hack to work around a bug in Love 0.10.2 not sending the last value

 341             table.insert(amplitudes, 1.)

 342             table.insert(offsets, 1.)

 343         end

 344         for _,shader in ipairs{noiseShader, gradShader} do

 345             shader:send("sigma",sigma)

 346             shader:send("resolution",resolution)

 347             shader:send("offsets",unpack(offsets))

 348             shader:send("amplitudes",unpack(amplitudes))

 349             shader:send("range_min", 0)

 350             shader:send("range_max", 0)

 351             if opts and opts.draw_outline then

 352                 shader:send("draw_outline", true)

 353             end

 354         end

 355         return noiseShader,gradShader

 356     end

 358     local shader2d = [[

 359     #define TAU 6.283185307179586476925286766559005768394338798750211641949

 360     #define PHI 1.618033988749894848204586834365638117720309179805762862135

 361     #define CDF(x) (.5 + .5*sign(x)*sqrt(1.-exp(-4./TAU * (x)*(x))))

 362     #define CDF_PRIME(x, dx) ((0.31831 * exp(-4./TAU * (x)*(x)) * abs(x)*(dx))/sqrt(1.0-exp(-4./TAU * (x)*(x))))

 363     #define SIGMOID(x) (1./(1.+exp(-2.*(x))))

 364     extern float sigma;

 365     extern float amplitudes[RESOLUTION];

 366     extern vec2 offsets[RESOLUTION];

 367     extern vec2 range_min, range_max;

 369     vec4 effect(vec4 color, Image texture, vec2 texture_coords, vec2 pixel_coords)

 371         vec2 pos = mix(range_min,range_max,texture_coords);

 372         float noise = 0.;

 373     #ifdef GRADIENT

 374         float dndx = 0., dndy = 0.;

 375     #endif

 376         for (int i=0; i < RESOLUTION; i++) {

 377             float angle = mod(float(i)*PHI, 1.)*TAU;

 378             float cosa = cos(angle), sina = sin(angle);

 379             float u = pos.x*cosa - pos.y*sina;

 380             float v = -pos.x*sina - pos.y*cosa;

 381             float a = amplitudes[i];

 382             float k = offsets[i].x, w = offsets[i].y;

 383             noise += a*(cos(u/a + k) + cos(v/a + w));

 384     #ifdef GRADIENT

 385             dndx += -cosa*sin(u/a + k) + sina*sin(v/a + w);

 386             dndy += sina*sin(u/a + k) + cosa*sin(v/a + w);

 387     #endif

 390         float _sigma = 2.*sigma;

 391         noise /= _sigma;

 392     #ifdef GRADIENT

 393         dndx = CDF_PRIME(noise, dndx/_sigma);

 394         dndx = SIGMOID(dndx*1.5);

 396         dndy = CDF_PRIME(noise, dndy/_sigma);

 397         dndy = SIGMOID(dndy*1.5);

 399         return vec4(dndx,dndy,0.,1.);

 400     #else

 401         noise = CDF(noise);

 402         return vec4(noise,noise,noise,1.);

 403     #endif

 407     noise.make2dShader = function(opts)

 408         local amplitudes, random = _defaultArgs(opts, "2dShader")

 409         local resolution = #amplitudes

 410         local sigma = calculateSigma(amplitudes)

 411         local offsets = makeOffsets(2*#amplitudes, random)

 412         local shaderCode = "#define RESOLUTION "..tostring(resolution).."\n"..shader2d

 413         local noiseShader = love.graphics.newShader(shaderCode)

 414         local gradShader = love.graphics.newShader("#define GRADIENT\n"..shaderCode)

 415         sigma = sigma/sqrt(2)

 416         local offsets2 = {}

 417         for i=1,#offsets-1,2 do

 418             table.insert(offsets2, {offsets[i],offsets[i+1]})

 419         end

 420         do -- Dumb hack to work around a bug in Love 0.10.2 not sending the last value

 421             table.insert(amplitudes, 1.)

 422             table.insert(offsets2, {1,1})

 423         end

 424         for _,shader in ipairs{noiseShader, gradShader} do

 425             shader:send("sigma",sigma)

 426             shader:send("offsets",unpack(offsets2))

 427             shader:send("amplitudes",unpack(amplitudes))

 428             shader:send("range_min", {0,0})

 429             shader:send("range_max", {1,1})

 430         end

 431         return noiseShader,gradShader

 432     end

 434     local shader3d = [[

 435     #define TAU 6.283185307179586476925286766559005768394338798750211641949

 436     #define PHI 1.618033988749894848204586834365638117720309179805762862135

 437     #define LOG_PHI 0.481211825059603447497758913424368423135184334385660519660

 438     #define SQRT5 2.236067977499789696409173668731276235440618359611525724270

 439     #define CDF(x) (.5 + .5*sign(x)*sqrt(1.-exp(-4./TAU * (x)*(x))))

 440     #define CDF_PRIME(x, dx) ((0.31831 * exp(-4./TAU * (x)*(x)) * abs(x)*(dx))/sqrt(1.0-exp(-4./TAU * (x)*(x))))

 441     #define SIGMOID(x) (1./(1.+exp(-2.*(x))))

 442     extern float sigma, z;

 443     extern float amplitudes[RESOLUTION];

 444     extern vec3 offsets[RESOLUTION];

 445     extern vec2 range_min, range_max;

 447     // https://www.graphics.rwth-aachen.de/media/papers/jgt.pdf

 448     vec4 effect(vec4 color, Image texture, vec2 texture_coords, vec2 pixel_coords)

 450         vec3 pos = vec3(mix(range_min,range_max,texture_coords), z);

 451         // Find the biggest fibonacci number F_n such that F_n < RESOLUTION

 452         int fib_n = int(log((float(RESOLUTION)-1.)*SQRT5 + .5)/LOG_PHI);

 453         int dec = int(.5 + pow(PHI,fib_n)/SQRT5); // F_n, using closed form Fibonacci

 454         int inc = int(.5 + dec/PHI); // F_(fib_n-1)

 456         float n = 0.;

 457     #ifdef GRADIENT

 458         float dndx = 0., dndy = 0., dndz = 0.;

 459     #endif

 460         for (int i=0, j=0; i<RESOLUTION; ++i) {

 461             if (j >= dec) {

 462                 j -= dec;

 463             } else {

 464                 j += inc;

 465                 if (j >= RESOLUTION)

 466                     j -= dec;

 468             // Convert golden ratio sequence into polar coordinate unit vector

 469             float phi = mod(float(i)*PHI,1.)*TAU;

 470             float theta = acos(mix(-1.,1.,mod(float(j)*PHI,1.)));

 471             // Make an orthonormal basis, where n1 is from polar phi/theta,

 472             // n2 is roated 90 degrees along phi, and n3 is the cross product of the two

 473             vec3 n1 = vec3(sin(phi)*cos(theta), sin(phi)*sin(theta), cos(phi));

 474             vec3 n2 = vec3(sin(phi+TAU/4.)*cos(theta), sin(phi+TAU/4.)*sin(theta), cos(phi+TAU/4.));

 475             vec3 n3 = cross(n1,n2);

 476             // Convert pos from x/y/z coordinates to n1/n2/n3 coordinates

 477             float u = dot(n1, pos);

 478             float v = dot(n2, pos);

 479             float w = dot(n3, pos);

 480             // Pull the amplitude from the shuffled array index ("j"), not "i",

 481             // otherwise neighboring unit vectors will have similar amplitudes!

 482             float a = amplitudes[j];

 483             //float a = pow(mod(float(i+1)*(PHI-1.), 1.), .3);

 484             // Noise is the average of cosine of distance along each axis, shifted by offsets and scaled by amplitude.

 485             n += a*(cos(u/a + offsets[i].x) + cos(v/a + offsets[i].y) + cos(w/a + offsets[i].z));

 486     #ifdef GRADIENT

 487             vec3 k = vec3(-sin(u/a+offsets[i].x),-sin(v/a+offsets[i].y),-sin(w/a + offsets[i].z));

 488             dndx += (n1.x*k.x + n2.x*k.y + n3.x*k.z)/3.;

 489             dndy += (n1.y*k.x + n2.y*k.y + n3.y*k.z)/3.;

 490             dndz += (n1.z*k.x + n2.z*k.y + n3.z*k.z)/3.;

 491     #endif

 493         float _sigma = 3.*sigma;

 494         n /= _sigma;

 495     #ifdef GRADIENT

 496         dndx = CDF_PRIME(n, dndx/_sigma);

 497         dndx = SIGMOID(dndx*4.);

 499         dndy = CDF_PRIME(n, dndy/_sigma);

 500         dndy = SIGMOID(dndy*4.);

 502         dndz = CDF_PRIME(n, dndz);

 503         dndz = SIGMOID(dndz*4.);

 505         return vec4(dndx,dndy,dndz, 1.);

 506     #else

 507         n = CDF(n);

 508         return vec4(n,n,n,1.);

 509     #endif

 513     noise.make3dShader = function(opts)

 514         local amplitudes, random = _defaultArgs(opts, "3dShader")

 515         local resolution = #amplitudes

 516         local sigma = calculateSigma(amplitudes)

 517         local offsets = makeOffsets(3*#amplitudes, random)

 518         local shaderCode = "#define RESOLUTION "..tostring(resolution).."\n"..shader3d

 519         local noiseShader = love.graphics.newShader(shaderCode)

 520         local gradShader = love.graphics.newShader("#define GRADIENT\n"..shaderCode)

 521         sigma = sigma/sqrt(3)

 522         local offsets2 = {}

 523         for i=1,#offsets-1,3 do

 524             table.insert(offsets2, {offsets[i],offsets[i+1],offsets[i+2]})

 525         end

 526         do -- Dumb hack to work around a bug in Love 0.10.2 not sending the last value

 527             table.insert(amplitudes, 1.)

 528             table.insert(offsets, {1,1,1})

 529         end

 530         for _,shader in ipairs{noiseShader, gradShader} do

 531             shader:send("sigma",sigma)

 532             shader:send("offsets",unpack(offsets2))

 533             shader:send("amplitudes",unpack(amplitudes))

 534             shader:send("range_min", {0,0})

 535             shader:send("range_max", {1,1})

 536         end

 537         return noiseShader, gradShader

 538     end

 539 end

 541 return noise