qvolray
view src/volume.c @ 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 |
line source
1 #include <stdio.h>
2 #include <stdlib.h>
3 #include <string.h>
4 #include <ctype.h>
6 #ifndef __APPLE__
7 #include <GL/gl.h>
8 #else
9 #include <OpenGL/gl.h>
10 #endif
12 #include <imago2.h>
13 #include "volume.h"
15 struct slice {
16 char *name;
17 struct slice *next;
18 };
20 static void calc_gradients(float *voxels, int xsz, int ysz, int zsz);
21 static struct slice *read_voldesc(const char *fname, struct volume *vol);
22 static char *trim(char *str);
24 struct volume *load_volume(const char *fname)
25 {
26 int i;
27 struct slice *slist;
28 struct volume *vol = 0;
29 float *voxels = 0, *vptr;
30 struct img_pixmap img;
32 if(!(vol = malloc(sizeof *vol))) {
33 perror("failed to allocate volume");
34 return 0;
35 }
36 memset(vol, 0, sizeof *vol);
37 vol->zaspect = 1;
39 if(!(slist = read_voldesc(fname, vol))) {
40 goto err;
41 }
43 /* load the first image to determine slice dimensions */
44 img_init(&img);
45 if(img_load(&img, slist->name) == -1) {
46 fprintf(stderr, "failed to load volume slice: %s\n", slist->name);
47 goto err;
48 }
50 vol->sz[0] = img.width;
51 vol->sz[1] = img.height;
52 printf("volume dimensions: %dx%dx%d\n", img.width, img.height, vol->sz[2]);
54 /* allocate the whole volume at once */
55 if(!(voxels = malloc(vol->sz[0] * vol->sz[1] * vol->sz[2] * 4 * sizeof *voxels))) {
56 img_destroy(&img);
57 goto err;
58 }
59 vptr = voxels;
60 img_destroy(&img);
62 /* put the volume data into the alpha component */
63 for(i=0; i<vol->sz[2]; i++) {
64 int x, y, xsz, ysz;
65 float *pixels, *src;
66 struct slice *slice = slist;
67 slist = slist->next;
69 if(!(pixels = img_load_pixels(slice->name, &xsz, &ysz, IMG_FMT_GREYF))) {
70 fprintf(stderr, "failed to load volume slice: %s\n", slice->name);
71 goto err;
72 }
74 if(xsz != vol->sz[0] || ysz != vol->sz[1]) {
75 fprintf(stderr, "inconsistent dimensions for slice %s: %dx%d (%dx%d expected)\n",
76 slice->name, xsz, ysz, vol->sz[0], vol->sz[1]);
77 goto err;
78 }
79 free(slice->name);
80 free(slice);
82 src = pixels;
83 for(y=0; y<ysz; y++) {
84 for(x=0; x<xsz; x++) {
85 vptr[0] = vptr[1] = vptr[2] = 0.0f;
86 vptr[3] = *src++;
87 vptr += 4;
88 }
89 }
90 img_free_pixels(pixels);
91 }
93 /* calculate gradients */
94 /*calc_gradients(voxels, vol->sz[0], vol->sz[1], vol->sz[2]);*/
96 /* create the volume texture */
97 glGenTextures(1, &vol->tex_vol);
98 glBindTexture(GL_TEXTURE_3D, vol->tex_vol);
99 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MIN_FILTER, GL_LINEAR);
100 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MAG_FILTER, GL_LINEAR);
101 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE);
102 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE);
103 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_R, GL_CLAMP_TO_EDGE);
104 glTexImage3D(GL_TEXTURE_3D, 0, GL_RGBA32F, vol->sz[0], vol->sz[1], vol->sz[2], 0,
105 GL_RGBA, GL_FLOAT, voxels);
107 /* ... and destroy our copy */
108 free(voxels);
109 return vol;
111 err:
112 while(slist) {
113 struct slice *tmp = slist;
114 slist = slist->next;
115 free(tmp->name);
116 free(tmp);
117 }
118 free_volume(vol);
119 free(voxels);
120 return 0;
121 }
123 void free_volume(struct volume *vol)
124 {
125 if(vol) {
126 if(vol->tex_vol) {
127 glDeleteTextures(1, &vol->tex_vol);
128 }
129 free(vol);
130 }
131 }
133 struct vec4 { float x, y, z, w; };
135 static void calc_gradients(float *voxels, int xsz, int ysz, int zsz)
136 {
137 int x, y, z, slice_pixels = xsz * ysz;
138 struct vec4 *vptr = (struct vec4*)voxels;
140 for(z=0; z<zsz; z++) {
141 for(y=0; y<ysz; y++) {
142 for(x=0; x<xsz; x++) {
143 float x0, x1, y0, y1, z0, z1;
145 x0 = x > 0 ? (vptr - 1)->w : vptr->w;
146 y0 = y > 0 ? (vptr - xsz)->w : vptr->w;
147 z0 = z > 0 ? (vptr - slice_pixels)->w : vptr->w;
148 x1 = x < xsz - 1 ? (vptr + 1)->w : vptr->w;
149 y1 = y < ysz - 1 ? (vptr + xsz)->w : vptr->w;
150 z1 = z < zsz - 1 ? (vptr + slice_pixels)->w : vptr->w;
152 vptr->x = x1 - x0;
153 vptr->y = y1 - y0;
154 vptr->z = z1 - z0;
156 vptr++;
157 }
158 }
159 }
160 }
162 static struct slice *read_voldesc(const char *fname, struct volume *vol)
163 {
164 FILE *fp = 0;
165 char buf[512], *slash, *prefix = "";
166 int mode_slices = 0;
167 struct slice *head = 0, *tail;
169 if(!(fp = fopen(fname, "r"))) {
170 fprintf(stderr, "failed to open volume descriptor: %s\n", fname);
171 return 0;
172 }
173 fgets(buf, sizeof buf, fp);
174 if(strstr(buf, "VOLDESC") != buf) {
175 fprintf(stderr, "invalid file format while trying to read volume descriptor: %s\n", fname);
176 goto err;
177 }
179 if((slash = strrchr(fname, '/'))) {
180 int len = slash - fname + 1;
181 prefix = alloca(len + 1);
182 memcpy(prefix, fname, len);
183 prefix[len] = 0;
184 }
186 while(fgets(buf, sizeof buf, fp)) {
187 char *line = trim(buf);
189 if(!*line || *line == '#')
190 continue;
192 if(mode_slices) {
193 /* we're in the part that contains filenames... append to the list */
194 struct slice *node = malloc(sizeof *node);
195 if(!node || !(node->name = malloc(strlen(prefix) + strlen(line) + 1))) {
196 perror("failed to allocate list node");
197 free(node);
198 goto err;
199 }
200 sprintf(node->name, "%s%s", prefix, line);
201 node->next = 0;
203 if(head) {
204 tail->next = node;
205 tail = node;
206 } else {
207 head = tail = node;
208 }
209 vol->sz[2]++;
210 } else {
211 /* TODO we're in the options part... parse */
213 if(strcmp(line, "SLICES") == 0) {
214 mode_slices = 1;
215 }
216 }
217 }
218 fclose(fp);
220 if(!vol->sz[2]) {
221 fprintf(stderr, "file: %s contains no slices\n", fname);
222 goto err;
223 }
225 return head;
227 err:
228 while(head) {
229 struct slice *tmp = head;
230 head = head->next;
231 free(tmp->name);
232 free(tmp);
233 }
234 if(fp) {
235 fclose(fp);
236 }
237 return 0;
238 }
240 static char *trim(char *str)
241 {
242 char *tmp = str + strlen(str) - 1;
244 while(isspace(*tmp))
245 tmp--;
246 tmp[1] = 0;
248 tmp = str;
249 while(isspace(*tmp))
250 tmp++;
251 return tmp;
252 }