X-Git-Url: http://git.mutantstargoat.com/user/nuclear/?p=dosrtxon;a=blobdiff_plain;f=src%2Fnoise.c;fp=src%2Fnoise.c;h=7c1debd87d93f448398d916fa3c65debb6867a0f;hp=0000000000000000000000000000000000000000;hb=6952fc3d7981a73ffbabeca5dca741fd790949cc;hpb=584c6161bb6f3f548d5e27bb7d6895a81375881e

diff --git a/src/noise.c b/src/noise.c
new file mode 100644
index 0000000..7c1debd
--- /dev/null
+++ b/src/noise.c
@@ -0,0 +1,459 @@
+#include <math.h>
+#include <stdlib.h>
+#include "noise.h"
+
+/* ---- Ken Perlin's implementation of noise ---- */
+#define B	0x100
+#define BM	0xff
+#define N	0x1000
+#define NP	12   /* 2^N */
+#define NM	0xfff
+
+#define s_curve(t) (t * t * (3.0f - 2.0f * t))
+
+#define setup(elem, b0, b1, r0, r1) \
+	do {							\
+		float t = elem + N;		\
+		b0 = ((int)t) & BM;			\
+		b1 = (b0 + 1) & BM;			\
+		r0 = t - (int)t;			\
+		r1 = r0 - 1.0f;				\
+	} while(0)
+
+#define setup_p(elem, b0, b1, r0, r1, p) \
+	do {							\
+		float t = elem + N;		\
+		b0 = (((int)t) & BM) % p;	\
+		b1 = ((b0 + 1) & BM) % p;	\
+		r0 = t - (int)t;			\
+		r1 = r0 - 1.0f;				\
+	} while(0)
+
+
+static int perm[B + B + 2];			/* permuted index from g_n onto themselves */
+static float grad3[B + B + 2][3];	/* 3D random gradients */
+static float grad2[B + B + 2][2];	/* 2D random gradients */
+static float grad1[B + B + 2];		/* 1D random ... slopes */
+static int tables_valid;
+
+#define init_once()	if(!tables_valid) init_noise()
+
+#define SQ(x)	((x) * (x))
+#define LERP(a, b, t)	((a) + ((b) - (a)) * (t))
+
+static void init_noise()
+{
+	int i;
+	float len;
+
+	/* calculate random gradients */
+	for(i=0; i<B; i++) {
+		perm[i] = i;	/* .. and initialize permutation mapping to identity */
+
+		grad1[i] = (float)((rand() % (B + B)) - B) / B;
+
+		grad2[i][0] = (float)((rand() % (B + B)) - B) / B;
+		grad2[i][1] = (float)((rand() % (B + B)) - B) / B;
+		len = sqrt(SQ(grad2[i][0]) + SQ(grad2[i][1]));
+		if(len != 0.0f) {
+			grad2[i][0] /= len;
+			grad2[i][1] /= len;
+		}
+
+		grad3[i][0] = (float)((rand() % (B + B)) - B) / B;
+		grad3[i][1] = (float)((rand() % (B + B)) - B) / B;
+		grad3[i][2] = (float)((rand() % (B + B)) - B) / B;
+		len = sqrt(SQ(grad3[i][0]) + SQ(grad3[i][1]) + SQ(grad3[i][2]));
+		if(len != 0.0f) {
+			grad3[i][0] /= len;
+			grad3[i][1] /= len;
+			grad3[i][2] /= len;
+		}
+	}
+
+	/* permute indices by swapping them randomly */
+	for(i=0; i<B; i++) {
+		int rand_idx = rand() % B;
+
+		int tmp = perm[i];
+		perm[i] = perm[rand_idx];
+		perm[rand_idx] = tmp;
+	}
+
+	/* fill up the rest of the arrays by duplicating the existing gradients */
+	/* and permutations */
+	for(i=0; i<B+2; i++) {
+		perm[B + i] = perm[i];
+		grad1[B + i] = grad1[i];
+		grad2[B + i][0] = grad2[i][0];
+		grad2[B + i][1] = grad2[i][1];
+		grad3[B + i][0] = grad3[i][0];
+		grad3[B + i][1] = grad3[i][1];
+		grad3[B + i][2] = grad3[i][2];
+	}
+
+	tables_valid = 1;
+}
+
+
+float noise1(float x)
+{
+	int bx0, bx1;
+	float rx0, rx1, sx, u, v;
+
+	init_once();
+
+	setup(x, bx0, bx1, rx0, rx1);
+	sx = s_curve(rx0);
+	u = rx0 * grad1[perm[bx0]];
+	v = rx1 * grad1[perm[bx1]];
+	return LERP(u, v, sx);
+}
+
+float noise2(float x, float y)
+{
+	int i, j, b00, b10, b01, b11;
+	int bx0, bx1, by0, by1;
+	float rx0, rx1, ry0, ry1;
+	float sx, sy, u, v, a, b;
+
+	init_once();
+
+	setup(x, bx0, bx1, rx0, rx1);
+	setup(y, by0, by1, ry0, ry1);
+
+	i = perm[bx0];
+	j = perm[bx1];
+
+	b00 = perm[i + by0];
+	b10 = perm[j + by0];
+	b01 = perm[i + by1];
+	b11 = perm[j + by1];
+
+	/* calculate hermite inteprolating factors */
+	sx = s_curve(rx0);
+	sy = s_curve(ry0);
+
+	/* interpolate along the left edge */
+	u = grad2[b00][0] * rx0 + grad2[b00][1] * ry0;
+	v = grad2[b10][0] * rx1 + grad2[b10][1] * ry0;
+	a = LERP(u, v, sx);
+
+	/* interpolate along the right edge */
+	u = grad2[b01][0] * rx0 + grad2[b01][1] * ry1;
+	v = grad2[b11][0] * rx1 + grad2[b11][1] * ry1;
+	b = LERP(u, v, sx);
+
+	/* interpolate between them */
+	return LERP(a, b, sy);
+}
+
+float noise3(float x, float y, float z)
+{
+	int i, j;
+	int bx0, bx1, by0, by1, bz0, bz1;
+	int b00, b10, b01, b11;
+	float rx0, rx1, ry0, ry1, rz0, rz1;
+	float sx, sy, sz;
+	float u, v, a, b, c, d;
+
+	init_once();
+
+	setup(x, bx0, bx1, rx0, rx1);
+	setup(y, by0, by1, ry0, ry1);
+	setup(z, bz0, bz1, rz0, rz1);
+
+	i = perm[bx0];
+	j = perm[bx1];
+
+	b00 = perm[i + by0];
+	b10 = perm[j + by0];
+	b01 = perm[i + by1];
+	b11 = perm[j + by1];
+
+	/* calculate hermite interpolating factors */
+	sx = s_curve(rx0);
+	sy = s_curve(ry0);
+	sz = s_curve(rz0);
+
+	/* interpolate along the top slice of the cell */
+	u = grad3[b00 + bz0][0] * rx0 + grad3[b00 + bz0][1] * ry0 + grad3[b00 + bz0][2] * rz0;
+	v = grad3[b10 + bz0][0] * rx1 + grad3[b10 + bz0][1] * ry0 + grad3[b10 + bz0][2] * rz0;
+	a = LERP(u, v, sx);
+
+	u = grad3[b01 + bz0][0] * rx0 + grad3[b01 + bz0][1] * ry1 + grad3[b01 + bz0][2] * rz0;
+	v = grad3[b11 + bz0][0] * rx1 + grad3[b11 + bz0][1] * ry1 + grad3[b11 + bz0][2] * rz0;
+	b = LERP(u, v, sx);
+
+	c = LERP(a, b, sy);
+
+	/* interpolate along the bottom slice of the cell */
+	u = grad3[b00 + bz0][0] * rx0 + grad3[b00 + bz0][1] * ry0 + grad3[b00 + bz0][2] * rz1;
+	v = grad3[b10 + bz0][0] * rx1 + grad3[b10 + bz0][1] * ry0 + grad3[b10 + bz0][2] * rz1;
+	a = LERP(u, v, sx);
+
+	u = grad3[b01 + bz0][0] * rx0 + grad3[b01 + bz0][1] * ry1 + grad3[b01 + bz0][2] * rz1;
+	v = grad3[b11 + bz0][0] * rx1 + grad3[b11 + bz0][1] * ry1 + grad3[b11 + bz0][2] * rz1;
+	b = LERP(u, v, sx);
+
+	d = LERP(a, b, sy);
+
+	/* interpolate between slices */
+	return LERP(c, d, sz);
+}
+
+
+float pnoise1(float x, int period)
+{
+	int bx0, bx1;
+	float rx0, rx1, sx, u, v;
+
+	init_once();
+
+	setup_p(x, bx0, bx1, rx0, rx1, period);
+	sx = s_curve(rx0);
+	u = rx0 * grad1[perm[bx0]];
+	v = rx1 * grad1[perm[bx1]];
+	return LERP(u, v, sx);
+}
+
+float pnoise2(float x, float y, int per_x, int per_y)
+{
+	int i, j, b00, b10, b01, b11;
+	int bx0, bx1, by0, by1;
+	float rx0, rx1, ry0, ry1;
+	float sx, sy, u, v, a, b;
+
+	init_once();
+
+	setup_p(x, bx0, bx1, rx0, rx1, per_x);
+	setup_p(y, by0, by1, ry0, ry1, per_y);
+
+	i = perm[bx0];
+	j = perm[bx1];
+
+	b00 = perm[i + by0];
+	b10 = perm[j + by0];
+	b01 = perm[i + by1];
+	b11 = perm[j + by1];
+
+	/* calculate hermite inteprolating factors */
+	sx = s_curve(rx0);
+	sy = s_curve(ry0);
+
+	/* interpolate along the left edge */
+	u = grad2[b00][0] * rx0 + grad2[b00][1] * ry0;
+	v = grad2[b10][0] * rx1 + grad2[b10][1] * ry0;
+	a = LERP(u, v, sx);
+
+	/* interpolate along the right edge */
+	u = grad2[b01][0] * rx0 + grad2[b01][1] * ry1;
+	v = grad2[b11][0] * rx1 + grad2[b11][1] * ry1;
+	b = LERP(u, v, sx);
+
+	/* interpolate between them */
+	return LERP(a, b, sy);
+}
+
+float pnoise3(float x, float y, float z, int per_x, int per_y, int per_z)
+{
+	int i, j;
+	int bx0, bx1, by0, by1, bz0, bz1;
+	int b00, b10, b01, b11;
+	float rx0, rx1, ry0, ry1, rz0, rz1;
+	float sx, sy, sz;
+	float u, v, a, b, c, d;
+
+	init_once();
+
+	setup_p(x, bx0, bx1, rx0, rx1, per_x);
+	setup_p(y, by0, by1, ry0, ry1, per_y);
+	setup_p(z, bz0, bz1, rz0, rz1, per_z);
+
+	i = perm[bx0];
+	j = perm[bx1];
+
+	b00 = perm[i + by0];
+	b10 = perm[j + by0];
+	b01 = perm[i + by1];
+	b11 = perm[j + by1];
+
+	/* calculate hermite interpolating factors */
+	sx = s_curve(rx0);
+	sy = s_curve(ry0);
+	sz = s_curve(rz0);
+
+	/* interpolate along the top slice of the cell */
+	u = grad3[b00 + bz0][0] * rx0 + grad3[b00 + bz0][1] * ry0 + grad3[b00 + bz0][2] * rz0;
+	v = grad3[b10 + bz0][0] * rx1 + grad3[b10 + bz0][1] * ry0 + grad3[b10 + bz0][2] * rz0;
+	a = LERP(u, v, sx);
+
+	u = grad3[b01 + bz0][0] * rx0 + grad3[b01 + bz0][1] * ry1 + grad3[b01 + bz0][2] * rz0;
+	v = grad3[b11 + bz0][0] * rx1 + grad3[b11 + bz0][1] * ry1 + grad3[b11 + bz0][2] * rz0;
+	b = LERP(u, v, sx);
+
+	c = LERP(a, b, sy);
+
+	/* interpolate along the bottom slice of the cell */
+	u = grad3[b00 + bz0][0] * rx0 + grad3[b00 + bz0][1] * ry0 + grad3[b00 + bz0][2] * rz1;
+	v = grad3[b10 + bz0][0] * rx1 + grad3[b10 + bz0][1] * ry0 + grad3[b10 + bz0][2] * rz1;
+	a = LERP(u, v, sx);
+
+	u = grad3[b01 + bz0][0] * rx0 + grad3[b01 + bz0][1] * ry1 + grad3[b01 + bz0][2] * rz1;
+	v = grad3[b11 + bz0][0] * rx1 + grad3[b11 + bz0][1] * ry1 + grad3[b11 + bz0][2] * rz1;
+	b = LERP(u, v, sx);
+
+	d = LERP(a, b, sy);
+
+	/* interpolate between slices */
+	return LERP(c, d, sz);
+}
+
+
+float fbm1(float x, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(i=0; i<octaves; i++) {
+		res += noise1(x * freq) / freq;
+		freq *= 2.0f;
+	}
+	return res;
+}
+
+float fbm2(float x, float y, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(i=0; i<octaves; i++) {
+		res += noise2(x * freq, y * freq) / freq;
+		freq *= 2.0f;
+	}
+	return res;
+}
+
+float fbm3(float x, float y, float z, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(i=0; i<octaves; i++) {
+		res += noise3(x * freq, y * freq, z * freq) / freq;
+		freq *= 2.0f;
+	}
+	return res;
+
+}
+
+
+float pfbm1(float x, int per, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(i=0; i<octaves; i++) {
+		res += pnoise1(x * freq, per) / freq;
+		freq *= 2.0f;
+		per *= 2;
+	}
+	return res;
+}
+
+float pfbm2(float x, float y, int per_x, int per_y, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(i=0; i<octaves; i++) {
+		res += pnoise2(x * freq, y * freq, per_x, per_y) / freq;
+		freq *= 2.0f;
+		per_x *= 2;
+		per_y *= 2;
+	}
+	return res;
+}
+
+float pfbm3(float x, float y, float z, int per_x, int per_y, int per_z, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(int i=0; i<octaves; i++) {
+		res += pnoise3(x * freq, y * freq, z * freq, per_x, per_y, per_z) / freq;
+		freq *= 2.0f;
+		per_x *= 2;
+		per_y *= 2;
+		per_z *= 2;
+	}
+	return res;
+}
+
+
+float turbulence1(float x, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(i=0; i<octaves; i++) {
+		res += fabs(noise1(x * freq) / freq);
+		freq *= 2.0f;
+	}
+	return res;
+}
+
+float turbulence2(float x, float y, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(i=0; i<octaves; i++) {
+		res += fabs(noise2(x * freq, y * freq) / freq);
+		freq *= 2.0f;
+	}
+	return res;
+}
+
+float turbulence3(float x, float y, float z, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(i=0; i<octaves; i++) {
+		res += fabs(noise3(x * freq, y * freq, z * freq) / freq);
+		freq *= 2.0f;
+	}
+	return res;
+}
+
+
+float pturbulence1(float x, int per, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(i=0; i<octaves; i++) {
+		res += fabs(pnoise1(x * freq, per) / freq);
+		freq *= 2.0f;
+		per *= 2;
+	}
+	return res;
+}
+
+float pturbulence2(float x, float y, int per_x, int per_y, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(i=0; i<octaves; i++) {
+		res += fabs(pnoise2(x * freq, y * freq, per_x, per_y) / freq);
+		freq *= 2.0f;
+		per_x *= 2;
+		per_y *= 2;
+	}
+	return res;
+}
+
+float pturbulence3(float x, float y, float z, int per_x, int per_y, int per_z, int octaves)
+{
+	int i;
+	float res = 0.0f, freq = 1.0f;
+	for(i=0; i<octaves; i++) {
+		res += fabs(pnoise3(x * freq, y * freq, z * freq, per_x, per_y, per_z) / freq);
+		freq *= 2.0f;
+		per_x *= 2;
+		per_y *= 2;
+		per_z *= 2;
+	}
+	return res;
+}