qvolray

changeset 9:a6765984e057

moved the volume loading to volume.c
author John Tsiombikas <nuclear@member.fsf.org>
date Sun, 08 Apr 2012 07:11:54 +0300
parents ae10631bb11b
children 9d2396738b60
files sdr/slice.p.glsl sdr/volray.p.glsl sdr/volray.v.glsl slice.p.glsl src/volray.c src/volume.c src/volume.h volray.p.glsl volray.v.glsl
diffstat 9 files changed, 311 insertions(+), 272 deletions(-) [+]
line diff
     1.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     1.2 +++ b/sdr/slice.p.glsl	Sun Apr 08 07:11:54 2012 +0300
     1.3 @@ -0,0 +1,8 @@
     1.4 +uniform sampler3D volume;
     1.5 +uniform sampler1D xfer_tex;
     1.6 +
     1.7 +void main()
     1.8 +{
     1.9 +	float val = texture3D(volume, gl_TexCoord[0].xyz).a;
    1.10 +	gl_FragColor = vec4(texture1D(xfer_tex, val).xxx, 1.0);
    1.11 +}
     2.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     2.2 +++ b/sdr/volray.p.glsl	Sun Apr 08 07:11:54 2012 +0300
     2.3 @@ -0,0 +1,155 @@
     2.4 +uniform sampler3D volume;
     2.5 +uniform sampler2D ray_tex;
     2.6 +uniform sampler1D xfer_tex;
     2.7 +uniform float ray_step;
     2.8 +uniform float zclip;
     2.9 +
    2.10 +struct Ray {
    2.11 +	vec3 origin, dir;
    2.12 +};
    2.13 +
    2.14 +struct AABBox {
    2.15 +	vec3 min, max;
    2.16 +};
    2.17 +
    2.18 +struct ISect {
    2.19 +	bool hit;
    2.20 +	float t0, t1;
    2.21 +	vec3 pos;
    2.22 +	vec3 normal;
    2.23 +};
    2.24 +
    2.25 +vec3 sky(Ray ray);
    2.26 +vec3 ray_march(Ray ray, float t0, float t1);
    2.27 +vec3 shade(Ray ray, vec3 pos, vec3 norm);
    2.28 +Ray get_primary_ray();
    2.29 +ISect intersect_aabb(Ray ray, AABBox aabb);
    2.30 +
    2.31 +void main()
    2.32 +{
    2.33 +	const AABBox aabb = AABBox(vec3(-1.0, -1.0, -1.0), vec3(1.0, 1.0, 1.0));
    2.34 +	Ray ray = get_primary_ray();
    2.35 +
    2.36 +	vec3 color = vec3(0.0, 0.0, 0.0);
    2.37 +
    2.38 +	ISect res = intersect_aabb(ray, aabb);
    2.39 +	if(res.hit) {
    2.40 +		color = ray_march(ray, res.t0, res.t1);
    2.41 +	}
    2.42 +
    2.43 +	gl_FragColor = vec4(color, 1.0);
    2.44 +}
    2.45 +
    2.46 +float eval(vec3 pos)
    2.47 +{
    2.48 +	vec3 tc = pos * 0.5 + 0.5;
    2.49 +
    2.50 +	if(tc.x < 0.0 || tc.y < 0.0 || tc.z < zclip || tc.x > 1.0 || tc.y > 1.0 || tc.z > 1.0) {
    2.51 +		return 0.0;
    2.52 +	}
    2.53 +
    2.54 +	return texture1D(xfer_tex, texture3D(volume, tc).a).x;
    2.55 +}
    2.56 +
    2.57 +#define OFFS	0.01
    2.58 +vec3 ray_march(Ray ray, float t0, float t1)
    2.59 +{
    2.60 +	float energy = 1.0;
    2.61 +	float t = t0;
    2.62 +	vec3 col = vec3(0.0, 0.0, 0.0);
    2.63 +
    2.64 +
    2.65 +	while(t < t1) {
    2.66 +		vec3 pos = ray.origin + ray.dir * t;
    2.67 +		t += ray_step;
    2.68 +
    2.69 +		float val = eval(pos);
    2.70 +		vec3 norm;
    2.71 +
    2.72 +		norm.x = eval(pos + vec3(OFFS, 0.0, 0.0)) - val;
    2.73 +		norm.y = eval(pos + vec3(0.0, OFFS, 0.0)) - val;
    2.74 +		norm.z = eval(pos + vec3(0.0, 0.0, OFFS)) - val;
    2.75 +
    2.76 +		col += shade(ray, pos, normalize(norm)) * val * energy;
    2.77 +		energy -= val;
    2.78 +		if(energy < 0.001) {
    2.79 +			break;
    2.80 +		}
    2.81 +		pos += ray.dir * ray_step;
    2.82 +	}
    2.83 +
    2.84 +	return col;
    2.85 +}
    2.86 +
    2.87 +vec3 shade(Ray ray, vec3 pos, vec3 norm)
    2.88 +{
    2.89 +	vec3 ldir = -pos;//normalize(vec3(10.0, 10.0, -10.0) - pos);
    2.90 +	vec3 vdir = -ray.dir;
    2.91 +	vec3 hdir = normalize(ldir + vdir);
    2.92 +
    2.93 +	float ndotl = dot(ldir, norm);
    2.94 +	float ndoth = dot(hdir, norm);
    2.95 +
    2.96 +	vec3 dcol = vec3(0.9, 0.9, 0.9) * max(ndotl, 0.0);
    2.97 +	vec3 scol = vec3(0.5, 0.5, 0.5) * pow(max(ndoth, 0.0), 50.0);
    2.98 +
    2.99 +	return vec3(0.01, 0.01, 0.01) + dcol + scol;
   2.100 +}
   2.101 +
   2.102 +Ray get_primary_ray()
   2.103 +{
   2.104 +	Ray ray;
   2.105 +	vec2 tc = gl_TexCoord[0].xy;
   2.106 +	ray.dir = gl_NormalMatrix * normalize(texture2D(ray_tex, tc).xyz);
   2.107 +	ray.origin = (gl_ModelViewMatrix * vec4(0.0, 0.0, 0.0, 1.0)).xyz;
   2.108 +	return ray;
   2.109 +}
   2.110 +
   2.111 +ISect intersect_aabb(Ray ray, AABBox aabb)
   2.112 +{
   2.113 +	ISect res;
   2.114 +	res.hit = false;
   2.115 +	res.t0 = res.t1 = 0.0;
   2.116 +
   2.117 +	if(ray.origin.x >= aabb.min.x && ray.origin.y >= aabb.min.y && ray.origin.z >= aabb.min.z &&
   2.118 +			ray.origin.x < aabb.max.x && ray.origin.y < aabb.max.y && ray.origin.z < aabb.max.z) {
   2.119 +		res.hit = true;
   2.120 +		return res;
   2.121 +	}
   2.122 +
   2.123 +	vec4 bbox[2];
   2.124 +	bbox[0] = vec4(aabb.min.x, aabb.min.y, aabb.min.z, 0);
   2.125 +	bbox[1] = vec4(aabb.max.x, aabb.max.y, aabb.max.z, 0);
   2.126 +
   2.127 +	int xsign = int(ray.dir.x < 0.0);
   2.128 +	float invdirx = 1.0 / ray.dir.x;
   2.129 +	float tmin = (bbox[xsign].x - ray.origin.x) * invdirx;
   2.130 +	float tmax = (bbox[1 - xsign].x - ray.origin.x) * invdirx;
   2.131 +
   2.132 +	int ysign = int(ray.dir.y < 0.0);
   2.133 +	float invdiry = 1.0 / ray.dir.y;
   2.134 +	float tymin = (bbox[ysign].y - ray.origin.y) * invdiry;
   2.135 +	float tymax = (bbox[1 - ysign].y - ray.origin.y) * invdiry;
   2.136 +
   2.137 +	if(tmin > tymax || tymin > tmax) {
   2.138 +		return res;
   2.139 +	}
   2.140 +
   2.141 +	if(tymin > tmin) tmin = tymin;
   2.142 +	if(tymax < tmax) tmax = tymax;
   2.143 +
   2.144 +	int zsign = int(ray.dir.z < 0.0);
   2.145 +	float invdirz = 1.0 / ray.dir.z;
   2.146 +	float tzmin = (bbox[zsign].z - ray.origin.z) * invdirz;
   2.147 +	float tzmax = (bbox[1 - zsign].z - ray.origin.z) * invdirz;
   2.148 +
   2.149 +	if(tmin > tzmax || tzmin > tmax) {
   2.150 +		return res;
   2.151 +	}
   2.152 +
   2.153 +	res.t0 = tmin;
   2.154 +	res.t1 = tmax;
   2.155 +	res.hit = true;
   2.156 +	return res;
   2.157 +}
   2.158 +
     3.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     3.2 +++ b/sdr/volray.v.glsl	Sun Apr 08 07:11:54 2012 +0300
     3.3 @@ -0,0 +1,5 @@
     3.4 +void main()
     3.5 +{
     3.6 +	gl_Position = gl_Vertex;
     3.7 +	gl_TexCoord[0] = gl_TextureMatrix[0] * gl_MultiTexCoord0;
     3.8 +}
     4.1 --- a/slice.p.glsl	Sat Apr 07 16:07:12 2012 +0300
     4.2 +++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
     4.3 @@ -1,8 +0,0 @@
     4.4 -uniform sampler3D volume;
     4.5 -uniform sampler1D xfer_tex;
     4.6 -
     4.7 -void main()
     4.8 -{
     4.9 -	float val = texture3D(volume, gl_TexCoord[0].xyz).x;
    4.10 -	gl_FragColor = vec4(texture1D(xfer_tex, val).xxx, 1.0);
    4.11 -}
     5.1 --- a/src/volray.c	Sat Apr 07 16:07:12 2012 +0300
     5.2 +++ b/src/volray.c	Sun Apr 08 07:11:54 2012 +0300
     5.3 @@ -12,14 +12,10 @@
     5.4  #include <vmath/vmath.h>
     5.5  #include <imago2.h>
     5.6  #include "sdr.h"
     5.7 +#include "volume.h"
     5.8  
     5.9  #define XFER_MAP_SZ	512
    5.10  
    5.11 -struct slice_file {
    5.12 -	char *name;
    5.13 -	struct slice_file *next;
    5.14 -};
    5.15 -
    5.16  enum {
    5.17  	UIMODE_DEFAULT,
    5.18  	UIMODE_XFER,
    5.19 @@ -43,24 +39,26 @@
    5.20  static int round_pow2(int x);
    5.21  static void create_transfer_map(float mean, float sdev);
    5.22  
    5.23 -float cam_theta = 0, cam_phi = 0, cam_dist = 4.0;
    5.24 -float cam_x, cam_y, cam_z;
    5.25 +static float cam_theta = 0, cam_phi = 0, cam_dist = 4.0;
    5.26 +static float cam_x, cam_y, cam_z;
    5.27  
    5.28 -vec2_t tex_scale;
    5.29 -struct slice_file *flist;
    5.30 -int nslices;
    5.31 -unsigned int vol_sdr, slice_sdr, vol_tex, ray_tex;
    5.32 -int win_xsz, win_ysz;
    5.33 -int raytex_needs_recalc = 1;
    5.34 +static vec2_t tex_scale;
    5.35 +static struct volume *volume;
    5.36 +static unsigned int vol_sdr, slice_sdr, ray_tex;
    5.37 +static int win_xsz, win_ysz;
    5.38 +static int raytex_needs_recalc = 1;
    5.39  
    5.40 -unsigned int xfer_tex;
    5.41 -float xfer_mean = 0.7, xfer_sdev = 0.1;
    5.42 -int xfertex_needs_recalc = 1;
    5.43 +static unsigned int xfer_tex;
    5.44 +static float xfer_mean = 0.7, xfer_sdev = 0.1;
    5.45 +static int xfertex_needs_recalc = 1;
    5.46  
    5.47  static int uimode;
    5.48  static float cur_z = 0.0;
    5.49  static float ray_step = 0.01;
    5.50  
    5.51 +static const char *fname;
    5.52 +
    5.53 +
    5.54  int main(int argc, char **argv)
    5.55  {
    5.56  	glutInit(&argc, argv);
    5.57 @@ -92,9 +90,7 @@
    5.58  
    5.59  int init(void)
    5.60  {
    5.61 -	int i, vol_xsz, vol_ysz;
    5.62 -
    5.63 -	if(!(vol_sdr = create_program_load("volray.v.glsl", "volray.p.glsl"))) {
    5.64 +	if(!(vol_sdr = create_program_load("sdr/volray.v.glsl", "sdr/volray.p.glsl"))) {
    5.65  		return -1;
    5.66  	}
    5.67  	set_uniform_int(vol_sdr, "volume", 0);
    5.68 @@ -103,51 +99,16 @@
    5.69  	set_uniform_float(vol_sdr, "ray_step", ray_step);
    5.70  	set_uniform_float(vol_sdr, "zclip", cur_z);
    5.71  
    5.72 -	if(!(slice_sdr = create_program_load(0, "slice.p.glsl"))) {
    5.73 +	if(!(slice_sdr = create_program_load(0, "sdr/slice.p.glsl"))) {
    5.74  		return -1;
    5.75  	}
    5.76  	set_uniform_int(slice_sdr, "volume", 0);
    5.77  	set_uniform_int(slice_sdr, "xfer_tex", 1);
    5.78  
    5.79 -	glGenTextures(1, &vol_tex);
    5.80 -	glBindTexture(GL_TEXTURE_3D, vol_tex);
    5.81 -	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MIN_FILTER, GL_LINEAR);
    5.82 -	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MAG_FILTER, GL_LINEAR);
    5.83 -	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE);
    5.84 -	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE);
    5.85 -	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_R, GL_CLAMP_TO_EDGE);
    5.86 +	if(!(volume = load_volume(fname))) {
    5.87 +		return -1;
    5.88 +	}
    5.89  
    5.90 -	for(i=0; i<nslices; i++) {
    5.91 -		int xsz, ysz;
    5.92 -		void *pix;
    5.93 -		struct slice_file *sfile = flist;
    5.94 -		flist = flist->next;
    5.95 -
    5.96 -		if(!(pix = img_load_pixels(sfile->name, &xsz, &ysz, IMG_FMT_RGBA32))) {
    5.97 -			fprintf(stderr, "failed to load image: %s\n", sfile->name);
    5.98 -			return -1;
    5.99 -		}
   5.100 -
   5.101 -		if(i == 0) {
   5.102 -			/* allocate storage for the texture */
   5.103 -			glTexImage3D(GL_TEXTURE_3D, 0, GL_RGBA32F, xsz, ysz, nslices, 0, GL_RGBA, GL_UNSIGNED_BYTE, 0);
   5.104 -
   5.105 -			vol_xsz = xsz;
   5.106 -			vol_ysz = ysz;
   5.107 -
   5.108 -		} else {
   5.109 -			if(xsz != vol_xsz || ysz != vol_ysz) {
   5.110 -				fprintf(stderr, "%s: inconsistent slice size: %dx%d. expected: %dx%d\n",
   5.111 -						sfile->name, xsz, ysz, vol_xsz, vol_ysz);
   5.112 -				img_free_pixels(pix);
   5.113 -				return -1;
   5.114 -			}
   5.115 -		}
   5.116 -		free(sfile);
   5.117 -
   5.118 -		glTexSubImage3D(GL_TEXTURE_3D, 0, 0, 0, i, xsz, ysz, 1, GL_RGBA, GL_UNSIGNED_BYTE, pix);
   5.119 -		img_free_pixels(pix);
   5.120 -	}
   5.121  	return 0;
   5.122  }
   5.123  
   5.124 @@ -190,7 +151,7 @@
   5.125  
   5.126  	/* tex unit0: volume data 3D texture */
   5.127  	glActiveTexture(GL_TEXTURE0);
   5.128 -	glBindTexture(GL_TEXTURE_3D, vol_tex);
   5.129 +	glBindTexture(GL_TEXTURE_3D, volume->tex_vol);
   5.130  	glEnable(GL_TEXTURE_3D);
   5.131  
   5.132  	/* tex unit1: primary rays in view space */
   5.133 @@ -235,7 +196,7 @@
   5.134  	glTranslatef(-1, -1, 0);
   5.135  
   5.136  	glActiveTexture(GL_TEXTURE0);
   5.137 -	glBindTexture(GL_TEXTURE_3D, vol_tex);
   5.138 +	glBindTexture(GL_TEXTURE_3D, volume->tex_vol);
   5.139  	glEnable(GL_TEXTURE_3D);
   5.140  
   5.141  	glActiveTexture(GL_TEXTURE1);
   5.142 @@ -426,7 +387,6 @@
   5.143  int parse_args(int argc, char **argv)
   5.144  {
   5.145  	int i;
   5.146 -	struct slice_file *tail;
   5.147  	char *endp;
   5.148  
   5.149  	for(i=1; i<argc; i++) {
   5.150 @@ -453,27 +413,16 @@
   5.151  				return -1;
   5.152  			}
   5.153  		} else {
   5.154 -			struct slice_file *sfile;
   5.155 -
   5.156 -			if(!(sfile = malloc(sizeof *sfile))) {
   5.157 -				perror("failed to allocate memory");
   5.158 +			if(fname) {
   5.159 +				fprintf(stderr, "unexpected argument: %s\n", argv[i]);
   5.160  				return -1;
   5.161  			}
   5.162 -			sfile->name = argv[i];
   5.163 -			sfile->next = 0;
   5.164 -
   5.165 -			if(!flist) {
   5.166 -				flist = tail = sfile;
   5.167 -			} else {
   5.168 -				tail->next = sfile;
   5.169 -				tail = sfile;
   5.170 -			}
   5.171 -			nslices++;
   5.172 +			fname = argv[i];
   5.173  		}
   5.174  	}
   5.175  
   5.176 -	if(!nslices) {
   5.177 -		fprintf(stderr, "pass the slice filenames\n");
   5.178 +	if(!fname) {
   5.179 +		fprintf(stderr, "pass the volume descriptor filename\n");
   5.180  		return -1;
   5.181  	}
   5.182  	return 0;
     6.1 --- a/src/volume.c	Sat Apr 07 16:07:12 2012 +0300
     6.2 +++ b/src/volume.c	Sun Apr 08 07:11:54 2012 +0300
     6.3 @@ -10,19 +10,24 @@
     6.4  #endif
     6.5  
     6.6  #include <imago2.h>
     6.7 +#include "volume.h"
     6.8  
     6.9  struct slice {
    6.10  	char *name;
    6.11  	struct slice *next;
    6.12  };
    6.13  
    6.14 +static void calc_gradients(float *voxels, int xsz, int ysz, int zsz);
    6.15  static struct slice *read_voldesc(const char *fname, struct volume *vol);
    6.16  static char *trim(char *str);
    6.17  
    6.18  struct volume *load_volume(const char *fname)
    6.19  {
    6.20 +	int i;
    6.21  	struct slice *slist;
    6.22  	struct volume *vol = 0;
    6.23 +	float *voxels = 0, *vptr;
    6.24 +	struct img_pixmap img;
    6.25  
    6.26  	if(!(vol = malloc(sizeof *vol))) {
    6.27  		perror("failed to allocate volume");
    6.28 @@ -35,6 +40,60 @@
    6.29  		goto err;
    6.30  	}
    6.31  
    6.32 +	/* load the first image to determine slice dimensions */
    6.33 +	img_init(&img);
    6.34 +	if(img_load(&img, slist->name) == -1) {
    6.35 +		fprintf(stderr, "failed to load volume slice: %s\n", slist->name);
    6.36 +		goto err;
    6.37 +	}
    6.38 +
    6.39 +	vol->sz[0] = img.width;
    6.40 +	vol->sz[1] = img.height;
    6.41 +	printf("volume dimensions: %dx%dx%d\n", img.width, img.height, vol->sz[2]);
    6.42 +
    6.43 +	/* allocate the whole volume at once */
    6.44 +	if(!(voxels = malloc(vol->sz[0] * vol->sz[1] * vol->sz[2] * 4 * sizeof *voxels))) {
    6.45 +		img_destroy(&img);
    6.46 +		goto err;
    6.47 +	}
    6.48 +	vptr = voxels;
    6.49 +	img_destroy(&img);
    6.50 +
    6.51 +	/* put the volume data into the alpha component */
    6.52 +	for(i=0; i<vol->sz[2]; i++) {
    6.53 +		int x, y, xsz, ysz;
    6.54 +		float *pixels, *src;
    6.55 +		struct slice *slice = slist;
    6.56 +		slist = slist->next;
    6.57 +
    6.58 +		if(!(pixels = img_load_pixels(slice->name, &xsz, &ysz, IMG_FMT_GREYF))) {
    6.59 +			fprintf(stderr, "failed to load volume slice: %s\n", slice->name);
    6.60 +			goto err;
    6.61 +		}
    6.62 +
    6.63 +		if(xsz != vol->sz[0] || ysz != vol->sz[1]) {
    6.64 +			fprintf(stderr, "inconsistent dimensions for slice %s: %dx%d (%dx%d expected)\n",
    6.65 +					slice->name, xsz, ysz, vol->sz[0], vol->sz[1]);
    6.66 +			goto err;
    6.67 +		}
    6.68 +		free(slice->name);
    6.69 +		free(slice);
    6.70 +
    6.71 +		src = pixels;
    6.72 +		for(y=0; y<ysz; y++) {
    6.73 +			for(x=0; x<xsz; x++) {
    6.74 +				vptr[0] = vptr[1] = vptr[2] = 0.0f;
    6.75 +				vptr[3] = *src++;
    6.76 +				vptr += 4;
    6.77 +			}
    6.78 +		}
    6.79 +		img_free_pixels(pixels);
    6.80 +	}
    6.81 +
    6.82 +	/* calculate gradients */
    6.83 +	/*calc_gradients(voxels, vol->sz[0], vol->sz[1], vol->sz[2]);*/
    6.84 +
    6.85 +	/* create the volume texture */
    6.86  	glGenTextures(1, &vol->tex_vol);
    6.87  	glBindTexture(GL_TEXTURE_3D, vol->tex_vol);
    6.88  	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MIN_FILTER, GL_LINEAR);
    6.89 @@ -42,26 +101,12 @@
    6.90  	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE);
    6.91  	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE);
    6.92  	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_R, GL_CLAMP_TO_EDGE);
    6.93 +	glTexImage3D(GL_TEXTURE_3D, 0, GL_RGBA32F, vol->sz[0], vol->sz[1], vol->sz[2], 0,
    6.94 +			GL_RGBA, GL_FLOAT, voxels);
    6.95  
    6.96 -	for(i=0; i<vol->sz[2]; i++) {
    6.97 -		struct img_pixmap img;
    6.98 -		struct slice *slice = head;
    6.99 -		head = head->next;
   6.100 -
   6.101 -		img_init(&img);
   6.102 -		if(img_load(&img, slice->name) == -1) {
   6.103 -			fprintf(stderr, "failed to load volume slice: %s\n", slice->name);
   6.104 -			goto err;
   6.105 -		}
   6.106 -	}
   6.107 -
   6.108 -	glGenTextures(1, &vol->tex_grad);
   6.109 -	glBindTexture(GL_TEXTURE_3D, vol->tex_grad);
   6.110 -	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MIN_FILTER, GL_LINEAR);
   6.111 -	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MAG_FILTER, GL_LINEAR);
   6.112 -	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE);
   6.113 -	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE);
   6.114 -	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_R, GL_CLAMP_TO_EDGE);
   6.115 +	/* ... and destroy our copy */
   6.116 +	free(voxels);
   6.117 +	return vol;
   6.118  
   6.119  err:
   6.120  	while(slist) {
   6.121 @@ -70,20 +115,54 @@
   6.122  		free(tmp->name);
   6.123  		free(tmp);
   6.124  	}
   6.125 +	free_volume(vol);
   6.126 +	free(voxels);
   6.127 +	return 0;
   6.128 +}
   6.129 +
   6.130 +void free_volume(struct volume *vol)
   6.131 +{
   6.132  	if(vol) {
   6.133 -		if(vol->tex_vol)
   6.134 +		if(vol->tex_vol) {
   6.135  			glDeleteTextures(1, &vol->tex_vol);
   6.136 -		if(vol->tex_grad)
   6.137 -			glDeleteTextures(1, &vol->tex_grad);
   6.138 +		}
   6.139  		free(vol);
   6.140  	}
   6.141 -	return 0;
   6.142 +}
   6.143 +
   6.144 +struct vec4 { float x, y, z, w; };
   6.145 +
   6.146 +static void calc_gradients(float *voxels, int xsz, int ysz, int zsz)
   6.147 +{
   6.148 +	int x, y, z, slice_pixels = xsz * ysz;
   6.149 +	struct vec4 *vptr = (struct vec4*)voxels;
   6.150 +
   6.151 +	for(z=0; z<zsz; z++) {
   6.152 +		for(y=0; y<ysz; y++) {
   6.153 +			for(x=0; x<xsz; x++) {
   6.154 +				float x0, x1, y0, y1, z0, z1;
   6.155 +
   6.156 +				x0 = x > 0 ? (vptr - 1)->w : vptr->w;
   6.157 +				y0 = y > 0 ? (vptr - xsz)->w : vptr->w;
   6.158 +				z0 = z > 0 ? (vptr - slice_pixels)->w : vptr->w;
   6.159 +				x1 = x < xsz - 1 ? (vptr + 1)->w : vptr->w;
   6.160 +				y1 = y < ysz - 1 ? (vptr + xsz)->w : vptr->w;
   6.161 +				z1 = z < zsz - 1 ? (vptr + slice_pixels)->w : vptr->w;
   6.162 +
   6.163 +				vptr->x = x1 - x0;
   6.164 +				vptr->y = y1 - y0;
   6.165 +				vptr->z = z1 - z0;
   6.166 +
   6.167 +				vptr++;
   6.168 +			}
   6.169 +		}
   6.170 +	}
   6.171  }
   6.172  
   6.173  static struct slice *read_voldesc(const char *fname, struct volume *vol)
   6.174  {
   6.175  	FILE *fp = 0;
   6.176 -	char buf[512];
   6.177 +	char buf[512], *slash, *prefix = "";
   6.178  	int mode_slices = 0;
   6.179  	struct slice *head = 0, *tail;
   6.180  
   6.181 @@ -97,6 +176,13 @@
   6.182  		goto err;
   6.183  	}
   6.184  
   6.185 +	if((slash = strrchr(fname, '/'))) {
   6.186 +		int len = slash - fname + 1;
   6.187 +		prefix = alloca(len + 1);
   6.188 +		memcpy(prefix, fname, len);
   6.189 +		prefix[len] = 0;
   6.190 +	}
   6.191 +
   6.192  	while(fgets(buf, sizeof buf, fp)) {
   6.193  		char *line = trim(buf);
   6.194  
   6.195 @@ -106,12 +192,12 @@
   6.196  		if(mode_slices) {
   6.197  			/* we're in the part that contains filenames... append to the list */
   6.198  			struct slice *node = malloc(sizeof *node);
   6.199 -			if(!node || !(node->name = malloc(strlen(line) + 1))) {
   6.200 +			if(!node || !(node->name = malloc(strlen(prefix) + strlen(line) + 1))) {
   6.201  				perror("failed to allocate list node");
   6.202  				free(node);
   6.203  				goto err;
   6.204  			}
   6.205 -			strcpy(node->name, line);
   6.206 +			sprintf(node->name, "%s%s", prefix, line);
   6.207  			node->next = 0;
   6.208  
   6.209  			if(head) {
   6.210 @@ -123,6 +209,10 @@
   6.211  			vol->sz[2]++;
   6.212  		} else {
   6.213  			/* TODO we're in the options part... parse */
   6.214 +
   6.215 +			if(strcmp(line, "SLICES") == 0) {
   6.216 +				mode_slices = 1;
   6.217 +			}
   6.218  		}
   6.219  	}
   6.220  	fclose(fp);
     7.1 --- a/src/volume.h	Sat Apr 07 16:07:12 2012 +0300
     7.2 +++ b/src/volume.h	Sun Apr 08 07:11:54 2012 +0300
     7.3 @@ -6,7 +6,7 @@
     7.4  struct volume {
     7.5  	int sz[3];
     7.6  	float zaspect;
     7.7 -	unsigned int tex_vol, tex_grad;
     7.8 +	unsigned int tex_vol;
     7.9  };
    7.10  
    7.11  struct volume *load_volume(const char *fname);
     8.1 --- a/volray.p.glsl	Sat Apr 07 16:07:12 2012 +0300
     8.2 +++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
     8.3 @@ -1,155 +0,0 @@
     8.4 -uniform sampler3D volume;
     8.5 -uniform sampler2D ray_tex;
     8.6 -uniform sampler1D xfer_tex;
     8.7 -uniform float ray_step;
     8.8 -uniform float zclip;
     8.9 -
    8.10 -struct Ray {
    8.11 -	vec3 origin, dir;
    8.12 -};
    8.13 -
    8.14 -struct AABBox {
    8.15 -	vec3 min, max;
    8.16 -};
    8.17 -
    8.18 -struct ISect {
    8.19 -	bool hit;
    8.20 -	float t0, t1;
    8.21 -	vec3 pos;
    8.22 -	vec3 normal;
    8.23 -};
    8.24 -
    8.25 -vec3 sky(Ray ray);
    8.26 -vec3 ray_march(Ray ray, float t0, float t1);
    8.27 -vec3 shade(Ray ray, vec3 pos, vec3 norm);
    8.28 -Ray get_primary_ray();
    8.29 -ISect intersect_aabb(Ray ray, AABBox aabb);
    8.30 -
    8.31 -void main()
    8.32 -{
    8.33 -	const AABBox aabb = AABBox(vec3(-1.0, -1.0, -1.0), vec3(1.0, 1.0, 1.0));
    8.34 -	Ray ray = get_primary_ray();
    8.35 -
    8.36 -	vec3 color = vec3(0.0, 0.0, 0.0);
    8.37 -
    8.38 -	ISect res = intersect_aabb(ray, aabb);
    8.39 -	if(res.hit) {
    8.40 -		color = ray_march(ray, res.t0, res.t1);
    8.41 -	}
    8.42 -
    8.43 -	gl_FragColor = vec4(color, 1.0);
    8.44 -}
    8.45 -
    8.46 -float eval(vec3 pos)
    8.47 -{
    8.48 -	vec3 tc = pos * 0.5 + 0.5;
    8.49 -
    8.50 -	if(tc.x < 0.0 || tc.y < 0.0 || tc.z < zclip || tc.x > 1.0 || tc.y > 1.0 || tc.z > 1.0) {
    8.51 -		return 0.0;
    8.52 -	}
    8.53 -
    8.54 -	return texture1D(xfer_tex, texture3D(volume, tc).x).x;
    8.55 -}
    8.56 -
    8.57 -#define OFFS	0.01
    8.58 -vec3 ray_march(Ray ray, float t0, float t1)
    8.59 -{
    8.60 -	float energy = 1.0;
    8.61 -	float t = t0;
    8.62 -	vec3 col = vec3(0.0, 0.0, 0.0);
    8.63 -
    8.64 -
    8.65 -	while(t < t1) {
    8.66 -		vec3 pos = ray.origin + ray.dir * t;
    8.67 -		t += ray_step;
    8.68 -
    8.69 -		float val = eval(pos);
    8.70 -		vec3 norm;
    8.71 -
    8.72 -		norm.x = eval(pos + vec3(OFFS, 0.0, 0.0)) - val;
    8.73 -		norm.y = eval(pos + vec3(0.0, OFFS, 0.0)) - val;
    8.74 -		norm.z = eval(pos + vec3(0.0, 0.0, OFFS)) - val;
    8.75 -
    8.76 -		col += shade(ray, pos, normalize(norm)) * val * energy;
    8.77 -		energy -= val;
    8.78 -		if(energy < 0.001) {
    8.79 -			break;
    8.80 -		}
    8.81 -		pos += ray.dir * ray_step;
    8.82 -	}
    8.83 -
    8.84 -	return col;
    8.85 -}
    8.86 -
    8.87 -vec3 shade(Ray ray, vec3 pos, vec3 norm)
    8.88 -{
    8.89 -	vec3 ldir = -pos;//normalize(vec3(10.0, 10.0, -10.0) - pos);
    8.90 -	vec3 vdir = -ray.dir;
    8.91 -	vec3 hdir = normalize(ldir + vdir);
    8.92 -
    8.93 -	float ndotl = dot(ldir, norm);
    8.94 -	float ndoth = dot(hdir, norm);
    8.95 -
    8.96 -	vec3 dcol = vec3(0.9, 0.9, 0.9) * max(ndotl, 0.0);
    8.97 -	vec3 scol = vec3(0.5, 0.5, 0.5) * pow(max(ndoth, 0.0), 50.0);
    8.98 -
    8.99 -	return vec3(0.01, 0.01, 0.01) + dcol + scol;
   8.100 -}
   8.101 -
   8.102 -Ray get_primary_ray()
   8.103 -{
   8.104 -	Ray ray;
   8.105 -	vec2 tc = gl_TexCoord[0].xy;
   8.106 -	ray.dir = gl_NormalMatrix * normalize(texture2D(ray_tex, tc).xyz);
   8.107 -	ray.origin = (gl_ModelViewMatrix * vec4(0.0, 0.0, 0.0, 1.0)).xyz;
   8.108 -	return ray;
   8.109 -}
   8.110 -
   8.111 -ISect intersect_aabb(Ray ray, AABBox aabb)
   8.112 -{
   8.113 -	ISect res;
   8.114 -	res.hit = false;
   8.115 -	res.t0 = res.t1 = 0.0;
   8.116 -
   8.117 -	if(ray.origin.x >= aabb.min.x && ray.origin.y >= aabb.min.y && ray.origin.z >= aabb.min.z &&
   8.118 -			ray.origin.x < aabb.max.x && ray.origin.y < aabb.max.y && ray.origin.z < aabb.max.z) {
   8.119 -		res.hit = true;
   8.120 -		return res;
   8.121 -	}
   8.122 -
   8.123 -	vec4 bbox[2];
   8.124 -	bbox[0] = vec4(aabb.min.x, aabb.min.y, aabb.min.z, 0);
   8.125 -	bbox[1] = vec4(aabb.max.x, aabb.max.y, aabb.max.z, 0);
   8.126 -
   8.127 -	int xsign = int(ray.dir.x < 0.0);
   8.128 -	float invdirx = 1.0 / ray.dir.x;
   8.129 -	float tmin = (bbox[xsign].x - ray.origin.x) * invdirx;
   8.130 -	float tmax = (bbox[1 - xsign].x - ray.origin.x) * invdirx;
   8.131 -
   8.132 -	int ysign = int(ray.dir.y < 0.0);
   8.133 -	float invdiry = 1.0 / ray.dir.y;
   8.134 -	float tymin = (bbox[ysign].y - ray.origin.y) * invdiry;
   8.135 -	float tymax = (bbox[1 - ysign].y - ray.origin.y) * invdiry;
   8.136 -
   8.137 -	if(tmin > tymax || tymin > tmax) {
   8.138 -		return res;
   8.139 -	}
   8.140 -
   8.141 -	if(tymin > tmin) tmin = tymin;
   8.142 -	if(tymax < tmax) tmax = tymax;
   8.143 -
   8.144 -	int zsign = int(ray.dir.z < 0.0);
   8.145 -	float invdirz = 1.0 / ray.dir.z;
   8.146 -	float tzmin = (bbox[zsign].z - ray.origin.z) * invdirz;
   8.147 -	float tzmax = (bbox[1 - zsign].z - ray.origin.z) * invdirz;
   8.148 -
   8.149 -	if(tmin > tzmax || tzmin > tmax) {
   8.150 -		return res;
   8.151 -	}
   8.152 -
   8.153 -	res.t0 = tmin;
   8.154 -	res.t1 = tmax;
   8.155 -	res.hit = true;
   8.156 -	return res;
   8.157 -}
   8.158 -
     9.1 --- a/volray.v.glsl	Sat Apr 07 16:07:12 2012 +0300
     9.2 +++ /dev/null	Thu Jan 01 00:00:00 1970 +0000
     9.3 @@ -1,5 +0,0 @@
     9.4 -void main()
     9.5 -{
     9.6 -	gl_Position = gl_Vertex;
     9.7 -	gl_TexCoord[0] = gl_TextureMatrix[0] * gl_MultiTexCoord0;
     9.8 -}