qvolray: 3e53a16d4667 src/volray.c

qvolray

view src/volray.c @ 4:3e53a16d4667

foo

author	John Tsiombikas <nuclear@member.fsf.org>
date	Mon, 02 Apr 2012 17:59:46 +0300
parents	6f275934717b
children	0c3874aa717a

line source

1 #include <stdio.h>

2 #include <stdlib.h>

3 #include <assert.h>

5 #include <GL/glew.h>

6 #ifndef __APPLE__

7 #include <GL/glut.h>

8 #else

9 #include <GLUT/glut.h>

10 #endif

12 #include <vmath/vmath.h>

13 #include <imago2.h>

14 #include "sdr.h"

16 #define XFER_MAP_SZ 512

18 struct slice_file {

19 char *name;

20 struct slice_file *next;

21 };

23 int init(void);

24 void disp(void);

25 void render_volume(void);

26 void draw_xfer_func(void);

27 void reshape(int x, int y);

28 void keyb(unsigned char key, int x, int y);

29 void keyb_up(unsigned char key, int x, int y);

30 void mouse(int bn, int state, int x, int y);

31 void motion(int x, int y);

32 int parse_args(int argc, char **argv);

34 static void create_ray_texture(int xsz, int ysz, float vfov, vec2_t *tex_scale);

35 static vec3_t get_primary_ray_dir(int x, int y, int w, int h, float vfov_deg);

36 static int round_pow2(int x);

37 static void create_transfer_map(float mean, float sdev);

39 float cam_theta = 0, cam_phi = 0, cam_dist = 4.0;

40 float cam_x, cam_y, cam_z;

42 vec2_t tex_scale;

43 struct slice_file *flist;

44 int nslices;

45 unsigned int sdr, vol_tex, ray_tex;

46 int win_xsz, win_ysz;

47 int raytex_needs_recalc = 1;

49 unsigned int xfer_tex;

50 float xfer_mean = 0.5, xfer_sdev = 1.0;

51 int xfertex_needs_recalc = 1;

53 static int uimode_xfer;

55 int main(int argc, char **argv)

56 {

57 glutInit(&argc, argv);

58 glutInitWindowSize(1280, 720);

59 glutInitDisplayMode(GLUT_RGB | GLUT_DOUBLE);

61 if(parse_args(argc, argv) == -1) {

62 return 1;

63 }

65 glutCreateWindow("Volume Raytracer");

67 glutDisplayFunc(disp);

68 glutReshapeFunc(reshape);

69 glutKeyboardFunc(keyb);

70 glutKeyboardUpFunc(keyb_up);

71 glutMouseFunc(mouse);

72 glutMotionFunc(motion);

74 glewInit();

76 if(init() == -1) {

77 return 1;

78 }

80 glutMainLoop();

81 return 0;

82 }

84 int init(void)

85 {

86 int i, vol_xsz, vol_ysz;

88 if(!(sdr = create_program_load("volray.v.glsl", "volray.p.glsl"))) {

89 return -1;

90 }

91 set_uniform_int(sdr, "volume", 0);

92 set_uniform_int(sdr, "ray_tex", 1);

93 set_uniform_int(sdr, "xfer_tex", 2);

95 glGenTextures(1, &vol_tex);

96 glBindTexture(GL_TEXTURE_3D, vol_tex);

97 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MIN_FILTER, GL_LINEAR);

98 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MAG_FILTER, GL_LINEAR);

99 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE);

100 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE);

101 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_R, GL_CLAMP_TO_EDGE);

102

103 for(i=0; i<nslices; i++) {

104 int xsz, ysz;

105 void *pix;

106 struct slice_file *sfile = flist;

107 flist = flist->next;

108

109 if(!(pix = img_load_pixels(sfile->name, &xsz, &ysz, IMG_FMT_RGBA32))) {

110 fprintf(stderr, "failed to load image: %s\n", sfile->name);

111 return -1;

112 }

113

114 if(i == 0) {

115 /* allocate storage for the texture */

116 glTexImage3D(GL_TEXTURE_3D, 0, GL_RGBA, xsz, ysz, nslices, 0, GL_RGBA, GL_UNSIGNED_BYTE, 0);

117

118 vol_xsz = xsz;

119 vol_ysz = ysz;

120

121 } else {

122 if(xsz != vol_xsz || ysz != vol_ysz) {

123 fprintf(stderr, "%s: inconsistent slice size: %dx%d. expected: %dx%d\n",

124 sfile->name, xsz, ysz, vol_xsz, vol_ysz);

125 img_free_pixels(pix);

126 return -1;

127 }

128 }

129 free(sfile);

130

131 glTexSubImage3D(GL_TEXTURE_3D, 0, 0, 0, i, xsz, ysz, 1, GL_RGBA, GL_UNSIGNED_BYTE, pix);

132 img_free_pixels(pix);

133 }

