source: trunk/src/testing/app/volvis/vFileVolume.cpp @ 4

Revision 4, 8.7 KB checked in by ajaworski, 13 years ago (diff)

Added modified SAGE sources

Line 
1/********************************************************************************
2 * Volatile - Volume Visualization Software for SAGE
3 * Copyright (C) 2004 Electronic Visualization Laboratory,
4 * University of Illinois at Chicago
5 *
6 * All rights reserved.
7 *
8 * Redistribution and use in source and binary forms, with or without
9 * modification, are permitted provided that the following conditions are met:
10 *
11 *  * Redistributions of source code must retain the above copyright
12 *    notice, this list of conditions and the following disclaimer.
13 *  * Redistributions in binary form must reproduce the above
14 *    copyright notice, this list of conditions and the following disclaimer
15 *    in the documentation and/or other materials provided with the distribution.
16 *  * Neither the name of the University of Illinois at Chicago nor
17 *    the names of its contributors may be used to endorse or promote
18 *    products derived from this software without specific prior written permission.
19 *
20 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
23 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
24 * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
25 * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
26 * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
27 * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
28 * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
29 * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
30 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
31 *
32 * Direct questions, comments etc about Volatile to www.evl.uic.edu/cavern/forum
33 *********************************************************************************/
34
35//temporarily commented for the demo
36#if 0
37//Class vFileVolume: Defines the 3D volume thats loaded from a file
38//This class only manages the data(doesnt do the rendering)
39//Shalini Venkataraman, Luc Renambot
40//Sep 7 2003
41#include <math.h>
42#include <stdio.h>
43#include <stdlib.h>
44#include <string.h>
45#include "glUE.h"
46#include "vVolume.h"
47#include "VectorMath.h"
48#include "global.h"
49
50vFileVolume::vFileVolume(const char* fname, int c, int r, int d):vVolume(c,r,d) {
51        strcpy(this->filename,fname);
52        padY = padX = padZ =0;
53        filename[0] = '\0';
54        gradientFile = 0;
55        fileBuffer = 0;
56        //vVolume(c,r,d);
57}
58
59vFileVolume::~vFileVolume() {
60        if (gradientFile)
61                delete [] gradientFile;
62}
63
64
65void vFileVolume::load(int loadx, int loady, int loadz) {
66        setLoadDims(loadx,loady,loadz);
67        voxelData = new GLubyte[dimZ*dimY*dimX];
68        loadVolumeFromFile(filename);
69        if (gradientFile)
70                loadGradientFromFile(gradientFile);
71        else
72                calcGradient(); //create grad volume
73}
74
75void vFileVolume::loadVolumeFromFile(const char* filename){
76        int rc= mmap.setFilename(filename, 0);
77        cout << " == " << rc <<endl;
78        // Call mmap
79        fileBuffer= (unsigned char *)mmap.doMmap();
80        // Did the mapping work?
81        if((char*)fileBuffer== (char *)-1)
82        {
83                cerr << "mmap failed -- No memory :-( " << endl;
84                exit(1);
85        }
86
87        printf("Loading Volume %s\n",filename);
88        int filePlaneDims = fileX*fileY;//plane in a volume file
89        int planeDims = dimX*dimY; //plane in the volume that we want to load
90        GLubyte* plane = new GLubyte[planeDims];
91        padVolume(); //pad volume if necessary
92        //if we padded, just fill the dimesions in the file
93        int fillZ = padZ?fileZ:dimZ;
94        int fillY = padY?fileY:dimY;
95        int fillX = padX?fileX:dimX;
96        //if there is a dimZ to skip, skip it
97        unsigned char* curFilePtr = fileBuffer+filePlaneDims*offsetZ*sizeof(GLubyte);
98        //load the volume one plane at a time
99        for (int r=0;r<fillZ;r++) {
100                readASlice(plane,curFilePtr);
101                for (int t=0;t< fillY;t++) {
102                        for (int s=0;s< fillX;s++) {
103                                setValue(r+padZ,t+padY,s+padX,plane[(t*fillX+s)]);
104                        }
105                }
106        }
107        printf("Loading Textures\n");
108        delete [] plane;
109}
110
111
112//if the volume is bigger than the file size, this will pad the volume
113//else do nothing
114void vFileVolume::padVolume() {
115        unsigned int r, s, t;
116        int lastSlice;
117        //pad a little with black if fileZ < dimZ
118        if (dimZ>fileZ) {
119                padZ = (dimZ-fileZ)/2;
120                lastSlice = (dimZ-padZ);
121                for (r=0;r<padZ;r++){
122                        for (t=0;t<dimY;t++){
123                                for (s=0;s<dimX;s++){
124                                        setValue(r,t,s,0);
125                                        setValue(r+lastSlice,t,s,0);
126                                }
127                        }
128                }
129        }
130        if (dimY>fileY) {
131                //pad a little with black if fileY doesnt match dimY
132                padY = (dimY-fileY)/2;
133                lastSlice = (dimY-padY);
134                for (r=0;r<dimZ;r++){
135                        for (t=0;t<padY;t++){
136                                for (s=0;s<dimX;s++){
137                                        setValue(r,t,s,0);
138                                        setValue(r,t+lastSlice,s,0);
139                                }
140                        }
141                }
142        }
143
144        if (dimX>fileX) {
145                //pad a little with black if fileX doesnt match dimX
146                padX = (dimX-fileX)/2;
147                lastSlice = (dimX-padX);
148                for (r=0;r<dimZ;r++){
149                        for (t=0;t<dimY;t++){
150                                for (s=0;s<padX;s++){
151                                        setValue(r,t,s,0);
152                                        setValue(r,t,s+lastSlice,0);
153                                }
154                        }
155                }
156        }
157}
158
159//reads a slice skipping the Y and X offsets
160void vFileVolume::readASlice(GLubyte* plane, unsigned char* &curFilePtr) {
161        //calculating the remaining rows and colums to skip once a sub image is read in
162        int remainingX = fileX-offsetX-dimX;
163        int remainingY = fileY-offsetY-dimY;
164
165        //skip the rows on one slice to start the read
166        curFilePtr += fileX*offsetY*sizeof(GLubyte);
167        GLubyte* volDataPtr = plane;
168        for (int i=0;i<dimY;i++) {
169                //skip the offset X
170                curFilePtr += offsetX*sizeof(GLubyte);
171                //read the sub-cols
172                memcpy(volDataPtr,curFilePtr,dimX);
173                //skip the cols read
174                curFilePtr += dimX*sizeof(GLubyte);
175                //skip the remaining cols
176                curFilePtr += remainingX*sizeof(GLubyte);
177                volDataPtr += dimX;//skip in the volume
178        }
179        //skip the remaining rows
180        curFilePtr += remainingY*fileX*sizeof(GLubyte);
181        fprintf(stderr,"ok here\n");
182}
183
184
185
186void vFileVolume::loadGradientFromFile(const char* filename) {
187        FILE* fp;
188        if (!(fp=fopen(filename,"rb"))) {
189                printf("Error: unable to open volume file %s\n",filename);
190                exit(0);
191        }
192        else
193                fprintf(stderr,"Opening file %s\n",filename);
194        int volSize = dimY*dimX*dimZ;
195        gradData = new unsigned char[volSize];
196        int numRead = fread(gradData,sizeof(unsigned char), volSize,fp);
197        if (numRead != volSize)
198                fprintf(stderr,"Error: gradient data size %d doesnt match volume size %d\n",numRead,volSize);
199        fclose(fp);
200}
201
202//calculate the value, gradient components of the volume
203//- may be needed by the transfer functions
204void vFileVolume::calcGradient()
205{
206        unsigned char *vdata = (unsigned char*)voxelData;
207        float *gradV3 = new float[dimX*dimY*dimZ*3]; //individual components of the gradient
208        float *gmag = new float[dimX*dimY*dimZ]; //gradient magnitude
209        float gmmax = -100000000;
210        float gmmin = 100000000;
211        float dmax = -100000000;
212        float dmin = 100000000;
213
214
215        //compute the gradient
216        cerr << "   computing 1st derivative";
217        int i;
218        for(i = 0; i < dimZ; ++i){
219                for(int j = 0; j < dimY; ++j){
220                        for(int k = 0; k < dimX; ++k){
221                                if((k<1)||(k>dimX-2)||(j<1)||(j>dimY-2)||(i<1)||(i>dimZ-2)){
222                                        gradV3[i*dimX*dimY*3 + j*dimX*3 + k*3 + 0] = 0;
223                                        gradV3[i*dimX*dimY*3 + j*dimX*3 + k*3 + 1] = 0;
224                                        gradV3[i*dimX*dimY*3 + j*dimX*3 + k*3 + 2] = 0;
225                                        gmag[i*dimX*dimY + j*dimX + k] = 0;
226                                }
227                                else {
228                                        float dx = (float)(vdata[i*dimX*dimY + j*dimX + (k+1)]
229                                                                        - vdata[i*dimX*dimY + j*dimX + (k-1)]);
230                                        float dy = (float)(vdata[i*dimX*dimY + (j+1)*dimX + k]
231                                                                        - vdata[i*dimX*dimY + (j-1)*dimX + k]);
232                                        float dz = (float)(vdata[(i+1)*dimX*dimY + j*dimX + k]
233                                                                        - vdata[(i-1)*dimX*dimY + j*dimX + k]);
234
235                                        gradV3[i*dimX*dimY*3 + j*dimX*3 + k*3 + 0] = dx;
236                                        gradV3[i*dimX*dimY*3 + j*dimX*3 + k*3 + 1] = dy;
237                                        gradV3[i*dimX*dimY*3 + j*dimX*3 + k*3 + 2] = dz;
238                                        gmag[i*dimX*dimY + j*dimX + k] = (float)sqrt(dx*dx+dy*dy+dz*dz);
239                                        gmmax = MAX(gmag[i*dimX*dimY + j*dimX + k], gmmax);
240                                        gmmin = MIN(gmag[i*dimX*dimY + j*dimX + k], gmmin);
241                                        dmax = MAX(vdata[i*dimX*dimY + j*dimX +k], dmax);
242                                        dmin = MIN(vdata[i*dimX*dimY + j*dimX +k], dmin);
243                                }
244                        }
245                }
246        }
247        cerr << " - done" << endl;
248        cerr << "   quantizing";
249        gradData = new unsigned char[dimX*dimY*dimZ]; //stores the 1st derivative
250
251        for(i = 0; i < dimZ; ++i){
252                for(int j = 0; j < (dimY); ++j){
253                        for(int k = 0; k < (dimX); ++k){
254                                if((k<1)||(k>dimX-2)||(j<1)||(j>dimY-2)||(i<1)||(i>dimZ-2)){
255                                        gradData[i*dimX*dimY + j*dimX + k] = 0;
256                                }
257                                else {
258                                        gradData[i*dimX*dimY + j*dimX + k] = (unsigned char)affine(gmmin, gmag[i*dimX*dimY + j*dimX + k], gmmax, 0, 255);
259                                }
260                        }
261                }
262        }
263        cerr << " - done" << endl;
264        delete [] gradV3;
265        delete [] gmag;
266}
267#endif
Note: See TracBrowser for help on using the repository browser.