From 3b52c25b520e0e1b9143cafaee2670e57c89ef0d Mon Sep 17 00:00:00 2001
From: John Tsiombikas <nuclear@mutantstargoat.com>
Date: Sun, 25 Mar 2018 22:02:45 +0300
Subject: [PATCH] separated the smoke-text effect to reuse it in multiple
 parts if necessary

---
 src/greets.c   |  241 +--------------------------
 src/noise.c    |  505 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 src/noise.h    |   35 ++++
 src/smoketxt.c |  371 +++++++++++++++++++++++++++++++++++++++++
 src/smoketxt.h |   21 +++
 5 files changed, 939 insertions(+), 234 deletions(-)
 create mode 100644 src/noise.c
 create mode 100644 src/noise.h
 create mode 100644 src/smoketxt.c
 create mode 100644 src/smoketxt.h

diff --git a/src/greets.c b/src/greets.c
index 353d29c..8b128c6 100644
--- a/src/greets.c
+++ b/src/greets.c
@@ -11,6 +11,7 @@
 #include "util.h"
 #include "gfxutil.h"
 #include "timer.h"
+#include "smoketxt.h"
 
 #ifdef MSDOS
 #include "dos/gfx.h"	/* for wait_vsync assembly macro */
@@ -18,66 +19,14 @@
 void wait_vsync(void);
 #endif
 
-/* if defined, use bilinear interpolation for dispersion field vectors */
-#define BILERP_FIELD
-/* if defined randomize field vectors by RAND_FIELD_MAX */
-#define RANDOMIZE_FIELD
-
-#define RAND_FIELD_MAX	0.7
-
 #define BLUR_RAD	5
 
-#define PCOUNT		4000
-#define MAX_LIFE	7.0f
-#define PALPHA		1.0f
-#define ZBIAS		0.25
-#define DRAG		0.95
-#define FORCE		0.07
-#define FREQ		0.085
-static float wind[] = {-0.0, 0.0, 0.01};
-
-
-struct vec2 {
-	float x, y;
-};
-
-struct vec3 {
-	float x, y, z;
-};
-
-struct particle {
-	float x, y, z;
-	float vx, vy, vz;	/* velocity */
-	int r, g, b;
-	float life;
-};
-
-struct emitter {
-	struct particle *plist;
-	int pcount;
-};
-
-struct vfield {
-	struct vec2 pos, size;
-
-	int width, height;
-	int xshift;
-	struct vec2 *v;
-};
-
 
 static int init(void);
 static void destroy(void);
 static void start(long trans_time);
 static void draw(void);
 
