// Perlin Noise http://mrl.nyu.edu/~perlin/doc/oscar.html
module coneneko.noise;
import std.random, std.math, coneneko.math;

class FloatArray3D
{
	const int width;
	private float[] array;
	
	this(int width)
	{
		this.width = width;
		array.length = width * width * width;
	}

	float get(int x, int y, int z)
	in
	{
		assert(0 <= x && x < width);
		if (!(0 <= y && y < width)) printf("%d %d", y, width);
		assert(0 <= y && y < width);
		assert(0 <= z && z < width);
	}
	body
	{
		return array[indexOf(x, y, z)];
	}

	void set(int x, int y, int z, float value)
	in
	{
		assert(0 <= x && x < width);
		assert(0 <= y && y < width);
		assert(0 <= z && z < width);
	}
	body
	{
		array[indexOf(x, y, z)] = value;
	}

	private int indexOf(int x, int y, int z)
	in
	{
		assert(0 <= x && x < width);
		assert(0 <= y && y < width);
		assert(0 <= z && z < width);
	}
	out (result)
	{
		assert(0 <= result && result < array.length);
	}
	body
	{
		return width * width * x + width * y + z;
	}
}

class NoiseBox
{
	public this(int width)
	in
	{
		assert(2 <= width && width <= 256);
	}
	body
	{
		this.width = width;
		randomTable = new FloatArray3D(width);
		for (int z = 0; z < width; z++)
		{
			for (int y = 0; y < width; y++)
			{
				for (int x = 0; x < width; x++)
				{
					float a = cast(float)rand() / cast(float)uint.max;
					assert(0 <= a && a <= 1);
					randomTable.set(x, y, z, a);
				}
			}
		}
	}

	private FloatArray3D randomTable;
	private const int width;

	public float at(float x, float y, float z)
	in
	{
		assert(0 <= x && x < 1);
		assert(0 <= y && y < 1);
		assert(0 <= z && z < 1);
	}
	body
	{
		return catmullRom3D(
			indexOf(x), indexOf(y), indexOf(z),
			toRelative(x), toRelative(y), toRelative(z)
			);
	}

	float catmullRom3D(int indexOfX, int indexOfY, int indexOfZ, float x, float y, float z)
	in
	{
		assert(0 <= x && x < 1);
		assert(0 <= y && y < 1);
		assert(0 <= z && z < 1);
	}
	body
	{
		return catmullRom(
			catmullRom2D(indexOfX, indexOfY, indexOfZ - 1, x, y),
			catmullRom2D(indexOfX, indexOfY, indexOfZ    , x, y),
			catmullRom2D(indexOfX, indexOfY, indexOfZ + 1, x, y),
			catmullRom2D(indexOfX, indexOfY, indexOfZ + 2, x, y),
			z);
	}

	private float catmullRom2D(int indexOfX, int indexOfY, int indexOfZ, float x, float y)
	{
		return catmullRom(
			catmullRom1D(indexOfX, indexOfY - 1, indexOfZ, x),
			catmullRom1D(indexOfX, indexOfY    , indexOfZ, x),
			catmullRom1D(indexOfX, indexOfY + 1, indexOfZ, x),
			catmullRom1D(indexOfX, indexOfY + 2, indexOfZ, x),
			y);
	}

	private float catmullRom1D(int indexOfX, int indexOfY, int indexOfZ, float x)
	{
		return catmullRom(
			randomTableAt(indexOfX - 1, indexOfY, indexOfZ),
			randomTableAt(indexOfX    , indexOfY, indexOfZ),
			randomTableAt(indexOfX + 1, indexOfY, indexOfZ),
			randomTableAt(indexOfX + 2, indexOfY, indexOfZ),
			x);
	}

	private float randomTableAt(int indexOfX, int indexOfY, int indexOfZ)
	{
		if (indexOfX < 0) indexOfX += width;
		if (indexOfY < 0) indexOfY += width;
		if (indexOfZ < 0) indexOfZ += width;
		if (indexOfX >= width) indexOfX -= width;
		if (indexOfY >= width) indexOfY -= width;
		if (indexOfZ >= width) indexOfZ -= width;
		return randomTable.get(indexOfX, indexOfY, indexOfZ);
	}

	private int indexOf(float f)
	in
	{
		assert(0 <= f && f < 1);
	}
	out (result)
	{
		assert(0 <= result && result < width);
	}
	body
	{
		return cast(int)floor(width * f);
	}

	private float toRelative(float f)
	out (result)
	{
		assert(0 <= result && result < 1);
	}
	body
	{
		return width * f % 1; // 小数部
	}
}

float catmullRom(float x0, float x1, float x2, float x3, float t) // 一次元用
{
	const float[4][4] mm =
	[
		[ -1,  3, -3,  1],
		[  2, -5,  4, -1],
		[ -1,  0,  1,  0],
		[  0,  2,  0,  0],
	];
	Matrix m;
	m.m[][] = mm[][];
	float t2 = t * t;
	auto vt = vector(t2 * t, t2, t, 1) * m;
	return (x0 * vt.x + x1 * vt.y + x2 * vt.z + x3 * vt.w) / 2;
}

class PerlinNoise
{
	private this()
	{
		noiseBox0 = new NoiseBox(4);
		noiseBox1 = new NoiseBox(8);
		noiseBox2 = new NoiseBox(16);
		noiseBox3 = new NoiseBox(32);
	}
	
	public float at(float x, float y, float z)
	in
	{
		assert(0 <= x && x < 1);
		assert(0 <= y && y < 1);
		assert(0 <= z && z < 1);
	}
	out (result)
	{
		assert(-1 <= result && result <= 1);
	}
	body
	{
		float result = noiseBox0.at(x, y, z);
		result += noiseBox1.at(x, y, z) * pow(cast(real)0.25, 1);
		result += noiseBox2.at(x, y, z) * pow(cast(real)0.25, 2);
		result += noiseBox3.at(x, y, z) * pow(cast(real)0.25, 3);
		result /= 1 + pow(cast(real)0.25, 1) + pow(cast(real)0.25, 2) + pow(cast(real)0.25, 3);
		return clamp(result * 2.0 - 1.0, -1.0, 1.0);
	}

	private NoiseBox noiseBox0, noiseBox1, noiseBox2, noiseBox3;

	public static PerlinNoise getInstance()
	{
		if (!_instance) _instance = new PerlinNoise();
		return _instance;
	}

	private static PerlinNoise _instance;
}

float noise(float x, float y, float z)
in
{
	assert(0 <= x && x < 1);
	assert(0 <= y && y < 1);
	assert(0 <= z && z < 1);
}
out (result)
{
	assert(-1 <= result && result <= 1);
}
body
{
	return PerlinNoise.getInstance().at(x, y, z);
}

ubyte noise2(ubyte x, ubyte y, ubyte z)
{
	float fx = x / 256.0, fy = y / 256.0, fz = z / 256.0;
	auto v = noise(fx, fy, fz) * 0.5 + 0.5;
	if (v < 0.0 || 1.0 <= v) throw new Exception("noise2");
	return cast(ubyte)floor(v * 256);
}