134 return 0;

135 }

136

137 void disp(void)

138 {

139 /* recalculate primary ray texture if needed */

140 if(raytex_needs_recalc) {

141 create_ray_texture(win_xsz, win_ysz, 50.0, &tex_scale);

142 }

143 /* recalculate transfer function texture if needed */

144 if(xfertex_needs_recalc) {

145 create_transfer_map(xfer_mean, xfer_sdev);

146 }

147

148 render_volume();

149 draw_xfer_func();

150

151 glutSwapBuffers();

152 assert(glGetError() == GL_NO_ERROR);

153 }

154

155 void render_volume(void)

156 {

157 /* set the camera transformation */

158 glMatrixMode(GL_MODELVIEW);

159 glPushMatrix();

160 glLoadIdentity();

161 glRotatef(-90, 1, 0, 0);

162 glTranslatef(cam_x, cam_y, -cam_z);

163 glRotatef(cam_theta, 0, 1, 0);

164 glRotatef(cam_phi, 1, 0, 0);

165 glTranslatef(0, 0, -cam_dist);

166

167 /* setup the texture matrix to map the useful part of the ray texture to [0,1] */

168 glMatrixMode(GL_TEXTURE);

169 glPushMatrix();

170 glLoadIdentity();

171 glScalef(tex_scale.x, tex_scale.y, 1.0);

172

173 /* tex unit0: volume data 3D texture */

174 glActiveTexture(GL_TEXTURE0);

175 glBindTexture(GL_TEXTURE_3D, vol_tex);

176 glEnable(GL_TEXTURE_3D);

177

178 /* tex unit1: primary rays in view space */

179 glActiveTexture(GL_TEXTURE1);

180 glBindTexture(GL_TEXTURE_2D, ray_tex);

181 glEnable(GL_TEXTURE_2D);

182

183 /* tex unit2: transfer function (1d) */

184 glActiveTexture(GL_TEXTURE2);

185 glBindTexture(GL_TEXTURE_1D, xfer_tex);

186 glEnable(GL_TEXTURE_1D);

187

188 bind_program(sdr);

189 glBegin(GL_QUADS);

190 glColor3f(1, 1, 1);

191 glTexCoord2f(0, 1); glVertex2f(-1, -1);

192 glTexCoord2f(1, 1); glVertex2f(1, -1);

193 glTexCoord2f(1, 0); glVertex2f(1, 1);

194 glTexCoord2f(0, 0); glVertex2f(-1, 1);

195 glEnd();

196 bind_program(0);

197

198 glActiveTexture(GL_TEXTURE2);

199 glDisable(GL_TEXTURE_1D);

200 glActiveTexture(GL_TEXTURE1);

201 glDisable(GL_TEXTURE_2D);

202 glActiveTexture(GL_TEXTURE0);

203 glDisable(GL_TEXTURE_3D);

204

205 glMatrixMode(GL_TEXTURE);

206 glPopMatrix();

207 glMatrixMode(GL_MODELVIEW);

208 glPopMatrix();

209 }

210

211 void draw_xfer_func(void)

212 {

213 glMatrixMode(GL_MODELVIEW);

214 glPushMatrix();

215 glTranslatef(-0.9, -0.9, 0);

216 glScalef(0.5, 0.1, 1);

217

218 glBindTexture(GL_TEXTURE_1D, xfer_tex);

219 glEnable(GL_TEXTURE_1D);

220

221 glBegin(GL_QUADS);

222 glColor3f(1, 1, 1);

223 glTexCoord1f(1);

224 glVertex2f(1, 0);

225 glVertex2f(1, 1);

226 glTexCoord1f(0);

227 glVertex2f(0, 1);

228 glVertex2f(0, 0);

229 glEnd();

230

231 glDisable(GL_TEXTURE_1D);

232

233 glLineWidth(2.0);

234 glBegin(GL_LINE_LOOP);

235 glColor3f(uimode_xfer ? 1 : 0, 0, uimode_xfer ? 0 : 1);

236 glVertex2f(0, 0);

237 glVertex2f(1, 0);

238 glVertex2f(1, 1);

239 glVertex2f(0, 1);

240 glEnd();

241

242 glPopMatrix();

243 }