-int init_emitter(struct emitter *em, int num, unsigned char *map, int xsz, int ysz);
-void update_particles(struct emitter *em, float dt);
-void draw_particles(struct emitter *em);
-
-int init_vfield_load(struct vfield *vf, const char *fname);
-void vfield_eval(struct vfield *vf, float x, float y, struct vec2 *dir);
-
 
 static struct screen scr = {
 	"greets",
@@ -87,11 +36,10 @@ static struct screen scr = {
 	draw
 };
 
-static struct emitter em;
-static struct vfield vfield;
-static struct g3d_vertex *varr;
 static long start_time;
 
+static struct smktxt *stx;
+
 static uint16_t *cur_smokebuf, *prev_smokebuf;
 static int smokebuf_size;
 #define smokebuf_start	(cur_smokebuf < prev_smokebuf ? cur_smokebuf : prev_smokebuf)
@@ -113,23 +61,7 @@ struct screen *greets_screen(void)
 
 static int init(void)
 {
-	int xsz, ysz;
-	unsigned char *pixels;
-
-	if(!(pixels = img_load_pixels("data/greets1.png", &xsz, &ysz, IMG_FMT_GREY8))) {
-		fprintf(stderr, "failed to load particle spawn map\n");
-		return -1;
-	}
-
-	init_emitter(&em, PCOUNT, pixels, xsz, ysz);
-	img_free_pixels(pixels);
-
-	if(!(varr = malloc(PCOUNT * sizeof *varr))) {
-		perror("failed to allocate particle vertex buffer\n");
-		return -1;
-	}
-
-	if(init_vfield_load(&vfield, "data/vfield1") == -1) {
+	if(!(stx = create_smktxt("data/greets1.png", "data/vfield1"))) {
 		return -1;
 	}
 
@@ -145,9 +77,8 @@ static int init(void)
 
 static void destroy(void)
 {
-	free(varr);
-	free(vfield.v);
 	free(smokebuf_start);
+	destroy_smktxt(stx);
 }
 
 static void start(long trans_time)
@@ -196,7 +127,7 @@ static void update(void)
 	prev_my = mouse_y;
 	prev_bmask = mouse_bmask;
 
-	update_particles(&em, dt);
+	update_smktxt(stx, dt);
 }
 
 static void draw(void)
@@ -220,7 +151,7 @@ static void draw(void)
 	memcpy(cur_smokebuf, prev_smokebuf, smokebuf_size);
 
 	g3d_framebuffer(fb_width, fb_height, cur_smokebuf);
-	draw_particles(&em);
+	draw_smktxt(stx);
 	g3d_framebuffer(fb_width, fb_height, fb_pixels);
 
 	dest = fb_pixels;
@@ -251,161 +182,3 @@ static void draw(void)
 	swap_buffers(fb_pixels);
 	last_swap = get_msec();
 }
-
-
-
-int init_emitter(struct emitter *em, int num, unsigned char *map, int xsz, int ysz)
-{
-	int i, x, y;
-	float aspect = (float)xsz / (float)ysz;
-	struct particle *p;
-
-	if(!(em->plist = malloc(num * sizeof *em->plist))) {
-		return -1;
-	}
-	em->pcount = num;
-
-	p = em->plist;
-	for(i=0; i<num; i++) {
-		do {
-			x = rand() % xsz;
-			y = rand() % ysz;
-		} while(map[y * xsz + x] < 128);
-
-		p->x = (float)x / (float)xsz - 0.5;
-		p->y = -(float)y / (float)xsz + 0.5 / aspect;
-		p->z = ((float)i / (float)num * 2.0 - 1.0) * 0.005;
-		p->r = p->g = p->b = 255;
-		p->vx = p->vy = p->vz = 0.0f;
-		p->life = MAX_LIFE;
-		++p;
-	}
-	return 0;
-}
-
-void update_particles(struct emitter *em, float dt)
-{
-	int i;
-	struct vec2 accel;
-	struct particle *p = em->plist;
-	struct g3d_vertex *v = varr;
-
-	for(i=0; i<em->pcount; i++) {
-		vfield_eval(&vfield, p->x, p->y, &accel);
-		p->x += p->vx * DRAG * dt;
-		p->y += p->vy * DRAG * dt;
-		p->z += p->vz * DRAG * dt;
-		p->vx += (wind[0] + accel.x * FORCE) * dt;
-		p->vy += (wind[1] + accel.y * FORCE) * dt;
-		p->vz += (wind[2] + p->z * ZBIAS) * dt;
-		p->life -= dt;
-		if(p->life < 0.0f) p->life = 0.0f;
-
-		v->x = p->x;
-		v->y = p->y;
-		v->z = p->z;
-		v->w = 1.0f;
-		v->a = cround64(p->life * 255.0 / MAX_LIFE);
-		v->r = 0;
-		v->g = (v->a & 0xe0) >> 3;
-		v->b = (v->a & 0x1f) << 3;
-		++v;
-
-		++p;
-	}
-}
-
-void draw_particles(struct emitter *em)
-{
-	g3d_draw(G3D_POINTS, varr, PCOUNT);
-}
-
-
-int init_vfield_load(struct vfield *vf, const char *fname)
-{
-	FILE *fp;
-	int tmp;
-
-	if(!(fp = fopen(fname, "rb"))) {
-		fprintf(stderr, "failed to open vector field: %s\n", fname);
-		return -1;
-	}
-	if(fread(&vf->width, sizeof vf->width, 1, fp) < 1 ||
-			fread(&vf->height, sizeof vf->height, 1, fp) < 1) {
-		fprintf(stderr, "init_vfield_load: unexpected end of file while reading header\n");
-		fclose(fp);
-		return -1;
-	}
-
-	/* assume xsz is pow2 otherwise fuck you */
-	tmp = vf->width - 1;
-	vf->xshift = 0;
-	while(tmp) {
-		++vf->xshift;
-		tmp >>= 1;
-	}
-
-	if(!(vf->v = malloc(vf->width * vf->height * sizeof *vf->v))) {
-		fprintf(stderr, "failed to allocate %dx%d vector field\n", vf->width, vf->height);
-		fclose(fp);
-		return -1;
-	}
-	if(fread(vf->v, sizeof *vf->v, vf->width * vf->height, fp) < vf->width * vf->height) {
-		fprintf(stderr, "init_vfield_load: unexpected end of file while reading %dx%d vector field\n",
-				vf->width, vf->height);
-		fclose(fp);
-		return -1;
-	}
-	fclose(fp);
-
-	vf->pos.x = vf->pos.y = 0;
-	vf->size.x = vf->size.y = 1;
-	return 0;
-}
-
-void vfield_eval(struct vfield *vf, float x, float y, struct vec2 *dir)
-{
-	int px, py;
-	float tx, ty;
-	struct vec2 *p1, *p2;
-	struct vec2 left, right;
-
-	x = ((x - vf->pos.x) / vf->size.x + 0.5f) * vf->width;
-	y = ((y - vf->pos.y) / vf->size.y + 0.5f) * vf->height;
-	x = floor(x);
-	y = floor(y);
-
-	if(x < 0) x = 0;
-	if(y < 0) y = 0;
-	if(x > vf->width - 2) x = vf->width - 2;
-	if(y > vf->height - 2) y = vf->height - 2;
-
-	px = (int)x;
-	py = (int)y;
-
-	p1 = vf->v + (py << vf->xshift) + px;
-#ifdef BILERP_FIELD
-	p2 = p1 + vf->width;
-
-	tx = fmod(x, 1.0f);
-	ty = fmod(y, 1.0f);
-
-	left.x = p1->x + (p2->x - p1->x) * ty;
-	left.y = p1->y + (p2->y - p1->y) * ty;
-	++p1;
-	++p2;
-	right.x = p1->x + (p2->x - p1->x) * ty;
-	right.y = p1->y + (p2->y - p1->y) * ty;
-
-	dir->x = left.x + (right.x - left.x) * tx;
-	dir->y = left.y + (right.y - left.y) * ty;
-#else
-	dir->x = p1->x;
-	dir->y = p1->y;
-#endif
-
-#ifdef RANDOMIZE_FIELD
-	dir->x += ((float)rand() / RAND_MAX - 0.5) * RAND_FIELD_MAX;
-	dir->y += ((float)rand() / RAND_MAX - 0.5) * RAND_FIELD_MAX;
-#endif
-}
diff --git a/src/noise.c b/src/noise.c
new file mode 100644
index 0000000..b3b3f9c
--- /dev/null
+++ b/src/noise.c
@@ -0,0 +1,505 @@
+#include <stdlib.h>
+#include <math.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()
+
+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;
+		if((len = sqrt(grad2[i][0] * grad2[i][0] + grad2[i][1] * grad2[i][1])) != 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;
+		if((len = sqrt(grad3[i][0] * grad3[i][0] + grad3[i][1] * grad3[i][1]) + grad3[i][2] * grad3[i][2]) != 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;
+}
+
+#define lerp(a, b, t) ((a) + ((b) - (a)) * (t))
+#define dotgrad2(g, x, y) ((g)[0] * (x) + (g)[1] * (y))
+#define dotgrad3(g, x, y, z) ((g)[0] * (x) + (g)[1] * (y) + (g)[2] * (z))
+
+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 = dotgrad2(grad2[b00], rx0, ry0);
+	v = dotgrad2(grad2[b10], rx1, ry0);
+	a = lerp(u, v, sx);
+
+	/* interpolate along the right edge */
+	u = dotgrad2(grad2[b01], rx0, ry1);
+	v = dotgrad2(grad2[b11], rx1, 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 = dotgrad3(grad3[b00 + bz0], rx0, ry0, rz0);
+	v = dotgrad3(grad3[b10 + bz0], rx1, ry0, rz0);
+	a = lerp(u, v, sx);
+
+	u = dotgrad3(grad3[b01 + bz0], rx0, ry1, rz0);
+	v = dotgrad3(grad3[b11 + bz0], rx1, ry1, rz0);
+	b = lerp(u, v, sx);
+
+	c = lerp(a, b, sy);
+
+	/* interpolate along the bottom slice of the cell */
+	u = dotgrad3(grad3[b00 + bz0], rx0, ry0, rz1);
+	v = dotgrad3(grad3[b10 + bz0], rx1, ry0, rz1);
+	a = lerp(u, v, sx);
+
+	u = dotgrad3(grad3[b01 + bz0], rx0, ry1, rz1);
+	v = dotgrad3(grad3[b11 + bz0], rx1, ry1, rz1);
+	b = lerp(u, v, sx);
+
+	d = lerp(a, b, sy);
+
+	/* interpolate between slices */
+	return lerp(c, d, sz);
+}
+
+float noise4(float x, float y, float z, float w)
+{
+	return 0;	/* TODO */
+}
+
+
+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 = dotgrad2(grad2[b00], rx0, ry0);
+	v = dotgrad2(grad2[b10], rx1, ry0);
+	a = lerp(u, v, sx);
+
+	/* interpolate along the right edge */
+	u = dotgrad2(grad2[b01], rx0, ry1);
+	v = dotgrad2(grad2[b11], rx1, 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 = dotgrad3(grad3[b00 + bz0], rx0, ry0, rz0);
+	v = dotgrad3(grad3[b10 + bz0], rx1, ry0, rz0);
+	a = lerp(u, v, sx);
+
+	u = dotgrad3(grad3[b01 + bz0], rx0, ry1, rz0);
+	v = dotgrad3(grad3[b11 + bz0], rx1, ry1, rz0);
+	b = lerp(u, v, sx);
+
+	c = lerp(a, b, sy);
+
+	/* interpolate along the bottom slice of the cell */
+	u = dotgrad3(grad3[b00 + bz0], rx0, ry0, rz1);
+	v = dotgrad3(grad3[b10 + bz0], rx1, ry0, rz1);
+	a = lerp(u, v, sx);
+
+	u = dotgrad3(grad3[b01 + bz0], rx0, ry1, rz1);
+	v = dotgrad3(grad3[b11 + bz0], rx1, ry1, rz1);
+	b = lerp(u, v, sx);
+
+	d = lerp(a, b, sy);
+
+	/* interpolate between slices */
+	return lerp(c, d, sz);
+}
+
+float pnoise4(float x, float y, float z, float w, int per_x, int per_y, int per_z, int per_w)
+{
+	return 0;
+}
+
+
+float fbm1(float x, int octaves)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int i=0; i<octaves; i++) {
+		res += noise1(x * freq) / freq;
+		freq *= 2.0f;
+	}
+	return res;
+}
+
+float fbm2(float x, float y, int octaves)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int 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)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int i=0; i<octaves; i++) {
+		res += noise3(x * freq, y * freq, z * freq) / freq;
+		freq *= 2.0f;
+	}
+	return res;
+
+}
+
+float fbm4(float x, float y, float z, float w, int octaves)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int i=0; i<octaves; i++) {
+		res += noise4(x * freq, y * freq, z * freq, w * freq) / freq;
+		freq *= 2.0f;
+	}
+	return res;
+}
+
+
+float pfbm1(float x, int per, int octaves)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int 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)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int 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)
+{
+	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 pfbm4(float x, float y, float z, float w, int per_x, int per_y, int per_z, int per_w, int octaves)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int i=0; i<octaves; i++) {
+		res += pnoise4(x * freq, y * freq, z * freq, w * freq,
+				per_x, per_y, per_z, per_w) / freq;
+		freq *= 2.0f;
+		per_x *= 2;
+		per_y *= 2;
+		per_z *= 2;
+		per_w *= 2;
+	}
+	return res;
+}
+
+
+float turbulence1(float x, int octaves)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int i=0; i<octaves; i++) {
+		res += fabs(noise1(x * freq) / freq);
+		freq *= 2.0f;
+	}
+	return res;
+}
+
+float turbulence2(float x, float y, int octaves)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int 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)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int i=0; i<octaves; i++) {
+		res += fabs(noise3(x * freq, y * freq, z * freq) / freq);
+		freq *= 2.0f;
+	}
+	return res;
+}
+
+float turbulence4(float x, float y, float z, float w, int octaves)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int i=0; i<octaves; i++) {
+		res += fabs(noise4(x * freq, y * freq, z * freq, w * freq) / freq);
+		freq *= 2.0f;
+	}
+	return res;
+}
+
+
+float pturbulence1(float x, int per, int octaves)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int 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)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int 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)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int 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;
+}
+
+float pturbulence4(float x, float y, float z, float w, int per_x, int per_y, int per_z, int per_w, int octaves)
+{
+	float res = 0.0f, freq = 1.0f;
+	for(int i=0; i<octaves; i++) {
+		res += fabs(pnoise4(x * freq, y * freq, z * freq, w * freq,
+				per_x, per_y, per_z, per_w) / freq);
+		freq *= 2.0f;
+		per_x *= 2;
+		per_y *= 2;
+		per_z *= 2;
+		per_w *= 2;
+	}
+	return res;
+}
diff --git a/src/noise.h b/src/noise.h
new file mode 100644
index 0000000..7c9eaaa
--- /dev/null
+++ b/src/noise.h
@@ -0,0 +1,35 @@
+#ifndef NOISE_H_
+#define NOISE_H_
+
+float noise1(float x);
+float noise2(float x, float y);
+float noise3(float x, float y, float z);
+/*float noise4(float x, float y, float z, float w);*/
+
+float pnoise1(float x, int period);
+float pnoise2(float x, float y, int per_x, int per_y);
+float pnoise3(float x, float y, float z, int per_x, int per_y, int per_z);
+/*float pnoise4(float x, float y, float z, float w, int per_x, int per_y, int per_z, int per_w);*/
+
+float fbm1(float x, int octaves);
+float fbm2(float x, float y, int octaves);
+float fbm3(float x, float y, float z, int octaves);
+/*float fbm4(float x, float y, float z, float w, int octaves);*/
+
+float pfbm1(float x, int per, int octaves);
+float pfbm2(float x, float y, int per_x, int per_y, int octaves);
+float pfbm3(float x, float y, float z, int per_x, int per_y, int per_z, int octaves);
+/*float pfbm4(float x, float y, float z, float w, int per_x, int per_y, int per_z, int per_w, int octaves);*/
+
+float turbulence1(float x, int octaves);
+float turbulence2(float x, float y, int octaves);
+float turbulence3(float x, float y, float z, int octaves);
+/*float turbulence4(float x, float y, float z, float w, int octaves);*/
+
+float pturbulence1(float x, int per, int octaves);
+float pturbulence2(float x, float y, int per_x, int per_y, int octaves);
+float pturbulence3(float x, float y, float z, int per_x, int per_y, int per_z, int octaves);
+/*float pturbulence4(float x, float y, float z, float w, int per_x, int per_y, int per_z, int per_w, int octaves);*/
+
+
+#endif	/* NOISE_H_ */
diff --git a/src/smoketxt.c b/src/smoketxt.c
new file mode 100644
index 0000000..841ca28
--- /dev/null
+++ b/src/smoketxt.c
@@ -0,0 +1,371 @@
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <math.h>
+#include <errno.h>
+#include "imago2.h"
+#include "3dgfx.h"
+#include "smoketxt.h"
+#include "util.h"
+#include "noise.h"
+
+/* if defined, use bilinear interpolation for dispersion field vectors */
+#define BILERP_FIELD
+/* if defined randomize field vectors by RAND_FIELD_MAX */
+#define RANDOMIZE_FIELD
+
+#define RAND_FIELD_MAX	0.7
+
+#define DFL_PCOUNT		4000
+#define DFL_MAX_LIFE	7.0f
+#define DFL_PALPHA		1.0f
+#define DFL_ZBIAS		0.25
+#define DFL_DRAG		0.95
+#define DFL_FORCE		0.07
+#define DFL_FREQ		0.085
+
+struct vec2 {
+	float x, y;
+};
+
+struct vec3 {
+	float x, y, z;
+};
+
+struct particle {
+	float x, y, z;
+	float vx, vy, vz;	/* velocity */
+	int r, g, b;
+	float life;
+};
+
+struct emitter {
+	struct particle *plist;
+	int pcount;
+
+	struct vec3 wind;
+	float drag;
+	float max_life;
+
+	struct g3d_vertex *varr;
+};
+
+struct vfield {
+	struct vec2 pos, size;
+
+	int width, height;
+	int xshift;
+	struct vec2 *v;
+};
+
+
+struct smktxt {
+	struct emitter em;
+	struct vfield vfield;
+
+	unsigned char *img_pixels;
+	int img_xsz, img_ysz;
+};
+
+static int init_emitter(struct emitter *em, int num, unsigned char *map, int xsz, int ysz);
+static int load_vfield(struct vfield *vf, const char *fname);
+static void vfield_eval(struct vfield *vf, float x, float y, struct vec2 *dir);
+
+struct smktxt *create_smktxt(const char *imgname, const char *vfieldname)
+{
+	struct smktxt *stx;
+
+	if(!imgname) return 0;
+
+	if(!(stx = calloc(sizeof *stx, 1))) {
+		return 0;
+	}
+	if(!(stx->img_pixels = img_load_pixels(imgname, &stx->img_xsz, &stx->img_ysz, IMG_FMT_GREY8))) {
+		fprintf(stderr, "create_smktxt: failed to load particle spawnmap: %s\n", imgname);
+		free(stx);
+		return 0;
+	}
+
+	stx->em.pcount = DFL_PCOUNT;
+	stx->em.wind.z = 0.01;
+	stx->em.drag = DFL_DRAG;
+	stx->em.max_life = DFL_MAX_LIFE;
+
+	if(vfieldname) {
+		if(load_vfield(&stx->vfield, vfieldname) == -1) {
+			img_free_pixels(stx->img_pixels);
+			free(stx);
+			return 0;
+		}
+	}
+
+	return stx;
+}
+
+void destroy_smktxt(struct smktxt *stx)
+{
+	if(!stx) return;
+
+	img_free_pixels(stx->img_pixels);
+	free(stx);
+}
+
+int gen_smktxt_vfield(struct smktxt *stx, int xres, int yres, float xfreq, float yfreq)
+{
+	int i, j;
+	unsigned int tmp;
+	struct vec2 *vptr;
+	struct vfield *vf = &stx->vfield;
+
+	if(!(vf->v = malloc(xres * yres * sizeof *vf->v))) {
+		fprintf(stderr, "failed to allocate %dx%d vector field\n", xres, yres);
+		return -1;
+	}
+
+	vf->width = xres;
+	vf->height = yres;
+	vf->pos.x = vf->pos.y = 0.0f;
+	vf->size.x = vf->size.y = 1.0f;
+
+	/* assume xres is pow2 otherwise fuck you */
+	tmp = xres - 1;
+	vf->xshift = 0;
+	while(tmp) {
+		++vf->xshift;
+		tmp >>= 1;
+	}
+
+	vptr = vf->v;
+	for(i=0; i<yres; i++) {
+		for(j=0; j<xres; j++) {
+			float x = noise3(j * xfreq, i * yfreq, 0.5);
+			float y = noise3(j * xfreq, i * yfreq, 10.5);
+			float len = sqrt(x * x + y * y);
+			if(len == 0.0f) len = 1.0f;
+
+			vptr->x = x / len;
+			vptr->y = y / len;
+
+			++vptr;
+		}
+	}
+	return 0;
+}
+
+int dump_smktxt_vfield(struct smktxt *stx, const char *fname)
+{
+	FILE *fp;
+	int xsz, ysz;
+
+	if(!(fp = fopen(fname, "wb"))) {
+		fprintf(stderr, "failed to open %s for writing: %s\n", fname, strerror(errno));
+		return -1;
+	}
+
+	xsz = stx->vfield.width;
+	ysz = stx->vfield.height;
+
+	fwrite(&xsz, sizeof xsz, 1, fp);
+	fwrite(&ysz, sizeof ysz, 1, fp);
+	fwrite(stx->vfield.v, sizeof *stx->vfield.v, xsz * ysz, fp);
+	fclose(fp);
+	return 0;
+}
+
+void set_smktxt_wind(struct smktxt *stx, float x, float y, float z)
+{
+	stx->em.wind.x = x;
+	stx->em.wind.y = y;
+	stx->em.wind.z = z;
+}
+
+void set_smktxt_plife(struct smktxt *stx, float life)
+{
+	stx->em.max_life = life;
+}
+
+void set_smktxt_pcount(struct smktxt *stx, int count)
+{
+	free(stx->em.plist);
+	stx->em.plist = 0;
+	stx->em.pcount = count;
+}
+
+void set_smktxt_drag(struct smktxt *stx, float drag)
+{
+	stx->em.drag = drag;
+}
+
+void update_smktxt(struct smktxt *stx, float dt)
+{
+	int i;
+	struct vec2 accel;
+	struct particle *p;
+	struct g3d_vertex *v;
+
+	if(!stx->em.plist) {
+		if(init_emitter(&stx->em, stx->em.pcount, stx->img_pixels, stx->img_xsz, stx->img_ysz) == -1) {
+			fprintf(stderr, "failed to initialize emitter with %d particles\n", stx->em.pcount);
+			return;
+		}
+	}
+
+	p = stx->em.plist;
+	v = stx->em.varr;
+
+	for(i=0; i<stx->em.pcount; i++) {
+		vfield_eval(&stx->vfield, p->x, p->y, &accel);
+		p->x += p->vx * stx->em.drag * dt;
+		p->y += p->vy * stx->em.drag * dt;
+		p->z += p->vz * stx->em.drag * dt;
+		p->vx += (stx->em.wind.x + accel.x * DFL_FORCE) * dt;
+		p->vy += (stx->em.wind.y + accel.y * DFL_FORCE) * dt;
+		p->vz += (stx->em.wind.z + p->z * DFL_ZBIAS) * dt;
+		p->life -= dt;
+		if(p->life < 0.0f) p->life = 0.0f;
+
+		v->x = p->x;
+		v->y = p->y;
+		v->z = p->z;
+		v->w = 1.0f;
+		v->a = cround64(p->life * 255.0f / stx->em.max_life);
+		v->r = 0;
+		v->g = (v->a & 0xe0) >> 3;
+		v->b = (v->a & 0x1f) << 3;
+		++v;
+
+		++p;
+	}
+}
+
+void draw_smktxt(struct smktxt *stx)
+{
+	g3d_draw(G3D_POINTS, stx->em.varr, stx->em.pcount);
+}
+
+
+static int init_emitter(struct emitter *em, int num, unsigned char *map, int xsz, int ysz)
+{
+	int i, x, y;
+	float aspect = (float)xsz / (float)ysz;
+	struct particle *p;
+
+	free(em->varr);
+	if(!(em->varr = malloc(num * sizeof *em->varr))) {
+		fprintf(stderr, "failed to allocate particle vertex array (%d verts)\n", num);
+		return -1;
+	}
+
+	free(em->plist);
+	if(!(em->plist = malloc(num * sizeof *em->plist))) {
+		free(em->varr);
+		return -1;
+	}
+	em->pcount = num;
+
+	p = em->plist;
+	for(i=0; i<num; i++) {
+		do {
+			x = rand() % xsz;
+			y = rand() % ysz;
+		} while(map[y * xsz + x] < 128);
+
+		p->x = (float)x / (float)xsz - 0.5;
+		p->y = -(float)y / (float)xsz + 0.5 / aspect;
+		p->z = ((float)i / (float)num * 2.0 - 1.0) * 0.005;
+		p->r = p->g = p->b = 255;
+		p->vx = p->vy = p->vz = 0.0f;
+		p->life = em->max_life;
+		++p;
+	}
+	return 0;
+}
+
+static int load_vfield(struct vfield *vf, const char *fname)
+{
+	FILE *fp;
+	int tmp;
+
+	if(!(fp = fopen(fname, "rb"))) {
+		fprintf(stderr, "failed to open vector field: %s\n", fname);
+		return -1;
+	}
+	if(fread(&vf->width, sizeof vf->width, 1, fp) < 1 ||
+			fread(&vf->height, sizeof vf->height, 1, fp) < 1) {
+		fprintf(stderr, "load_vfield: unexpected end of file while reading header\n");
+		fclose(fp);
+		return -1;
+	}
+
+	/* assume xsz is pow2 otherwise fuck you */
+	tmp = vf->width - 1;
+	vf->xshift = 0;
+	while(tmp) {
+		++vf->xshift;
+		tmp >>= 1;
+	}
+
+	if(!(vf->v = malloc(vf->width * vf->height * sizeof *vf->v))) {
+		fprintf(stderr, "failed to allocate %dx%d vector field\n", vf->width, vf->height);
+		fclose(fp);
+		return -1;
+	}
+	if(fread(vf->v, sizeof *vf->v, vf->width * vf->height, fp) < vf->width * vf->height) {
+		fprintf(stderr, "load_vfield: unexpected end of file while reading %dx%d vector field\n",
+				vf->width, vf->height);
+		fclose(fp);
+		return -1;
+	}
+	fclose(fp);
+
+	vf->pos.x = vf->pos.y = 0;
+	vf->size.x = vf->size.y = 1;
+	return 0;
+}
+
+static void vfield_eval(struct vfield *vf, float x, float y, struct vec2 *dir)
+{
+	int px, py;
+	float tx, ty;
+	struct vec2 *p1, *p2;
+	struct vec2 left, right;
+
+	x = ((x - vf->pos.x) / vf->size.x + 0.5f) * vf->width;
+	y = ((y - vf->pos.y) / vf->size.y + 0.5f) * vf->height;
+	x = floor(x);
+	y = floor(y);
+
+	if(x < 0) x = 0;
+	if(y < 0) y = 0;
+	if(x > vf->width - 2) x = vf->width - 2;
+	if(y > vf->height - 2) y = vf->height - 2;
+
+	px = (int)x;
+	py = (int)y;
+
+	p1 = vf->v + (py << vf->xshift) + px;
+#ifdef BILERP_FIELD
+	p2 = p1 + vf->width;
+
+	tx = fmod(x, 1.0f);
+	ty = fmod(y, 1.0f);
+
+	left.x = p1->x + (p2->x - p1->x) * ty;
+	left.y = p1->y + (p2->y - p1->y) * ty;
+	++p1;
+	++p2;
+	right.x = p1->x + (p2->x - p1->x) * ty;
+	right.y = p1->y + (p2->y - p1->y) * ty;
+
+	dir->x = left.x + (right.x - left.x) * tx;
+	dir->y = left.y + (right.y - left.y) * ty;
+#else
+	dir->x = p1->x;
+	dir->y = p1->y;
+#endif
+
+#ifdef RANDOMIZE_FIELD
+	dir->x += ((float)rand() / RAND_MAX - 0.5) * RAND_FIELD_MAX;
+	dir->y += ((float)rand() / RAND_MAX - 0.5) * RAND_FIELD_MAX;
+#endif
+}
diff --git a/src/smoketxt.h b/src/smoketxt.h
new file mode 100644
index 0000000..42ba49d
--- /dev/null
+++ b/src/smoketxt.h
@@ -0,0 +1,21 @@
+#ifndef SMOKE_TEXT_H_
+#define SMOKE_TEXT_H_
+
+struct smktxt;
+
+struct smktxt *create_smktxt(const char *imgname, const char *vfieldname);
+void destroy_smktxt(struct smktxt *stx);
+
+int gen_smktxt_vfield(struct smktxt *stx, int xres, int yres, float xfreq, float yfreq);
+int dump_smktxt_vfield(struct smktxt *stx, const char *fname);
+
+void set_smktxt_wind(struct smktxt *stx, float x, float y, float z);
+void set_smktxt_plife(struct smktxt *stx, float life);
+void set_smktxt_pcount(struct smktxt *stx, int count);
+void set_smktxt_drag(struct smktxt *stx, float drag);
+
+void update_smktxt(struct smktxt *stx, float dt);
+void draw_smktxt(struct smktxt *stx);
+
+
+#endif	/* SMOKE_TEXT_ */
-- 
1.7.10.4