# Perlin noise plot with Matplotlib

from pylab import *
import random

XDIM, YDIM = 11, 11     # gradient mesh size
SCALE = 50              # noise grid size per gradient mesh unit size
XDP, YDP = (XDIM - 2) * SCALE + 1, (YDIM - 2) * SCALE + 1

phi = zeros((YDIM, XDIM))   # gradient phase angles
perlin = zeros((YDP, XDP))  # Perlin noise array

# make the plot reproducible by using a fixed random seed
random.seed("Perlin")

# choose random phase angles for gradients
for y in range(YDIM):
    for x in range(XDIM):
        phi[y,x] = 2 * pi * random.random()

def grad(x0, y0, dx, dy):
    "Compute the dot product between the gradient and position vector"
    gx, gy = cos(phi[y0,x0]), sin(phi[y0,x0])
    return gx * dx + gy * dy

def inter(a, b, w):
    "Interpolate between a and b using weight w"
    return (b - a) * ((w * (w * 6 - 15) + 10) * w * w * w) + a  # smootherstep
    #return (b - a) * (3 - 2 * w) * w * w + a   # smoothstep
    #return (b - a) * w + a                     # linear

# compute Perlin noise
for yi in range(YDP):
    for xi in range(XDP):
        x, y = xi / SCALE, yi / SCALE
        x0, y0 = int(x), int(y)
        dx, dy = x - x0, y - y0

        g1 = grad(x0, y0, dx, dy)
        g2 = grad(x0 + 1, y0, -1 + dx, dy)
        g3 = grad(x0, y0 + 1, dx, -1 + dy)
        g4 = grad(x0 + 1, y0 + 1, -1 + dx, -1 + dy)

        perlin[yi,xi] = inter(inter(g1, g2, dx), inter(g3, g4, dx), dy)

print(f"Max {perlin.max()}, min {perlin.min()}")

# use Matplotlib to create a plot
figure(figsize = (8, 9))
X = array([x / SCALE for x in range(XDP)])
Y = array([y / SCALE for y in range(YDP)])
imshow(perlin, interpolation = "bicubic", cmap = plt.cm.bwr, vmin = -.8, vmax = .8,
  extent = (X.min(), X.max(), Y.max(), Y.min()))
cb = colorbar(orientation = "horizontal")
contour(X, Y, perlin, (0,), linewidths = 2, colors = "green")
xticks(range(XDIM-1), labels = "" * XDIM)
yticks(range(YDIM-1), labels = "" * YDIM)
grid(lw = 1.2, color = "black", alpha = .8, ls = "dashed")
gca().set_position([.1, .2, .8, .8])
cb.ax.set_position([.1, -.68, .8, .8])
title("2-D Perlin noise with contour line at zero")
savefig("perlin_noise_with_contour.svg")
show()