244

245 void reshape(int x, int y)

246 {

247 glViewport(0, 0, x, y);

248

249 if(x != win_xsz || y != win_ysz) {

250 raytex_needs_recalc = 1;

251 win_xsz = x;

252 win_ysz = y;

253 }

254 }

255

256 void keyb(unsigned char key, int x, int y)

257 {

258 switch(key) {

259 case 27:

260 exit(0);

261

262 case 'x':

263 uimode_xfer = 1;

264 glutPostRedisplay();

265 break;

266 }

267 }

268

269 void keyb_up(unsigned char key, int x, int y)

270 {

271 switch(key) {

272 case 'x':

273 uimode_xfer = 0;

274 glutPostRedisplay();

275 break;

276 }

277 }

278

279 static int bnstate[32];

280 static int prev_x, prev_y;

281

282 void mouse(int bn, int state, int x, int y)

283 {

284 bnstate[bn - GLUT_LEFT_BUTTON] = state == GLUT_DOWN;

285 prev_x = x;

286 prev_y = y;

287 }

288

289 void motion(int x, int y)

290 {

291 int dx = x - prev_x;

292 int dy = y - prev_y;

293 prev_x = x;

294 prev_y = y;

295

296 if(uimode_xfer) {

297 if(dx || dy) {

298 xfer_mean += dx / (float)win_xsz;

299 xfer_sdev += 0.5 * dy / (float)win_ysz;

300

301 xfer_mean = xfer_mean < 0.0 ? 0.0 : (xfer_mean > 1.0 ? 1.0 : xfer_mean);

302 xfer_sdev = xfer_sdev < 0.0 ? 0.0 : (xfer_sdev > 1.0 ? 1.0 : xfer_sdev);

303

304 xfertex_needs_recalc = 1;

305 glutPostRedisplay();

306 }

307 } else {

308

309 if(bnstate[0]) {

310 cam_theta += dx * 0.5;

311 cam_phi += dy * 0.5;

312

313 if(cam_phi <= -90) cam_phi = -89;

314 if(cam_phi >= 90) cam_phi = 89;

315

316 glutPostRedisplay();

317 }

318

319 if(bnstate[1]) {

320 cam_x += dx * 0.025;

321 cam_y += dy * 0.025;

322 glutPostRedisplay();

323 }

324

325 if(bnstate[2]) {

326 cam_dist += dy * 0.025;

327 if(cam_dist < 0.0) cam_dist = 0.0;

328 glutPostRedisplay();

329 }

330 }

331 }

332

333

334 int parse_args(int argc, char **argv)

335 {

336 int i;

337 struct slice_file *tail;

338 char *endp;

339

340 for(i=1; i<argc; i++) {

341 if(argv[i][0] == '-' && argv[i][2] == 0) {

342 switch(argv[i][1]) {

343 case 'm':

344 xfer_mean = strtod(argv[++i], &endp);

345 if(endp == argv[i]) {

346 fprintf(stderr, "-m must be followed by the transfer function mean\n");

347 return -1;

348 }

349 break;

350

351 case 'd':

352 xfer_sdev = strtod(argv[++i], &endp);

353 if(endp == argv[i]) {

354 fprintf(stderr, "-d must be followed by the transfer function std.deviation\n");

355 return -1;

356 }

357 break;

358

359 default:

360 fprintf(stderr, "unrecognized option: %s\n", argv[i]);

361 return -1;

362 }

363 } else {

364 struct slice_file *sfile;

365

366 if(!(sfile = malloc(sizeof *sfile))) {

367 perror("failed to allocate memory");

368 return -1;

369 }

370 sfile->name = argv[i];

371 sfile->next = 0;

372

373 if(!flist) {

374 flist = tail = sfile;

375 } else {

376 tail->next = sfile;

377 tail = sfile;

378 }

379 nslices++;

380 }

381 }

382

383 if(!nslices) {

384 fprintf(stderr, "pass the slice filenames\n");

385 return -1;

386 }

387 return 0;

388 }

389

390

391 static void create_ray_texture(int xsz, int ysz, float vfov, vec2_t *tex_scale)

