volmetrics

view src/volume.cc @ 37:1d935677d3ab

Created new branch non-qt-version
author Eleni Maria Stea <elene.mst@gmail.com>
date Mon, 17 Apr 2017 13:26:24 +0000
parents eb87d4a12bd3
children
line source
1 #include "opengl.h"
3 #include <assert.h>
5 #include <ctype.h>
6 #include <errno.h>
7 #include <float.h>
8 #include <stdio.h>
9 #include <string.h>
11 #include <string>
12 #include <vector>
14 #include "mesh.h"
15 #include "volume.h"
17 static char *strip_whitespaces(char *buf);
19 Volume::Volume()
20 {
21 width = height = 0;
22 zaspect = 1;
23 x_rot = 0;
24 vol_tex_valid = false;
26 memset(histogram, 0, HIST_SIZE * sizeof(int));
27 max_histogram_value = 0;
29 glGenTextures(1, &vol_tex);
30 }
32 Volume::~Volume()
33 {
34 glDeleteTextures(1, &vol_tex);
35 }
37 void Volume::create_vol_tex() const
38 {
39 if(slices.empty())
40 return;
42 glBindTexture(GL_TEXTURE_3D, vol_tex);
43 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MIN_FILTER, GL_LINEAR);
44 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MAG_FILTER, GL_LINEAR);
45 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE);
46 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE);
47 glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_R, GL_CLAMP_TO_EDGE);
49 assert(glGetError() == GL_NO_ERROR);
50 glTexImage3D(GL_TEXTURE_3D, 0, GL_LUMINANCE, width, height,
51 slices.size(), 0, GL_LUMINANCE, GL_FLOAT, 0);
52 assert(glGetError() == GL_NO_ERROR);
54 for(size_t i=0; i<slices.size(); i++) {
55 glTexSubImage3D(GL_TEXTURE_3D, 0, 0, 0, i, width, height,
56 1, GL_LUMINANCE, GL_FLOAT, slices[i].get_pixels());
57 }
58 }
60 void Volume::create_histogram()
61 {
62 if (slices.empty())
63 return;
65 for (size_t i=0; i<slices.size(); i++)
66 {
67 float *slice_pixels = (float*)slices[i].get_pixels();
68 for (int j=0; j < width * height; j++)
69 {
70 float val = slice_pixels[j];
72 int ival = (int)(val * 256);
73 histogram[ival]++;
74 }
75 }
77 max_histogram_value = 0;
78 for (int i=0; i<HIST_SIZE; i++) {
79 if (i < 10)
80 continue;
81 max_histogram_value = std::max(max_histogram_value, histogram[i]);
82 }
83 }
85 bool Volume::load(const char *fname)
86 {
87 char up = 'y';
89 FILE *fp = fopen(fname, "r");
90 if(!fp) {
91 fprintf(stderr, "Failed to open file: %s: %s\n", fname, strerror(errno));
92 return false;
93 }
95 char buf[512];
96 if(!fgets(buf, sizeof buf, fp)) {
97 fprintf(stderr, "Empty file: %s.\n", fname);
98 return false;
99 }
100 if(strstr(buf, "VOLUME") != buf) {
101 fprintf(stderr, "Invalid volume file format: %s\n", fname);
102 return false;
103 }
105 bool reading_slices = false;
106 std::vector<std::string> img_fnames;
108 while(fgets(buf, sizeof buf, fp)) {
109 char *line = strip_whitespaces(buf);
110 if(!line || *line == 0 || *line == '#')
111 continue;
113 if(reading_slices == false) {
114 if(strcmp(line, "SLICES") == 0)
115 reading_slices = true;
117 sscanf(line, "zaspect = %f", &zaspect);
118 sscanf(line, "up = %c", &up);
119 }
120 else {
121 img_fnames.push_back(line);
122 }
123 }
125 fclose(fp);
127 for(size_t i=0; i<img_fnames.size(); i++) {
128 Image img;
129 if(img.load(img_fnames[i].c_str())) {
130 push_slice(std::move(img));
131 }
132 else {
133 fprintf(stderr, "Failed to load slice: %s\n", img_fnames[i].c_str());
134 }
135 }
137 if (up == 'z')
138 x_rot = 90;
140 create_histogram();
141 return true;
142 }
144 bool Volume::push_slice(Image &&slice)
145 {
146 if(slices.empty()) {
147 width = slice.get_width();
148 height = slice.get_height();
149 }
150 else {
151 if(width != slice.get_width() || height != slice.get_height()) {
152 fprintf(stderr, "failed to load slice no: %d\n", (int)slices.size() + 1);
153 return false;
154 }
155 }
157 slices.push_back(slice);
158 return true;
159 }
161 const Image *Volume::get_slice(int num_slice) const
162 {
163 if(num_slice >= (int)slices.size() || num_slice < 0)
164 return 0;
166 return &slices[num_slice];
167 }
169 const Image *Volume::get_slice_by_z(float z) const
170 {
171 int idx = get_slice_idx_by_z(z);
172 return get_slice(idx);
173 }
175 int Volume::get_slice_count() const
176 {
177 return (int)slices.size();
178 }
180 int Volume::get_slice_idx_by_z(float z) const
181 {
182 int last_slice_idx = (int)slices.size() - 1;
183 int idx = z * last_slice_idx;
184 if(idx < 0) {
185 idx = 0;
186 }
187 else if(idx > last_slice_idx) {
188 idx = last_slice_idx;
189 }
191 return idx;
192 }
194 unsigned int Volume::get_texture() const
195 {
196 if(!vol_tex_valid) {
197 create_vol_tex();
198 vol_tex_valid = true;
199 }
201 return vol_tex;
202 }
204 float Volume::get_volume_rotation() const
205 {
206 return x_rot;
207 }
209 static Mesh *cur_mesh;
210 static Volume *cur_vol;
211 static float low_thres, high_thres;
213 static float smoothstep(float x, float low, float high)
214 {
215 if(x < low) return 0.0;
216 if(x >= high) return 1.0;
218 x = (x - low) / (high - low);
219 return x * x * (3.0 - 2.0 * x);
220 }
222 float transfer_function(float x, float low_thres, float high_thres)
223 {
224 float dt = 0.25 * (high_thres - low_thres);
225 return smoothstep(x, low_thres - dt, low_thres + dt) *
226 (1 - smoothstep(x, high_thres - dt, high_thres + dt));
227 }
229 void Volume::create_mesh(Mesh *mesh, float tmin, float tmax, int xres, int yres, int zres)
230 {
231 if (tmin > tmax) {
232 float t = tmax;
233 tmax = tmin;
234 tmin = t;
235 }
236 low_thres = tmin; high_thres = tmax;
238 if(xres <= 0)
239 xres = width;
240 if(yres <= 0)
241 yres = height;
242 if(zres <= 0)
243 zres = slices.size();
245 cur_mesh = mesh;
246 cur_vol = this;
247 }
249 void Volume::draw() const
250 {
251 }
253 static char *strip_whitespaces(char *buf)
254 {
255 while(*buf && isspace(*buf))
256 buf++;
258 if(!*buf)
259 return 0;
261 char *end = buf + strlen(buf) - 1;
262 while(end > buf && isspace(*end))
263 end--;
264 end[1] = 0;
265 return buf;
266 }