Main Page | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members | Related Pages

octree.cpp

Go to the documentation of this file.
00001 
00002 //                                                                    //
00003 // Informatikpraktikum I "Implementation of the ShearWarp Algorithm"  //
00004 // Februar - October 2002                                             //
00005 //                                                                    //
00006 // author: Sebastian Zambal                                           //
00007 //         e9826978@student.tuwien.ac.at                              //
00008 //                                                                    //
00009 // file:   shearWarp.cpp (regular)                                    //
00010 //                                                                    //
00012 
00013 #include "octree.h"
00014 #include "math.h"
00015 
00016 //-----------------------------------------------------------------------
00017 //-- Constructor                                                      ---
00018 //-----------------------------------------------------------------------
00019 Octree::Octree(byte *volumeData, int dim1Size, int dim2Size, int dim3Size)
00020 {
00021     width  = dim1Size;
00022     height = dim2Size;
00023     depth  = dim3Size;
00024     data   = volumeData;
00025 
00026     step_dim3 = dim1Size * dim2Size;
00027 
00028     int max = width;
00029     if (height > max) max = height;
00030     if (depth  > max) max = depth;
00031 
00032     int treeDepth = (int)(log((double)max)/log(2.0f)) - 1;
00033     if (treeDepth < 1) treeDepth = 1;
00034 
00035     root = build(0,0,0,width-1,height-1,depth-1,treeDepth);
00036 }
00037 
00038 //-----------------------------------------------------------------------
00039 //-- Destructor                                                       ---
00040 //-----------------------------------------------------------------------
00041 Octree::~Octree() 
00042 {
00043     remove(root);
00044 }
00045 
00046 //-----------------------------------------------------------------------
00047 //-- Delete octree, free memory                                       ---
00048 //-----------------------------------------------------------------------
00049 void Octree::remove(Node *tree) 
00050 {
00051     for (int i = 0; i < 8; ++i) {
00052         if (tree->children[i] != 0) remove(tree->children[i]);
00053     }
00054 
00055     delete tree;
00056     tree = NULL;
00057 }
00058 
00059 //-----------------------------------------------------------------------
00060 //-- Build octree                                                     ---
00061 //-----------------------------------------------------------------------
00062 Node *Octree::build(int minX, int minY, int minZ,
00063                     int maxX, int maxY, int maxZ, int count)
00064 {
00065     Node *result;
00066 
00067     int partX = 1;
00068     int partY = 1;
00069     int partZ = 1;
00070 
00071     if (minX > maxX) minX = maxX;
00072     if (minY > maxY) minY = maxY;
00073     if (minZ > maxZ) minZ = maxZ;
00074 
00075     result = new Node();                 // create new node/subcube
00076 
00077     result->maxX = maxX;                 // save the coordinates of the subcube
00078     result->maxY = maxY;
00079     result->maxZ = maxZ;
00080     result->minX = minX;
00081     result->minY = minY;
00082     result->minZ = minZ;
00083 
00084     if (count == 0) {
00085 
00086         for (int i = 0; i < 8; ++i) {
00087             result->children[i] = 0;                           // no children needed
00088         }
00089         result->noChildren = 1;
00090 
00091         result->min = 255;
00092         result->max = 0;
00093 
00094         for (int x = minX; x <= maxX; ++x) {
00095             for (int y = minY; y <= maxY; ++y) {
00096                 for (int z = minZ; z <= maxZ; ++z) {
00097                     byte voxel = data[x+y*width+z*step_dim3];
00098                     if (voxel > result->max) result->max = voxel;
00099                     if (voxel < result->min) result->min = voxel;
00100                 }
00101             }
00102         }
00103 
00104     } else {
00105         if ((minX == maxX) && (minY == maxY) && (minZ == maxZ)) {
00106             // domain represents only one voxel
00107 
00108             for (int i = 0; i < 8; ++i) {
00109                 result->children[i] = 0;                           // no children needed
00110             }
00111             result->noChildren = 1;
00112 
00113             result->min = result->max = data[minX+minY*width+minZ*step_dim3]; // get the value of the voxel
00114 
00115         } else {
00116 
00117             int midX = (minX + maxX) / 2;
00118             int midY = (minY + maxY) / 2;
00119             int midZ = (minZ + maxZ) / 2;
00120 
00121             if (minX == maxX) partX = 0;
00122             if (minY == maxY) partY = 0;
00123             if (minZ == maxZ) partZ = 0;
00124 
00125             for (int i = 0; i < 8; ++i) result->children[i] = 0;
00126 
00127             result->children[0] = build(minX,   minY,   minZ,   midX, midY, midZ, count-1);
00128             if (partX == 1)
00129                 result->children[1] = build(midX+1, minY,   minZ,   maxX, midY, midZ, count-1);
00130             else result->children[1] = 0;
00131             if (partY == 1)
00132                 result->children[2] = build(minX,   midY+1, minZ,   midX, maxY, midZ, count-1);
00133             else result->children[2] = 0;
00134             if ((partX == 1) && (partY == 1)) 
00135                 result->children[3] = build(midX+1, midY+1, minZ,   maxX, maxY, midZ, count-1);
00136             else result->children[3] = 0;
00137             if (partZ == 1)
00138                 result->children[4] = build(minX,   minY,   midZ+1, midX, midY, maxZ, count-1);
00139             else result->children[4] = 0;
00140             if ((partX == 1) && (partZ == 1)) 
00141                 result->children[5] = build(midX+1, minY,   midZ+1, maxX, midY, maxZ, count-1);
00142             else result->children[5] = 0;
00143             if ((partY == 1) && (partZ == 1))
00144                 result->children[6] = build(minX,   midY+1, midZ+1, midX, maxY, maxZ, count-1);
00145             else result->children[6] = 0;
00146             if ((partX == 1) && (partY == 1) && (partZ == 1))
00147                 result->children[7] = build(midX+1, midY+1, midZ+1, maxX, maxY, maxZ, count-1);
00148             else result->children[7] = 0;
00149 
00150             result->noChildren = 0;                 // there is at least one child-node
00151 
00152             result->min = 255;
00153             result->max = 0;
00154 
00155             for (int i = 0; i < 8; ++i) {
00156                 if (result->children[i] != 0) {
00157                     if (result->children[i]->min < result->min) result->min = result->children[i]->min;
00158                     if (result->children[i]->max > result->max) result->max = result->children[i]->max;
00159                 }
00160             }
00161         }
00162     }
00163 
00164     return result;
00165 }
00166 
00167 //-----------------------------------------------------------------------
00168 //-- Do classification                                                ---
00169 //-----------------------------------------------------------------------
00170 void Octree::classify(SummedAreaTable *sat)
00171 {
00172     printf("Starting classification");
00173     classify(sat,root);
00174     printf("ready.");
00175 }
00176 
00177 //-----------------------------------------------------------------------
00178 //-- Do classification recursively                                    ---
00179 //-----------------------------------------------------------------------
00180 void Octree::classify(SummedAreaTable *sat, Node *tree)
00181 {
00182     if (sat->transparent(tree->min,tree->max) == 1) {   // there are only transparent voxels in the subcube
00183         tree->classification = ALL_TRANSPARENT;
00184     } else {                                            // there are also non-transparent voxels in the subcube
00185         tree->classification = ALL_NON_TRANSPARENT;
00186         for (int i = 0; i < 8; ++i) {
00187             if (tree->children[i] != 0) {
00188                 tree->classification = COMBINATION;
00189                 classify(sat,tree->children[i]);
00190             }
00191         }
00192     }
00193 }
00194 
00195 //-----------------------------------------------------------------------
00196 //-- Skip transparent voxels                                          ---
00197 //-----------------------------------------------------------------------
00198 int Octree::skip(int x, int y, int z, int mainViewingDir)
00199 {
00200     return skip(x, y, z, mainViewingDir, root);
00201 }
00202 
00203 //-----------------------------------------------------------------------
00204 //-- Skip transparent voxels recursively                              ---
00205 //-----------------------------------------------------------------------
00206 int Octree::skip(int x, int y, int z, int mainViewingDir, Node *tree)
00207 {
00208     int result = 1;
00209 
00210     if (tree->classification == ALL_TRANSPARENT) {    // only transparent voxels => skip the whole subtree
00211         switch(mainViewingDir) {
00212             case XDIR: {
00213                 result = tree->maxY - y;
00214                 break;
00215             }
00216             case YDIR: {
00217                 result = tree->maxZ - z;
00218                 break;
00219             }
00220             case ZDIR: {
00221                 result = tree->maxX - x;
00222                 break;
00223             }
00224         }
00225     } else {
00226 
00227         if (tree->classification == COMBINATION) {   // both, transparent and non-transparent voxels
00228 
00229             int midX = (tree->minX + tree->maxX) / 2;
00230             int midY = (tree->minY + tree->maxY) / 2;
00231             int midZ = (tree->minZ + tree->maxZ) / 2;
00232 
00233             if (x <= midX) {
00234                 if (y <= midY) {
00235                     if (z <= midZ) {
00236                         result = skip(x, y, z, mainViewingDir, tree->children[0]);
00237                     } else {
00238                         result = skip(x, y, z, mainViewingDir, tree->children[4]);
00239                     }
00240                 } else {
00241                     if (z <= midZ) {
00242                         result = skip(x, y, z, mainViewingDir, tree->children[2]);
00243                     } else {
00244                         result = skip(x, y, z, mainViewingDir, tree->children[6]);
00245                     }
00246                 }
00247             } else {
00248                 if (y <= midY) {
00249                     if (z <= midZ) {
00250                         result = skip(x, y, z, mainViewingDir, tree->children[1]);
00251                     } else {
00252                         result = skip(x, y, z, mainViewingDir, tree->children[5]);
00253                     }
00254                 } else {
00255                     if (z <= midZ) {
00256                         result = skip(x, y, z, mainViewingDir, tree->children[3]);
00257                     } else {
00258                         result = skip(x, y, z, mainViewingDir, tree->children[7]);
00259                     }
00260                 }
00261             }
00262         } else {                                         // only non-transparent voxels
00263             result = 1;
00264         }
00265     }
00266 
00267     if (result <= 0) result = 1;
00268 
00269     return result;
00270 }
00271 
00272 //-----------------------------------------------------------------------
00273 //-- Constructor                                                      ---
00274 //-----------------------------------------------------------------------
00275 SummedAreaTable::SummedAreaTable() 
00276 {
00277     table = new int[256];
00278 }
00279 
00280 //-----------------------------------------------------------------------
00281 //-- Destructor                                                       ---
00282 //-----------------------------------------------------------------------
00283 SummedAreaTable::~SummedAreaTable() 
00284 {
00285     if (table) delete [] table;
00286     table = NULL;
00287 }
00288 
00289 //-----------------------------------------------------------------------
00290 //-- Build the summed area table                                      ---
00291 //-----------------------------------------------------------------------
00292 void SummedAreaTable::build(float *transferFunction, float threshold_runlength)
00293 {
00294     byte current = 0;
00295     for (int i = 0; i < 256; ++i) {
00296         if (transferFunction[i] > threshold_runlength) {
00297             current++;
00298         }
00299         table[i] = current;
00300     }
00301 }
00302 
00303 //-----------------------------------------------------------------------
00304 //-- Given interval transparent?                                      ---
00305 //-----------------------------------------------------------------------
00306 int SummedAreaTable::transparent(byte min, byte max) 
00307 {
00308     int result;
00309 
00310     if (min == max) {
00311         if (max == 0) {
00312             if (table[max] >= 1) result = 0; else result = 1;
00313         } else {
00314             if (table[max] > table[max-1]) result = 0; else result = 1;
00315         }
00316     } else {
00317         if (table[max] == table[min]) {   // transparent domain
00318             result = 1;
00319         } else {
00320             result = 0;
00321         }
00322     }
00323 
00324     return result;
00325 }
00326 
00327 //-----------------------------------------------------------------------
00328 //-- Constructor                                                      ---
00329 //-----------------------------------------------------------------------
00330 FastClassification::FastClassification(byte *volumeData, int dim1Size, int dim2Size, int dim3Size)
00331 {
00332     sTable = new SummedAreaTable();
00333     octree = new Octree(volumeData, dim1Size, dim2Size, dim3Size);
00334 
00335     width  = dim1Size;
00336     height = dim2Size;
00337     depth  = dim3Size;
00338 
00339     data = volumeData;
00340 }
00341 
00342 //-----------------------------------------------------------------------
00343 //-- Destructor                                                       ---
00344 //-----------------------------------------------------------------------
00345 FastClassification::~FastClassification() 
00346 {
00347   if (octree) delete(octree); // octree->~Octree();
00348   if (sTable) delete(sTable); //  sTable->~SummedAreaTable();
00349   octree = NULL;
00350   sTable = NULL;
00351 }
00352 
00353 //-----------------------------------------------------------------------
00354 //-- Do classification                                                ---
00355 //-----------------------------------------------------------------------
00356 void FastClassification::classify() 
00357 {
00358     octree->classify(sTable);
00359 }
00360 
00361 //-----------------------------------------------------------------------
00362 //-- Skip transparent voxels                                          ---
00363 //-----------------------------------------------------------------------
00364 int FastClassification::skip(int x, int y, int z, int mainViewingDirection)
00365 {
00366     int result;
00367 
00368     result = octree->skip(x, y, z, mainViewingDirection);
00369 
00370     return result;
00371 }
00372 
00373 //-----------------------------------------------------------------------
00374 //-- Build summed area table                                          ---
00375 //-----------------------------------------------------------------------
00376 void FastClassification::buildSummedAreaTable(float *transferFunction, float threshold_runlength)
00377 {
00378     sTable->build(transferFunction,threshold_runlength);
00379 }
00380 
00381 
00382 
00383 
00384 
00385 
00386 
00387 
00388 
00389 
00390 
00391 
00392 
00393 

Generated on Wed Dec 15 21:20:30 2004 for vuVolume by  doxygen 1.3.9.1