1.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     1.2 +++ b/src/mballs.cc	Mon Feb 18 03:46:52 2013 +0200
     1.3 @@ -0,0 +1,140 @@
     1.4 +#include <vector>
     1.5 +#include "opengl.h"
     1.6 +#include "mballs.h"
     1.7 +#include "metasurf.h"
     1.8 +#include "vmath/vmath.h"
     1.9 +
    1.10 +struct MetaBall {
    1.11 +	Vector3 pos;
    1.12 +	float orbit;
    1.13 +	float energy;
    1.14 +	float phase_offs;
    1.15 +};
    1.16 +
    1.17 +#define VOL_SZ			2
    1.18 +#define MBALL_GRID_SZ	50
    1.19 +
    1.20 +static void update();
    1.21 +static float calc_field(float x, float y, float z);
    1.22 +static float eval(float x, float y, float z);
    1.23 +static void vertex(float x, float y, float z);
    1.24 +static void normal(float x, float y, float z);
    1.25 +
    1.26 +static float grid[MBALL_GRID_SZ][MBALL_GRID_SZ][MBALL_GRID_SZ];
    1.27 +static std::vector<MetaBall> balls;
    1.28 +static struct metasurface *msurf;
    1.29 +
    1.30 +bool mball_init()
    1.31 +{
    1.32 +	static const float bbmin = -VOL_SZ / 2.0;
    1.33 +	static const float bbmax = VOL_SZ / 2.0;
    1.34 +
    1.35 +	if(!(msurf = msurf_create())) {
    1.36 +		return false;
    1.37 +	}
    1.38 +	msurf_threshold(msurf, 10);
    1.39 +	msurf_resolution(msurf, MBALL_GRID_SZ, MBALL_GRID_SZ, MBALL_GRID_SZ);
    1.40 +	msurf_bounds(msurf, bbmin, bbmin, bbmin, bbmax, bbmax, bbmax);
    1.41 +	msurf_eval_func(msurf, eval);
    1.42 +	msurf_vertex_func(msurf, vertex);
    1.43 +	//msurf_normal_func(msurf, normal);
    1.44 +
    1.45 +	for(int i=0; i<10; i++) {
    1.46 +		MetaBall mb;
    1.47 +		mb.orbit = 0.25 * rand() / (float)RAND_MAX + 0.35;
    1.48 +		mb.energy = 0.2 * rand() / (float)RAND_MAX + 0.15;
    1.49 +		mb.phase_offs = rand() / (float)RAND_MAX * M_PI * 2.0;
    1.50 +		balls.push_back(mb);
    1.51 +	}
    1.52 +
    1.53 +	return true;
    1.54 +}
    1.55 +
    1.56 +void mball_render()
    1.57 +{
    1.58 +	update();
    1.59 +
    1.60 +	glMatrixMode(GL_MODELVIEW);
    1.61 +	glPushMatrix();
    1.62 +	glScalef(4.0, 4.0, 4.0);
    1.63 +
    1.64 +	glBegin(GL_TRIANGLES);
    1.65 +	msurf_polygonize(msurf);
    1.66 +	glEnd();
    1.67 +
    1.68 +	glPopMatrix();
    1.69 +}
    1.70 +
    1.71 +
    1.72 +static void update()
    1.73 +{
    1.74 +	unsigned int msec = glutGet(GLUT_ELAPSED_TIME);
    1.75 +	float sec = msec / 1000.0;
    1.76 +
    1.77 +	for(size_t i=0; i<balls.size(); i++) {
    1.78 +		float t = sec + balls[i].phase_offs;
    1.79 +		balls[i].pos.x = cos(t * 1.8) * balls[i].orbit;
    1.80 +		balls[i].pos.z = sin(t * 1.2) * balls[i].orbit;
    1.81 +		balls[i].pos.y = (sin(t) + cos(t * 2.0) / 2.0 + sin(t * 3.0) / 3.0) * 0.45;
    1.82 +	}
    1.83 +
    1.84 +	for(int i=0; i<MBALL_GRID_SZ; i++) {
    1.85 +		float x = (float)i / (float)MBALL_GRID_SZ * 2.0 - 1.0;
    1.86 +		for(int j=0; j<MBALL_GRID_SZ; j++) {
    1.87 +			float y = (float)j / (float)MBALL_GRID_SZ * 2.0 - 1.0;
    1.88 +			for(int k=0; k<MBALL_GRID_SZ; k++) {
    1.89 +				float z = (float)k / (float)MBALL_GRID_SZ * 2.0 - 1.0;
    1.90 +				grid[i][j][k] = calc_field(x, y, z);
    1.91 +			}
    1.92 +		}
    1.93 +	}
    1.94 +}
    1.95 +
    1.96 +static float calc_field(float x, float y, float z)
    1.97 +{
    1.98 +	Vector3 pt(x, y, z);
    1.99 +
   1.100 +	float sum = 0.0f;
   1.101 +	for(size_t i=0; i<balls.size(); i++) {
   1.102 +		float dist_sq = (balls[i].pos - pt).length_sq();
   1.103 +		if(dist_sq > 1e-6) {
   1.104 +			sum += balls[i].energy / dist_sq;
   1.105 +		} else {
   1.106 +			sum += 1000.0;
   1.107 +		}
   1.108 +	}
   1.109 +	return sum;
   1.110 +}
   1.111 +
   1.112 +static inline int clamp(int x, int a, int b)
   1.113 +{
   1.114 +	return x < a ? a : (x > b ? b : x);
   1.115 +}
   1.116 +
   1.117 +static float eval(float x, float y, float z)
   1.118 +{
   1.119 +	int cell_x = clamp((int)((x / VOL_SZ + 0.5) * MBALL_GRID_SZ), 0, MBALL_GRID_SZ - 1);
   1.120 +	int cell_y = clamp((int)((y / VOL_SZ + 0.5) * MBALL_GRID_SZ), 0, MBALL_GRID_SZ - 1);
   1.121 +	int cell_z = clamp((int)((z / VOL_SZ + 0.5) * MBALL_GRID_SZ), 0, MBALL_GRID_SZ - 1);
   1.122 +
   1.123 +	return grid[cell_x][cell_y][cell_z];
   1.124 +}
   1.125 +
   1.126 +static void vertex(float x, float y, float z)
   1.127 +{
   1.128 +	float delta = (float)VOL_SZ / (float)MBALL_GRID_SZ;
   1.129 +
   1.130 +	float dfdx = calc_field(x - delta, y, z) - calc_field(x + delta, y, z);
   1.131 +	float dfdy = calc_field(x, y - delta, z) - calc_field(x, y + delta, z);
   1.132 +	float dfdz = calc_field(x, y, z - delta) - calc_field(x, y, z + delta);
   1.133 +
   1.134 +	float len = sqrt(dfdx * dfdx + dfdy * dfdy + dfdz * dfdz);
   1.135 +
   1.136 +	glNormal3f(dfdx / len, dfdy / len, dfdz / len);
   1.137 +	glVertex3f(x, y, z);
   1.138 +}
   1.139 +
   1.140 +static void normal(float x, float y, float z)
   1.141 +{
   1.142 +	glNormal3f(x, y, z);
   1.143 +}