nuclear@3: #include <stdio.h>
nuclear@3: #include <stdlib.h>
nuclear@3: #include <math.h>
nuclear@3: #include <assert.h>
nuclear@3: 
nuclear@3: #ifndef NO_SHADERS
nuclear@3: #include <GL/glew.h>
nuclear@3: #include "sdr.h"
nuclear@3: #endif
nuclear@3: 
nuclear@3: #ifndef __APPLE__
nuclear@3: #include <GL/glut.h>
nuclear@3: #else
nuclear@3: #include <GLUT/glut.h>
nuclear@3: #endif
nuclear@3: 
nuclear@3: #include <imago2.h>
nuclear@3: 
nuclear@3: #include "cam.h"
nuclear@3: #include "metasurf.h"
nuclear@3: 
nuclear@3: float eval(float x, float y, float z);
nuclear@3: void vertex(float x, float y, float z);
nuclear@3: void render(void);
nuclear@3: void disp(void);
nuclear@3: void reshape(int x, int y);
nuclear@3: void keyb(unsigned char key, int x, int y);
nuclear@11: void keyb_up(unsigned char key, int x, int y);
nuclear@3: void mouse(int bn, int state, int x, int y);
nuclear@3: void motion(int x, int y);
nuclear@3: int parse_args(int argc, char **argv);
nuclear@3: 
nuclear@3: int stereo, fullscreen;
nuclear@3: int orig_xsz, orig_ysz;
nuclear@3: 
nuclear@3: struct metasurface *msurf;
nuclear@3: float threshold = 0.5;
nuclear@3: #ifndef NO_SHADERS
nuclear@3: unsigned int sdr;
nuclear@3: #endif
nuclear@3: 
nuclear@11: float yscale = 1.0;
nuclear@11: 
nuclear@3: struct img_pixmap *volume;
nuclear@3: int xres, yres, num_slices;
nuclear@3: 
nuclear@3: int dlist, need_update = 1;
nuclear@3: 
nuclear@3: int main(int argc, char **argv)
nuclear@3: {
nuclear@3: 	float amb[] = {0, 0, 0, 0};
nuclear@3: 	float lpos[] = {-0.2, 0.2, 1, 0};
nuclear@3: 
nuclear@3: 	glutInitWindowSize(1280, 720);
nuclear@3: 	glutInit(&argc, argv);
nuclear@3: 
nuclear@3: 	if(parse_args(argc, argv) == -1) {
nuclear@3: 		return 1;
nuclear@3: 	}
nuclear@3: 	glutInitDisplayMode(GLUT_RGB | GLUT_DEPTH | GLUT_DOUBLE | (stereo ? GLUT_STEREO : 0));
nuclear@3: 	glutCreateWindow("metasurf - volume rendering");
nuclear@3: 
nuclear@3: 	orig_xsz = glutGet(GLUT_WINDOW_WIDTH);
nuclear@3: 	orig_ysz = glutGet(GLUT_WINDOW_HEIGHT);
nuclear@3: 
nuclear@3: 	if(fullscreen) {
nuclear@3: 		glutFullScreen();
nuclear@3: 	}
nuclear@3: 
nuclear@3: 	glutDisplayFunc(disp);
nuclear@3: 	glutReshapeFunc(reshape);
nuclear@3: 	glutKeyboardFunc(keyb);
nuclear@11: 	glutKeyboardUpFunc(keyb_up);
nuclear@3: 	glutMouseFunc(mouse);
nuclear@3: 	glutMotionFunc(motion);
nuclear@3: 
nuclear@3: #ifndef NO_SHADERS
nuclear@3: 	glewInit();
nuclear@3: 	if(!(sdr = create_program_load("sdr/vert.glsl", "sdr/frag.glsl"))) {
nuclear@3: 		return 1;
nuclear@3: 	}
nuclear@3: #endif
nuclear@3: 
nuclear@3: 	glEnable(GL_CULL_FACE);
nuclear@3: 	glEnable(GL_DEPTH_TEST);
nuclear@3: 
nuclear@3: 	glLightModelfv(GL_LIGHT_MODEL_AMBIENT, amb);
nuclear@3: 
nuclear@3: 	glEnable(GL_LIGHTING);
nuclear@3: 	glEnable(GL_LIGHT0);
nuclear@3: 	glLightfv(GL_LIGHT0, GL_POSITION, lpos);
nuclear@3: 
nuclear@3: 	glEnable(GL_NORMALIZE);
nuclear@3: 
nuclear@3: 	cam_focus_dist(3.0);
nuclear@3: 	cam_clip(0.1, 200.0);
nuclear@3: 	cam_rotate(0, 0);
nuclear@3: 	cam_dolly(2);
nuclear@3: 
nuclear@3: 	msurf = msurf_create();
nuclear@3: 	msurf_eval_func(msurf, eval);
nuclear@3: 	msurf_vertex_func(msurf, vertex);
nuclear@3: 	msurf_threshold(msurf, threshold);
nuclear@3: 	msurf_resolution(msurf, xres, yres, num_slices);
nuclear@3: 	msurf_bounds(msurf, -1, -1, -1, 1, 1, 1);
nuclear@3: 
nuclear@3: 	glClearColor(0.6, 0.6, 0.6, 1.0);
nuclear@3: 
nuclear@3: 	dlist = glGenLists(1);
nuclear@3: 
nuclear@3: 	glutMainLoop();
nuclear@3: 	return 0;
nuclear@3: }
nuclear@3: 
nuclear@3: float eval(float x, float y, float z)
nuclear@3: {
nuclear@3: 	int px, py, slice;
nuclear@3: 	struct img_pixmap *img;
nuclear@3: 
nuclear@3: 	px = round((x * 0.5 + 0.5) * xres);
nuclear@3: 	py = round((y * 0.5 + 0.5) * yres);
nuclear@3: 	slice = round((z * 0.5 + 0.5) * num_slices);
nuclear@3: 
nuclear@3: 	if(px < 0) px = 0;
nuclear@3: 	if(px >= xres) px = xres - 1;
nuclear@3: 
nuclear@3: 	if(py < 0) py = 0;
nuclear@3: 	if(py >= yres) py = yres - 1;
nuclear@3: 
nuclear@3: 	if(slice < 0) slice = 0;
nuclear@3: 	if(slice >= num_slices) slice = num_slices - 1;
nuclear@3: 
nuclear@3: 	img = volume + slice;
nuclear@3: 	return *((unsigned char*)img->pixels + py * img->width + px) / 255.0;
nuclear@3: }
nuclear@3: 
nuclear@3: void vertex(float x, float y, float z)
nuclear@3: {
nuclear@3: 	float dx = 1.0 / xres;
nuclear@3: 	float dy = 1.0 / yres;
nuclear@3: 	float dz = 1.0 / num_slices;
nuclear@3: 	float dfdx = eval(x - dx, y, z) - eval(x + dx, y, z);
nuclear@3: 	float dfdy = eval(x, y - dy, z) - eval(x, y + dy, z);
nuclear@3: 	float dfdz = eval(x, y, z - dz) - eval(x, y, z + dz);
nuclear@3: 
nuclear@3: 	glNormal3f(dfdx, dfdy, dfdz);
nuclear@3: 	glVertex3f(x, y, z);
nuclear@3: }
nuclear@3: 
nuclear@3: void render(void)
nuclear@3: {
nuclear@3: 	float kd[] = {0.87, 0.82, 0.74, 1.0};
nuclear@3: 	float ks[] = {0.9, 0.9, 0.9, 1.0};
nuclear@3: 
nuclear@3: 	glMaterialfv(GL_FRONT_AND_BACK, GL_AMBIENT_AND_DIFFUSE, kd);
nuclear@3: 	glMaterialfv(GL_FRONT_AND_BACK, GL_SPECULAR, ks);
nuclear@3: 	glMaterialf(GL_FRONT_AND_BACK, GL_SHININESS, 60.0);
nuclear@3: 
nuclear@3: #ifndef NO_SHADERS
nuclear@3: 	bind_program(sdr);
nuclear@3: #endif
nuclear@3: 
nuclear@3: 	glMatrixMode(GL_MODELVIEW);
nuclear@3: 	glPushMatrix();
nuclear@11: 	glScalef(1.0, yscale, 1.0);
nuclear@3: 	glRotatef(90, 1, 0, 0);
nuclear@3: 
nuclear@3: 	if(need_update) {
nuclear@3: 		glNewList(dlist, GL_COMPILE);
nuclear@3: 		glBegin(GL_TRIANGLES);
nuclear@3: 		printf("generating mesh... ");
nuclear@3: 		fflush(stdout);
nuclear@3: 		msurf_polygonize(msurf);
nuclear@3: 		glEnd();
nuclear@3: 		glEndList();
nuclear@3: 		need_update = 0;
nuclear@3: 		printf("done\n");
nuclear@3: 	}
nuclear@3: 	glCallList(dlist);
nuclear@3: 
nuclear@3: 	glPopMatrix();
nuclear@3: 
nuclear@3: 	assert(glGetError() == GL_NO_ERROR);
nuclear@3: }
nuclear@3: 
nuclear@3: void disp(void)
nuclear@3: {
nuclear@3: 	if(stereo) {
nuclear@3: 		glDrawBuffer(GL_BACK_LEFT);
nuclear@3: 	}
nuclear@3: 
nuclear@3: 	glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
nuclear@3: 
nuclear@3: 	glMatrixMode(GL_PROJECTION);
nuclear@3: 	glLoadIdentity();
nuclear@3: 	cam_stereo_proj_matrix(stereo ? CAM_LEFT : CAM_CENTER);
nuclear@3: 
nuclear@3: 	glMatrixMode(GL_MODELVIEW);
nuclear@3: 	glLoadIdentity();
nuclear@3: 	cam_stereo_view_matrix(stereo ? CAM_LEFT : CAM_CENTER);
nuclear@3: 
nuclear@3: 	render();
nuclear@3: 
nuclear@3: 	if(stereo) {
nuclear@3: 		glDrawBuffer(GL_BACK_RIGHT);
nuclear@3: 		glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
nuclear@3: 
nuclear@3: 		glMatrixMode(GL_PROJECTION);
nuclear@3: 		glLoadIdentity();
nuclear@3: 		cam_stereo_proj_matrix(CAM_RIGHT);
nuclear@3: 
nuclear@3: 		glMatrixMode(GL_MODELVIEW);
nuclear@3: 		glLoadIdentity();
nuclear@3: 		cam_stereo_view_matrix(CAM_RIGHT);
nuclear@3: 
nuclear@3: 		render();
nuclear@3: 	}
nuclear@3: 	glutSwapBuffers();
nuclear@3: }
nuclear@3: 
nuclear@3: void reshape(int x, int y)
nuclear@3: {
nuclear@3: 	glViewport(0, 0, x, y);
nuclear@3: 	cam_aspect((float)x / (float)y);
nuclear@3: }
nuclear@3: 
nuclear@11: int mode_scale;
nuclear@11: 
nuclear@3: void keyb(unsigned char key, int x, int y)
nuclear@3: {
nuclear@3: 	static int wire;
nuclear@3: 
nuclear@3: 	switch(key) {
nuclear@3: 	case 27:
nuclear@3: 		exit(0);
nuclear@3: 
nuclear@3: 	case 'f':
nuclear@3: 		fullscreen = !fullscreen;
nuclear@3: 		if(fullscreen) {
nuclear@3: 			glutFullScreen();
nuclear@3: 		} else {
nuclear@3: 			glutReshapeWindow(orig_xsz, orig_ysz);
nuclear@3: 		}
nuclear@3: 		break;
nuclear@3: 
nuclear@3: 	case 's':
nuclear@3: 		stereo = !stereo;
nuclear@3: 		glutPostRedisplay();
nuclear@3: 		break;
nuclear@3: 
nuclear@3: 	case 'w':
nuclear@3: 		wire = !wire;
nuclear@3: 		glPolygonMode(GL_FRONT_AND_BACK, wire ? GL_LINE : GL_FILL);
nuclear@3: 		glutPostRedisplay();
nuclear@3: 		break;
nuclear@3: 
nuclear@3: 	case '=':
nuclear@3: 		threshold += 0.05;
nuclear@3: 		msurf_threshold(msurf, threshold);
nuclear@3: 		printf("threshold: %f\n", threshold);
nuclear@3: 		glutPostRedisplay();
nuclear@3: 		need_update = 1;
nuclear@3: 		break;
nuclear@3: 
nuclear@3: 	case '-':
nuclear@3: 		threshold -= 0.05;
nuclear@3: 		msurf_threshold(msurf, threshold);
nuclear@3: 		printf("threshold: %f\n", threshold);
nuclear@3: 		glutPostRedisplay();
nuclear@3: 		need_update = 1;
nuclear@3: 		break;
nuclear@3: 
nuclear@11: 	case 'y':
nuclear@11: 		mode_scale = 1;
nuclear@11: 		break;
nuclear@11: 
nuclear@3: 	default:
nuclear@3: 		break;
nuclear@3: 	}
nuclear@3: }
nuclear@3: 
nuclear@11: void keyb_up(unsigned char key, int x, int y)
nuclear@11: {
nuclear@11: 	switch(key) {
nuclear@11: 	case 'y':
nuclear@11: 		mode_scale = 0;
nuclear@11: 		break;
nuclear@11: 	}
nuclear@11: }
nuclear@11: 
nuclear@3: int bnstate[32];
nuclear@3: int prev_x, prev_y;
nuclear@3: 
nuclear@3: void mouse(int bn, int state, int x, int y)
nuclear@3: {
nuclear@3: 	bnstate[bn] = state == GLUT_DOWN;
nuclear@3: 	prev_x = x;
nuclear@3: 	prev_y = y;
nuclear@3: }
nuclear@3: 
nuclear@3: void motion(int x, int y)
nuclear@3: {
nuclear@3: 	int dx, dy;
nuclear@3: 
nuclear@3: 	dx = x - prev_x;
nuclear@3: 	dy = y - prev_y;
nuclear@3: 	prev_x = x;
nuclear@3: 	prev_y = y;
nuclear@3: 
nuclear@11: 	if(mode_scale) {
nuclear@11: 		yscale += dy * 0.001;
nuclear@11: 	} else {
nuclear@11: 		if(bnstate[GLUT_LEFT_BUTTON]) {
nuclear@11: 			cam_inp_rotate(dx, dy);
nuclear@11: 		}
nuclear@11: 		if(bnstate[GLUT_RIGHT_BUTTON]) {
nuclear@11: 			cam_inp_zoom(dy);
nuclear@11: 		}
nuclear@3: 	}
nuclear@3: 	glutPostRedisplay();
nuclear@3: }
nuclear@3: 
nuclear@3: struct list_node {
nuclear@3: 	struct img_pixmap img;
nuclear@3: 	struct list_node *next;
nuclear@3: };
nuclear@3: 
nuclear@3: int parse_args(int argc, char **argv)
nuclear@3: {
nuclear@3: 	int i;
nuclear@6: 	char *endp;
nuclear@3: 	struct list_node *head = 0, *tail = 0;
nuclear@3: 
nuclear@3: 	for(i=1; i<argc; i++) {
nuclear@3: 		if(argv[i][0] == '-' && argv[i][2] == 0) {
nuclear@3: 			switch(argv[i][1]) {
nuclear@3: 			case 'f':
nuclear@3: 				fullscreen = !fullscreen;
nuclear@3: 				break;
nuclear@3: 
nuclear@3: 			case 's':
nuclear@3: 				stereo = !stereo;
nuclear@3: 				break;
nuclear@3: 
nuclear@6: 			case 't':
nuclear@6: 				threshold = strtod(argv[++i], &endp);
nuclear@6: 				if(endp == argv[i]) {
nuclear@6: 					fprintf(stderr, "-t must be followed by a number\n");
nuclear@6: 					return -1;
nuclear@6: 				}
nuclear@6: 				break;
nuclear@6: 
nuclear@3: 			case 'h':
nuclear@3: 				printf("usage: %s [opt]\n", argv[0]);
nuclear@3: 				printf("options:\n");
nuclear@3: 				printf("  -f    start in fullscreen\n");
nuclear@3: 				printf("  -s    enable stereoscopic rendering\n");
nuclear@3: 				printf("  -h    print usage and exit\n");
nuclear@3: 				exit(0);
nuclear@3: 
nuclear@3: 			default:
nuclear@3: 				fprintf(stderr, "unrecognized option: %s\n", argv[i]);
nuclear@3: 				return -1;
nuclear@3: 			}
nuclear@3: 		} else {
nuclear@3: 			struct list_node *slice;
nuclear@3: 
nuclear@3: 			if(!(slice = malloc(sizeof *slice))) {
nuclear@3: 				fprintf(stderr, "failed to allocate volume slice: %d\n", num_slices);
nuclear@3: 				return -1;
nuclear@3: 			}
nuclear@3: 			slice->next = 0;
nuclear@3: 
nuclear@3: 			img_init(&slice->img);
nuclear@3: 			if(img_load(&slice->img, argv[i]) == -1) {
nuclear@3: 				fprintf(stderr, "failed to load volume slice %d: %s\n", num_slices, argv[i]);
nuclear@3: 				free(slice);
nuclear@3: 				return -1;
nuclear@3: 			}
nuclear@3: 			img_convert(&slice->img, IMG_FMT_GREY8);
nuclear@3: 
nuclear@3: 			if(num_slices > 0 && (xres != slice->img.width || yres != slice->img.height)) {
nuclear@3: 				fprintf(stderr, "error: slice %d (%s) is %dx%d, up to now we had %dx%d images\n", num_slices, argv[i],
nuclear@3: 						slice->img.width, slice->img.height, xres, yres);
nuclear@3: 				img_destroy(&slice->img);
nuclear@3: 				free(slice);
nuclear@3: 				return -1;
nuclear@3: 			}
nuclear@3: 			xres = slice->img.width;
nuclear@3: 			yres = slice->img.height;
nuclear@3: 
nuclear@3: 			if(head) {
nuclear@3: 				tail->next = slice;
nuclear@3: 				tail = slice;
nuclear@3: 			} else {
nuclear@3: 				head = tail = slice;
nuclear@3: 			}
nuclear@3: 			printf("loaded volume slice %d: %s\n", num_slices++, argv[i]);
nuclear@3: 		}
nuclear@3: 	}
nuclear@3: 
nuclear@3: 	if(!head) {
nuclear@3: 		fprintf(stderr, "you must specify a list of images for the volume data slices\n");
nuclear@3: 		return -1;
nuclear@3: 	}
nuclear@3: 
nuclear@3: 	if(!(volume = malloc(num_slices * sizeof *volume))) {
nuclear@3: 		fprintf(stderr, "failed to allocate volume data (%d slices)\n", num_slices);
nuclear@3: 		return -1;
nuclear@3: 	}
nuclear@3: 
nuclear@3: 	for(i=0; i<num_slices; i++) {
nuclear@3: 		void *tmp;
nuclear@3: 
nuclear@3: 		assert(head);
nuclear@3: 		volume[i] = head->img;
nuclear@3: 
nuclear@3: 		tmp = head;
nuclear@3: 		head = head->next;
nuclear@3: 		free(tmp);
nuclear@3: 	}
nuclear@3: 
nuclear@3: 	return 0;
nuclear@3: }