qvolray

view src/volume.cc @ 11:8990b5d2c7fe

find changesets by author, revision, files, or words in the commit message
moving to qt
author John Tsiombikas <nuclear@member.fsf.org>
date Mon, 09 Apr 2012 23:42:57 +0300
parents src/volume.c@9d2396738b60
children 17d9dc2edc91
line source
1 #include <stdio.h>
2 #include <stdlib.h>
3 #include <string.h>
4 #include <ctype.h>
5 #include <list>
6 #include <string>
7
8 #ifndef __APPLE__
9 #include <GL/gl.h>
10 #else
11 #include <OpenGL/gl.h>
12 #endif
13
14 #include <imago2.h>
15 #include "volume.h"
16
17 static void calc_gradients(float *voxels, int xsz, int ysz, int zsz);
18 static struct slice *read_voldesc(const char *fname, struct volume *vol);
19 static char *trim(char *str);
20
21 Volume::Volume()
22 {
23 sz[0] = sz[1] = sz[2] = 0;
24 zaspect = 1.0;
25 tex = 0;
26 }
27
28 Volume::~Volume()
29 {
30 if(tex) {
31 glDeleteTextures(1, &tex);
32 }
33 }
34
35 bool Volume::load(const char *fname)
36 {
37 std::list<std::string> slist;
38 float *voxels = 0, *vptr;
39 struct img_pixmap img;
40
41 if(!(read_voldesc(fname, &slist))) {
42 return false;
43 }
44
45 /* load the first image to determine slice dimensions */
46 img_init(&img);
47 if(img_load(&img, slist.begin().name.c_str()) == -1) {
48 fprintf(stderr, "failed to load volume slice: %s\n", slist->name);
49 return false;
50 }
51
52 sz[0] = img.width;
53 sz[1] = img.height;
54 printf("volume dimensions: %dx%dx%d\n", img.width, img.height, vol->sz[2]);
55
56 /* allocate the whole volume at once */
57 voxels = new float[sz[0] * sz[1] * sz[2] * 4];
58 vptr = voxels;
59 img_destroy(&img);
60
61 /* put the volume data into the alpha component */
62 for(auto slice : slist) {
63 int x, y, xsz, ysz;
64 float *pixels, *src;
65
66 if(!(pixels = img_load_pixels(slice.name.c_str(), &xsz, &ysz, IMG_FMT_GREYF))) {
67 fprintf(stderr, "failed to load volume slice: %s\n", slice.name.c_str());
68 delete [] voxels;
69 return false;
70 }
71
72 if(xsz != sz[0] || ysz != sz[1]) {
73 fprintf(stderr, "inconsistent dimensions for slice %s: %dx%d (%dx%d expected)\n",
74 slice.name.c_str(), xsz, ysz, sz[0], sz[1]);
75 delete [] voxels;
76 return false;
77 }
78
79 src = pixels;
80 for(y=0; y<ysz; y++) {
81 for(x=0; x<xsz; x++) {
82 vptr[0] = vptr[1] = vptr[2] = 0.0f;
83 vptr[3] = *src++;
84 vptr += 4;
85 }
86 }
87 img_free_pixels(pixels);
88 }
89
90 /* calculate gradients */
91 calc_gradients(voxels, sz[0], sz[1], sz[2]);
92
93 /* create the volume texture */
94 glGenTextures(1, &tex);
95 glBindTexture(GL_TEXTURE_3D, tex);
96 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MIN_FILTER, GL_LINEAR);
97 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MAG_FILTER, GL_LINEAR);
98 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE);
99 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE);
100 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_R, GL_CLAMP_TO_EDGE);
101 glTexImage3D(GL_TEXTURE_3D, 0, GL_RGBA32F_ARB, sz[0], sz[1], sz[2], 0, GL_RGBA, GL_FLOAT, voxels);
102
103 /* ... and destroy our copy */
104 delete [] voxels;
105 return true;
106 }
107
108 unsigned int Volume::get_texture() const
109 {
110 return tex;
111 }
112
113 void Volume::bind(int tex_unit) const
114 {
115 glActiveTexture(GL_TEXTURE0 + tex_unit);
116 glBindTexture(GL_TEXTURE_3D, tex);
117 }
118
119 struct vec4 { float x, y, z, w; };
120
121 static void calc_gradients(float *voxels, int xsz, int ysz, int zsz)
122 {
123 int x, y, z, slice_pixels = xsz * ysz;
124 vec4 *vptr = (vec4*)voxels;
125
126 for(z=0; z<zsz; z++) {
127 for(y=0; y<ysz; y++) {
128 for(x=0; x<xsz; x++) {
129 float x0, x1, y0, y1, z0, z1;
130
131 x0 = x > 0 ? (vptr - 1)->w : vptr->w;
132 y0 = y > 0 ? (vptr - xsz)->w : vptr->w;
133 z0 = z > 0 ? (vptr - slice_pixels)->w : vptr->w;
134 x1 = x < xsz - 1 ? (vptr + 1)->w : vptr->w;
135 y1 = y < ysz - 1 ? (vptr + xsz)->w : vptr->w;
136 z1 = z < zsz - 1 ? (vptr + slice_pixels)->w : vptr->w;
137
138 vptr->x = x1 - x0;
139 vptr->y = y1 - y0;
140 vptr->z = z1 - z0;
141
142 vptr++;
143 }
144 }
145 }
146 }
147
148 bool read_voldesc(const char *fname, std::list<std::string> *slist) const
149 {
150 FILE *fp = 0;
151 char buf[512], *slash, *prefix = "";
152 int mode_slices = 0;
153
154 slist->clear();
155
156 if(!(fp = fopen(fname, "r"))) {
157 fprintf(stderr, "failed to open volume descriptor: %s\n", fname);
158 return false;
159 }
160 fgets(buf, sizeof buf, fp);
161 if(strstr(buf, "VOLDESC") != buf) {
162 fprintf(stderr, "invalid file format while trying to read volume descriptor: %s\n", fname);
163 fclose(fp);
164 return false;
165 }
166
167 if((slash = strrchr(fname, '/'))) {
168 int len = slash - fname + 1;
169 prefix = alloca(len + 1);
170 memcpy(prefix, fname, len);
171 prefix[len] = 0;
172 }
173
174 while(fgets(buf, sizeof buf, fp)) {
175 char *line = trim(buf);
176
177 if(!*line || *line == '#')
178 continue;
179
180 if(mode_slices) {
181 /* we're in the part that contains filenames... append to the list */
182 slist->push_back(std::string(prefix) + std::string(line));
183 sz[2]++;
184 } else {
185 /* TODO we're in the options part... parse */
186 if(strcmp(line, "SLICES") == 0) {
187 mode_slices = 1;
188 }
189 }
190 }
191 fclose(fp);
192
193 if(slist->empty()) {
194 fprintf(stderr, "file: %s contains no slices\n", fname);
195 return false;
196 }
197 return true;
198 }
199
200 static char *trim(char *str)
201 {
202 char *tmp = str + strlen(str) - 1;
203
204 while(isspace(*tmp))
205 tmp--;
206 tmp[1] = 0;
207
208 tmp = str;
209 while(isspace(*tmp))
210 tmp++;
211 return tmp;
212 }