392 {

393 int i, j;

394 int cur_tex_xsz, cur_tex_ysz;

395 int tex_xsz = round_pow2(xsz);

396 int tex_ysz = round_pow2(ysz);

397 float *teximg, *dir;

398

399 if(!(teximg = malloc(3 * xsz * ysz * sizeof *teximg))) {

400 return;

401 }

402 dir = teximg;

403

404 for(i=0; i<ysz; i++) {

405 for(j=0; j<xsz; j++) {

406 vec3_t rdir = get_primary_ray_dir(j, i, xsz, ysz, vfov);

407 *dir++ = rdir.x;

408 *dir++ = rdir.y;

409 *dir++ = rdir.z;

410 }

411 }

412

413 if(!ray_tex) {

414 glGenTextures(1, &ray_tex);

415 }

416

417 glGetTexLevelParameteriv(GL_TEXTURE_2D, 0, GL_TEXTURE_WIDTH, &cur_tex_xsz);

418 glGetTexLevelParameteriv(GL_TEXTURE_2D, 0, GL_TEXTURE_HEIGHT, &cur_tex_ysz);

419

420 if(tex_xsz > cur_tex_xsz || tex_ysz > cur_tex_ysz) {

421 glBindTexture(GL_TEXTURE_2D, ray_tex);

422 glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MIN_FILTER, GL_NEAREST);

423 glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_MAG_FILTER, GL_NEAREST);

424 glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_WRAP_S, GL_CLAMP);

425 glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_WRAP_T, GL_CLAMP);

426 glTexImage2D(GL_TEXTURE_2D, 0, GL_RGB32F_ARB, tex_xsz, tex_ysz, 0, GL_RGB, GL_FLOAT, 0);

427 }

428

429 glTexSubImage2D(GL_TEXTURE_2D, 0, 0, 0, xsz, ysz, GL_RGB, GL_FLOAT, teximg);

430 free(teximg);

431

432 if(tex_scale) {

433 tex_scale->x = (float)xsz / (float)tex_xsz;

434 tex_scale->y = (float)ysz / (float)tex_ysz;

435 }

436 raytex_needs_recalc = 0;

437 }

438

439 static vec3_t get_primary_ray_dir(int x, int y, int w, int h, float vfov_deg)

440 {

441 float vfov = M_PI * vfov_deg / 180.0;

442 float aspect = (float)w / (float)h;

443

444 float ysz = 2.0;

445 float xsz = aspect * ysz;

446

447 float px = ((float)x / (float)w) * xsz - xsz / 2.0;

448 float py = 1.0 - ((float)y / (float)h) * ysz;

449 float pz = 1.0 / tan(0.5 * vfov);

450

451 float mag = sqrt(px * px + py * py + pz * pz);

452

453 return v3_cons(px / mag, py / mag, pz / mag);

454 }

455

456 static int round_pow2(int x)

457 {

458 x--;

459 x = (x >> 1) | x;

460 x = (x >> 2) | x;

461 x = (x >> 4) | x;

462 x = (x >> 8) | x;

463 x = (x >> 16) | x;

464 return x + 1;

465 }

466

467 static void create_transfer_map(float mean, float sdev)

468 {

469 static float map[XFER_MAP_SZ];

470 int i;

471

472 if(!xfer_tex) {

473 glGenTextures(1, &xfer_tex);

474 glBindTexture(GL_TEXTURE_1D, xfer_tex);

475 glTexParameteri(GL_TEXTURE_1D, GL_TEXTURE_MIN_FILTER, GL_LINEAR);

476 glTexParameteri(GL_TEXTURE_1D, GL_TEXTURE_MAG_FILTER, GL_LINEAR);

477 glTexParameteri(GL_TEXTURE_1D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE);

478 glTexParameteri(GL_TEXTURE_1D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE);

479 glTexImage1D(GL_TEXTURE_1D, 0, GL_LUMINANCE32F_ARB, XFER_MAP_SZ, 0, GL_LUMINANCE, GL_FLOAT, 0);

480 }

481

482 for(i=0; i<XFER_MAP_SZ; i++) {

483 float x = (float)i / (float)(XFER_MAP_SZ - 1);

484 map[i] = gaussian(x, mean, sdev) - 1.0;

485 }

486 putchar('\n');

487

488 glTexSubImage1D(GL_TEXTURE_1D, 0, 0, XFER_MAP_SZ, GL_LUMINANCE, GL_FLOAT, map);

489 xfertex_needs_recalc = 0;

490 }