marchingcubes.cpp

Go to the documentation of this file.
00001 /*
00002  * Author: Tao Ju, 5/16/2007 <taoju@cs.wustl.edu>, code ported by Grant Tang
00003  * code extensively modified and optimized by David Woolford
00004  * Copyright (c) 2000-2006 Baylor College of Medicine
00005  *
00006  * This software is issued under a joint BSD/GNU license. You may use the
00007  * source code in this file under either license. However, note that the
00008  * complete EMAN2 and SPARX software packages have some GPL dependencies,
00009  * so you are responsible for compliance with the licenses of these packages
00010  * if you opt to use BSD licensing. The warranty disclaimer below holds
00011  * in either instance.
00012  *
00013  * This complete copyright notice must be included in any revised version of the
00014  * source code. Additional authorship citations may be added, but existing
00015  * author citations must be preserved.
00016  *
00017  * This program is free software; you can redistribute it and/or modify
00018  * it under the terms of the GNU General Public License as published by
00019  * the Free Software Foundation; either version 2 of the License, or
00020  * (at your option) any later version.
00021  *
00022  * This program is distributed in the hope that it will be useful,
00023  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00024  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
00025  * GNU General Public License for more details.
00026  *
00027  * You should have received a copy of the GNU General Public License
00028  * along with this program; if not, write to the Free Software
00029  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
00030  *
00031  * */
00032 #ifdef EMAN2_USING_OPENGL
00033 
00034 #ifdef _WIN32
00035         #include <windows.h>
00036 #endif  //_WIN32
00037 
00038 #ifndef GL_GLEXT_PROTOTYPES
00039         #define GL_GLEXT_PROTOTYPES
00040 #endif  //GL_GLEXT_PROTOTYPES
00041 
00042 
00043 #include "marchingcubes.h"
00044 
00045 #include <time.h>
00046 #include <math.h>
00047 
00048 using namespace EMAN;
00049 #include "transform.h"
00050 #include "emobject.h"
00051 #include "vec3.h"
00052 
00053 #define min(a,b)(((a) < (b)) ? (a) : (b))
00054 #define max(a,b)(((a) > (b)) ? (a) : (b))
00055 
00056 
00057 //a2fVertexOffset lists the positions, relative to vertex0, of each of the 8 vertices of a cube
00058 static const int a2fVertexOffset[8][3] =
00059 {
00060                 {0, 0, 0},{1, 0, 0},{1, 1, 0},{0, 1, 0},
00061                 {0, 0, 1},{1, 0, 1},{1, 1, 1},{0, 1, 1}
00062 };
00063 
00064 
00065 static const int a2fPosXOffset[4][3] =
00066 {
00067         {2, 0, 0},{2, 1, 0},
00068         {2, 0, 1},{2, 1, 1}
00069 };
00070 
00071 static const int a2fPosYOffset[4][3] =
00072 {
00073         {1, 2, 0},{0, 2, 0},
00074         {1, 2, 1},{0, 2, 1}
00075 };
00076 
00077 static const int a2fPosZOffset[4][3] =
00078 {
00079         {0, 0, 2},{1, 0, 2},
00080         {1, 1, 2},{0, 1, 2}
00081 };
00082 
00083 static const int a2fPosXPosYOffset[2][3] =
00084 {
00085         {2, 2, 0},{2, 2, 1}
00086 };
00087 
00088 static const int a2fPosXPosZOffset[2][3] =
00089 {
00090         {2, 0, 2},{2, 1, 2}
00091 };
00092 
00093 static const int a2fPosYPosZOffset[2][3] =
00094 {
00095         {1, 2, 2},{0, 2, 2}
00096 };
00097 
00098 
00099 static const int a2fPosXPosZPosYOffset[3] =
00100 {
00101         2, 2, 2
00102 };
00103 
00104 //a2fVertexOffset lists the positions, relative to vertex0, of each of the 8 vertices of a cube
00105 static const int a2OddXOffset[8] =
00106 {
00107         1,0,0,1,
00108         1,0,0,1
00109 };
00110 
00111 static const int a2OddYOffset[8] =
00112 {
00113         1,1,0,0,
00114         1,1,0,0
00115 };
00116 
00117 static const int a2OddZOffset[8] =
00118 {
00119         1,1,1,1,
00120         0,0,0,0
00121 };
00122 
00123 //a2iEdgeConnection lists the index of the endpoint vertices for each of the 12 edges of the cube
00124 static const int a2iEdgeConnection[12][2] =
00125 {
00126                 {0,1}, {1,2}, {2,3}, {3,0},
00127                 {4,5}, {5,6}, {6,7}, {7,4},
00128                 {0,4}, {1,5}, {2,6}, {3,7}
00129 };
00130 
00131 //a2fEdgeDirection lists the direction vector (vertex1-vertex0) for each edge in the cube
00132 static const float a2fEdgeDirection[12][3] =
00133 {
00134                 {1.0, 0.0, 0.0},{0.0, 1.0, 0.0},{-1.0, 0.0, 0.0},{0.0, -1.0, 0.0},
00135                 {1.0, 0.0, 0.0},{0.0, 1.0, 0.0},{-1.0, 0.0, 0.0},{0.0, -1.0, 0.0},
00136                 {0.0, 0.0, 1.0},{0.0, 0.0, 1.0},{ 0.0, 0.0, 1.0},{0.0,  0.0, 1.0}
00137 };
00138 
00139 // right = 0, up = 1, back = 2
00140 static const int edgeLookUp[12][4] =
00141 {
00142         {0,0,0,0},{1,0,0,1},{0,1,0,0},{0,0,0,1},
00143         {0,0,1,0},{1,0,1,1},{0,1,1,0},{0,0,1,1},
00144         {0,0,0,2},{1,0,0,2},{1,1,0,2},{0,1,0,2}
00145 };
00146 
00147 // For any edge, if one vertex is inside of the surface and the other is outside of the surface
00148 // then the edge intersects the surface
00149 // For each of the 8 vertices of the cube can be two possible states : either inside or outside of the surface
00150 // For any cube the are 2^8=256 possible sets of vertex states
00151 // This table lists the edges intersected by the surface for all 256 possible vertex states
00152 // There are 12 edges.  For each entry in the table, if edge #n is intersected, then bit #n is set to 1
00153 
00154 int aiCubeEdgeFlags[256]=
00155 {
00156         0x000, 0x109, 0x203, 0x30a, 0x406, 0x50f, 0x605, 0x70c, 0x80c, 0x905, 0xa0f, 0xb06, 0xc0a, 0xd03, 0xe09, 0xf00,
00157         0x190, 0x099, 0x393, 0x29a, 0x596, 0x49f, 0x795, 0x69c, 0x99c, 0x895, 0xb9f, 0xa96, 0xd9a, 0xc93, 0xf99, 0xe90,
00158         0x230, 0x339, 0x033, 0x13a, 0x636, 0x73f, 0x435, 0x53c, 0xa3c, 0xb35, 0x83f, 0x936, 0xe3a, 0xf33, 0xc39, 0xd30,
00159         0x3a0, 0x2a9, 0x1a3, 0x0aa, 0x7a6, 0x6af, 0x5a5, 0x4ac, 0xbac, 0xaa5, 0x9af, 0x8a6, 0xfaa, 0xea3, 0xda9, 0xca0,
00160         0x460, 0x569, 0x663, 0x76a, 0x066, 0x16f, 0x265, 0x36c, 0xc6c, 0xd65, 0xe6f, 0xf66, 0x86a, 0x963, 0xa69, 0xb60,
00161         0x5f0, 0x4f9, 0x7f3, 0x6fa, 0x1f6, 0x0ff, 0x3f5, 0x2fc, 0xdfc, 0xcf5, 0xfff, 0xef6, 0x9fa, 0x8f3, 0xbf9, 0xaf0,
00162         0x650, 0x759, 0x453, 0x55a, 0x256, 0x35f, 0x055, 0x15c, 0xe5c, 0xf55, 0xc5f, 0xd56, 0xa5a, 0xb53, 0x859, 0x950,
00163         0x7c0, 0x6c9, 0x5c3, 0x4ca, 0x3c6, 0x2cf, 0x1c5, 0x0cc, 0xfcc, 0xec5, 0xdcf, 0xcc6, 0xbca, 0xac3, 0x9c9, 0x8c0,
00164         0x8c0, 0x9c9, 0xac3, 0xbca, 0xcc6, 0xdcf, 0xec5, 0xfcc, 0x0cc, 0x1c5, 0x2cf, 0x3c6, 0x4ca, 0x5c3, 0x6c9, 0x7c0,
00165         0x950, 0x859, 0xb53, 0xa5a, 0xd56, 0xc5f, 0xf55, 0xe5c, 0x15c, 0x055, 0x35f, 0x256, 0x55a, 0x453, 0x759, 0x650,
00166         0xaf0, 0xbf9, 0x8f3, 0x9fa, 0xef6, 0xfff, 0xcf5, 0xdfc, 0x2fc, 0x3f5, 0x0ff, 0x1f6, 0x6fa, 0x7f3, 0x4f9, 0x5f0,
00167         0xb60, 0xa69, 0x963, 0x86a, 0xf66, 0xe6f, 0xd65, 0xc6c, 0x36c, 0x265, 0x16f, 0x066, 0x76a, 0x663, 0x569, 0x460,
00168         0xca0, 0xda9, 0xea3, 0xfaa, 0x8a6, 0x9af, 0xaa5, 0xbac, 0x4ac, 0x5a5, 0x6af, 0x7a6, 0x0aa, 0x1a3, 0x2a9, 0x3a0,
00169         0xd30, 0xc39, 0xf33, 0xe3a, 0x936, 0x83f, 0xb35, 0xa3c, 0x53c, 0x435, 0x73f, 0x636, 0x13a, 0x033, 0x339, 0x230,
00170         0xe90, 0xf99, 0xc93, 0xd9a, 0xa96, 0xb9f, 0x895, 0x99c, 0x69c, 0x795, 0x49f, 0x596, 0x29a, 0x393, 0x099, 0x190,
00171         0xf00, 0xe09, 0xd03, 0xc0a, 0xb06, 0xa0f, 0x905, 0x80c, 0x70c, 0x605, 0x50f, 0x406, 0x30a, 0x203, 0x109, 0x000
00172 };
00173 
00174 //  For each of the possible vertex states listed in aiCubeEdgeFlags there is a specific triangulation
00175 //  of the edge intersection points.  a2iTriangleConnectionTable lists all of thminvals[cur_level]-> in the form of
00176 //  0-5 edge triples with the list terminated by the invalid value -1.
00177 //  For example: a2iTriangleConnectionTable[3] list the 2 triangles formed when corner[0]
00178 //  and corner[1] are inside of the surface, but the rest of the cube is not.
00179 //
00180 //  I found this table in an example program someone wrote long ago.  It was probably generated by hand
00181 
00182 int a2iTriangleConnectionTable[256][16] =
00183 {
00184         {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00185         {0, 8, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00186         {0, 1, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00187         {1, 8, 3, 9, 8, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00188         {1, 2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00189         {0, 8, 3, 1, 2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00190         {9, 2, 10, 0, 2, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00191         {2, 8, 3, 2, 10, 8, 10, 9, 8, -1, -1, -1, -1, -1, -1, -1},
00192         {3, 11, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00193         {0, 11, 2, 8, 11, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00194         {1, 9, 0, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00195         {1, 11, 2, 1, 9, 11, 9, 8, 11, -1, -1, -1, -1, -1, -1, -1},
00196         {3, 10, 1, 11, 10, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00197         {0, 10, 1, 0, 8, 10, 8, 11, 10, -1, -1, -1, -1, -1, -1, -1},
00198         {3, 9, 0, 3, 11, 9, 11, 10, 9, -1, -1, -1, -1, -1, -1, -1},
00199         {9, 8, 10, 10, 8, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00200         {4, 7, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00201         {4, 3, 0, 7, 3, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00202         {0, 1, 9, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00203         {4, 1, 9, 4, 7, 1, 7, 3, 1, -1, -1, -1, -1, -1, -1, -1},
00204         {1, 2, 10, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00205         {3, 4, 7, 3, 0, 4, 1, 2, 10, -1, -1, -1, -1, -1, -1, -1},
00206         {9, 2, 10, 9, 0, 2, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1},
00207         {2, 10, 9, 2, 9, 7, 2, 7, 3, 7, 9, 4, -1, -1, -1, -1},
00208         {8, 4, 7, 3, 11, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00209         {11, 4, 7, 11, 2, 4, 2, 0, 4, -1, -1, -1, -1, -1, -1, -1},
00210         {9, 0, 1, 8, 4, 7, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1},
00211         {4, 7, 11, 9, 4, 11, 9, 11, 2, 9, 2, 1, -1, -1, -1, -1},
00212         {3, 10, 1, 3, 11, 10, 7, 8, 4, -1, -1, -1, -1, -1, -1, -1},
00213         {1, 11, 10, 1, 4, 11, 1, 0, 4, 7, 11, 4, -1, -1, -1, -1},
00214         {4, 7, 8, 9, 0, 11, 9, 11, 10, 11, 0, 3, -1, -1, -1, -1},
00215         {4, 7, 11, 4, 11, 9, 9, 11, 10, -1, -1, -1, -1, -1, -1, -1},
00216         {9, 5, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00217         {9, 5, 4, 0, 8, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00218         {0, 5, 4, 1, 5, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00219         {8, 5, 4, 8, 3, 5, 3, 1, 5, -1, -1, -1, -1, -1, -1, -1},
00220         {1, 2, 10, 9, 5, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00221         {3, 0, 8, 1, 2, 10, 4, 9, 5, -1, -1, -1, -1, -1, -1, -1},
00222         {5, 2, 10, 5, 4, 2, 4, 0, 2, -1, -1, -1, -1, -1, -1, -1},
00223         {2, 10, 5, 3, 2, 5, 3, 5, 4, 3, 4, 8, -1, -1, -1, -1},
00224         {9, 5, 4, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00225         {0, 11, 2, 0, 8, 11, 4, 9, 5, -1, -1, -1, -1, -1, -1, -1},
00226         {0, 5, 4, 0, 1, 5, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1},
00227         {2, 1, 5, 2, 5, 8, 2, 8, 11, 4, 8, 5, -1, -1, -1, -1},
00228         {10, 3, 11, 10, 1, 3, 9, 5, 4, -1, -1, -1, -1, -1, -1, -1},
00229         {4, 9, 5, 0, 8, 1, 8, 10, 1, 8, 11, 10, -1, -1, -1, -1},
00230         {5, 4, 0, 5, 0, 11, 5, 11, 10, 11, 0, 3, -1, -1, -1, -1},
00231         {5, 4, 8, 5, 8, 10, 10, 8, 11, -1, -1, -1, -1, -1, -1, -1},
00232         {9, 7, 8, 5, 7, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00233         {9, 3, 0, 9, 5, 3, 5, 7, 3, -1, -1, -1, -1, -1, -1, -1},
00234         {0, 7, 8, 0, 1, 7, 1, 5, 7, -1, -1, -1, -1, -1, -1, -1},
00235         {1, 5, 3, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00236         {9, 7, 8, 9, 5, 7, 10, 1, 2, -1, -1, -1, -1, -1, -1, -1},
00237         {10, 1, 2, 9, 5, 0, 5, 3, 0, 5, 7, 3, -1, -1, -1, -1},
00238         {8, 0, 2, 8, 2, 5, 8, 5, 7, 10, 5, 2, -1, -1, -1, -1},
00239         {2, 10, 5, 2, 5, 3, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1},
00240         {7, 9, 5, 7, 8, 9, 3, 11, 2, -1, -1, -1, -1, -1, -1, -1},
00241         {9, 5, 7, 9, 7, 2, 9, 2, 0, 2, 7, 11, -1, -1, -1, -1},
00242         {2, 3, 11, 0, 1, 8, 1, 7, 8, 1, 5, 7, -1, -1, -1, -1},
00243         {11, 2, 1, 11, 1, 7, 7, 1, 5, -1, -1, -1, -1, -1, -1, -1},
00244         {9, 5, 8, 8, 5, 7, 10, 1, 3, 10, 3, 11, -1, -1, -1, -1},
00245         {5, 7, 0, 5, 0, 9, 7, 11, 0, 1, 0, 10, 11, 10, 0, -1},
00246         {11, 10, 0, 11, 0, 3, 10, 5, 0, 8, 0, 7, 5, 7, 0, -1},
00247         {11, 10, 5, 7, 11, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00248         {10, 6, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00249         {0, 8, 3, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00250         {9, 0, 1, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00251         {1, 8, 3, 1, 9, 8, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1},
00252         {1, 6, 5, 2, 6, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00253         {1, 6, 5, 1, 2, 6, 3, 0, 8, -1, -1, -1, -1, -1, -1, -1},
00254         {9, 6, 5, 9, 0, 6, 0, 2, 6, -1, -1, -1, -1, -1, -1, -1},
00255         {5, 9, 8, 5, 8, 2, 5, 2, 6, 3, 2, 8, -1, -1, -1, -1},
00256         {2, 3, 11, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00257         {11, 0, 8, 11, 2, 0, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1},
00258         {0, 1, 9, 2, 3, 11, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1},
00259         {5, 10, 6, 1, 9, 2, 9, 11, 2, 9, 8, 11, -1, -1, -1, -1},
00260         {6, 3, 11, 6, 5, 3, 5, 1, 3, -1, -1, -1, -1, -1, -1, -1},
00261         {0, 8, 11, 0, 11, 5, 0, 5, 1, 5, 11, 6, -1, -1, -1, -1},
00262         {3, 11, 6, 0, 3, 6, 0, 6, 5, 0, 5, 9, -1, -1, -1, -1},
00263         {6, 5, 9, 6, 9, 11, 11, 9, 8, -1, -1, -1, -1, -1, -1, -1},
00264         {5, 10, 6, 4, 7, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00265         {4, 3, 0, 4, 7, 3, 6, 5, 10, -1, -1, -1, -1, -1, -1, -1},
00266         {1, 9, 0, 5, 10, 6, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1},
00267         {10, 6, 5, 1, 9, 7, 1, 7, 3, 7, 9, 4, -1, -1, -1, -1},
00268         {6, 1, 2, 6, 5, 1, 4, 7, 8, -1, -1, -1, -1, -1, -1, -1},
00269         {1, 2, 5, 5, 2, 6, 3, 0, 4, 3, 4, 7, -1, -1, -1, -1},
00270         {8, 4, 7, 9, 0, 5, 0, 6, 5, 0, 2, 6, -1, -1, -1, -1},
00271         {7, 3, 9, 7, 9, 4, 3, 2, 9, 5, 9, 6, 2, 6, 9, -1},
00272         {3, 11, 2, 7, 8, 4, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1},
00273         {5, 10, 6, 4, 7, 2, 4, 2, 0, 2, 7, 11, -1, -1, -1, -1},
00274         {0, 1, 9, 4, 7, 8, 2, 3, 11, 5, 10, 6, -1, -1, -1, -1},
00275         {9, 2, 1, 9, 11, 2, 9, 4, 11, 7, 11, 4, 5, 10, 6, -1},
00276         {8, 4, 7, 3, 11, 5, 3, 5, 1, 5, 11, 6, -1, -1, -1, -1},
00277         {5, 1, 11, 5, 11, 6, 1, 0, 11, 7, 11, 4, 0, 4, 11, -1},
00278         {0, 5, 9, 0, 6, 5, 0, 3, 6, 11, 6, 3, 8, 4, 7, -1},
00279         {6, 5, 9, 6, 9, 11, 4, 7, 9, 7, 11, 9, -1, -1, -1, -1},
00280         {10, 4, 9, 6, 4, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00281         {4, 10, 6, 4, 9, 10, 0, 8, 3, -1, -1, -1, -1, -1, -1, -1},
00282         {10, 0, 1, 10, 6, 0, 6, 4, 0, -1, -1, -1, -1, -1, -1, -1},
00283         {8, 3, 1, 8, 1, 6, 8, 6, 4, 6, 1, 10, -1, -1, -1, -1},
00284         {1, 4, 9, 1, 2, 4, 2, 6, 4, -1, -1, -1, -1, -1, -1, -1},
00285         {3, 0, 8, 1, 2, 9, 2, 4, 9, 2, 6, 4, -1, -1, -1, -1},
00286         {0, 2, 4, 4, 2, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00287         {8, 3, 2, 8, 2, 4, 4, 2, 6, -1, -1, -1, -1, -1, -1, -1},
00288         {10, 4, 9, 10, 6, 4, 11, 2, 3, -1, -1, -1, -1, -1, -1, -1},
00289         {0, 8, 2, 2, 8, 11, 4, 9, 10, 4, 10, 6, -1, -1, -1, -1},
00290         {3, 11, 2, 0, 1, 6, 0, 6, 4, 6, 1, 10, -1, -1, -1, -1},
00291         {6, 4, 1, 6, 1, 10, 4, 8, 1, 2, 1, 11, 8, 11, 1, -1},
00292         {9, 6, 4, 9, 3, 6, 9, 1, 3, 11, 6, 3, -1, -1, -1, -1},
00293         {8, 11, 1, 8, 1, 0, 11, 6, 1, 9, 1, 4, 6, 4, 1, -1},
00294         {3, 11, 6, 3, 6, 0, 0, 6, 4, -1, -1, -1, -1, -1, -1, -1},
00295         {6, 4, 8, 11, 6, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00296         {7, 10, 6, 7, 8, 10, 8, 9, 10, -1, -1, -1, -1, -1, -1, -1},
00297         {0, 7, 3, 0, 10, 7, 0, 9, 10, 6, 7, 10, -1, -1, -1, -1},
00298         {10, 6, 7, 1, 10, 7, 1, 7, 8, 1, 8, 0, -1, -1, -1, -1},
00299         {10, 6, 7, 10, 7, 1, 1, 7, 3, -1, -1, -1, -1, -1, -1, -1},
00300         {1, 2, 6, 1, 6, 8, 1, 8, 9, 8, 6, 7, -1, -1, -1, -1},
00301         {2, 6, 9, 2, 9, 1, 6, 7, 9, 0, 9, 3, 7, 3, 9, -1},
00302         {7, 8, 0, 7, 0, 6, 6, 0, 2, -1, -1, -1, -1, -1, -1, -1},
00303         {7, 3, 2, 6, 7, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00304         {2, 3, 11, 10, 6, 8, 10, 8, 9, 8, 6, 7, -1, -1, -1, -1},
00305         {2, 0, 7, 2, 7, 11, 0, 9, 7, 6, 7, 10, 9, 10, 7, -1},
00306         {1, 8, 0, 1, 7, 8, 1, 10, 7, 6, 7, 10, 2, 3, 11, -1},
00307         {11, 2, 1, 11, 1, 7, 10, 6, 1, 6, 7, 1, -1, -1, -1, -1},
00308         {8, 9, 6, 8, 6, 7, 9, 1, 6, 11, 6, 3, 1, 3, 6, -1},
00309         {0, 9, 1, 11, 6, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00310         {7, 8, 0, 7, 0, 6, 3, 11, 0, 11, 6, 0, -1, -1, -1, -1},
00311         {7, 11, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00312         {7, 6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00313         {3, 0, 8, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00314         {0, 1, 9, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00315         {8, 1, 9, 8, 3, 1, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1},
00316         {10, 1, 2, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00317         {1, 2, 10, 3, 0, 8, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1},
00318         {2, 9, 0, 2, 10, 9, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1},
00319         {6, 11, 7, 2, 10, 3, 10, 8, 3, 10, 9, 8, -1, -1, -1, -1},
00320         {7, 2, 3, 6, 2, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00321         {7, 0, 8, 7, 6, 0, 6, 2, 0, -1, -1, -1, -1, -1, -1, -1},
00322         {2, 7, 6, 2, 3, 7, 0, 1, 9, -1, -1, -1, -1, -1, -1, -1},
00323         {1, 6, 2, 1, 8, 6, 1, 9, 8, 8, 7, 6, -1, -1, -1, -1},
00324         {10, 7, 6, 10, 1, 7, 1, 3, 7, -1, -1, -1, -1, -1, -1, -1},
00325         {10, 7, 6, 1, 7, 10, 1, 8, 7, 1, 0, 8, -1, -1, -1, -1},
00326         {0, 3, 7, 0, 7, 10, 0, 10, 9, 6, 10, 7, -1, -1, -1, -1},
00327         {7, 6, 10, 7, 10, 8, 8, 10, 9, -1, -1, -1, -1, -1, -1, -1},
00328         {6, 8, 4, 11, 8, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00329         {3, 6, 11, 3, 0, 6, 0, 4, 6, -1, -1, -1, -1, -1, -1, -1},
00330         {8, 6, 11, 8, 4, 6, 9, 0, 1, -1, -1, -1, -1, -1, -1, -1},
00331         {9, 4, 6, 9, 6, 3, 9, 3, 1, 11, 3, 6, -1, -1, -1, -1},
00332         {6, 8, 4, 6, 11, 8, 2, 10, 1, -1, -1, -1, -1, -1, -1, -1},
00333         {1, 2, 10, 3, 0, 11, 0, 6, 11, 0, 4, 6, -1, -1, -1, -1},
00334         {4, 11, 8, 4, 6, 11, 0, 2, 9, 2, 10, 9, -1, -1, -1, -1},
00335         {10, 9, 3, 10, 3, 2, 9, 4, 3, 11, 3, 6, 4, 6, 3, -1},
00336         {8, 2, 3, 8, 4, 2, 4, 6, 2, -1, -1, -1, -1, -1, -1, -1},
00337         {0, 4, 2, 4, 6, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00338         {1, 9, 0, 2, 3, 4, 2, 4, 6, 4, 3, 8, -1, -1, -1, -1},
00339         {1, 9, 4, 1, 4, 2, 2, 4, 6, -1, -1, -1, -1, -1, -1, -1},
00340         {8, 1, 3, 8, 6, 1, 8, 4, 6, 6, 10, 1, -1, -1, -1, -1},
00341         {10, 1, 0, 10, 0, 6, 6, 0, 4, -1, -1, -1, -1, -1, -1, -1},
00342         {4, 6, 3, 4, 3, 8, 6, 10, 3, 0, 3, 9, 10, 9, 3, -1},
00343         {10, 9, 4, 6, 10, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00344         {4, 9, 5, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00345         {0, 8, 3, 4, 9, 5, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1},
00346         {5, 0, 1, 5, 4, 0, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1},
00347         {11, 7, 6, 8, 3, 4, 3, 5, 4, 3, 1, 5, -1, -1, -1, -1},
00348         {9, 5, 4, 10, 1, 2, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1},
00349         {6, 11, 7, 1, 2, 10, 0, 8, 3, 4, 9, 5, -1, -1, -1, -1},
00350         {7, 6, 11, 5, 4, 10, 4, 2, 10, 4, 0, 2, -1, -1, -1, -1},
00351         {3, 4, 8, 3, 5, 4, 3, 2, 5, 10, 5, 2, 11, 7, 6, -1},
00352         {7, 2, 3, 7, 6, 2, 5, 4, 9, -1, -1, -1, -1, -1, -1, -1},
00353         {9, 5, 4, 0, 8, 6, 0, 6, 2, 6, 8, 7, -1, -1, -1, -1},
00354         {3, 6, 2, 3, 7, 6, 1, 5, 0, 5, 4, 0, -1, -1, -1, -1},
00355         {6, 2, 8, 6, 8, 7, 2, 1, 8, 4, 8, 5, 1, 5, 8, -1},
00356         {9, 5, 4, 10, 1, 6, 1, 7, 6, 1, 3, 7, -1, -1, -1, -1},
00357         {1, 6, 10, 1, 7, 6, 1, 0, 7, 8, 7, 0, 9, 5, 4, -1},
00358         {4, 0, 10, 4, 10, 5, 0, 3, 10, 6, 10, 7, 3, 7, 10, -1},
00359         {7, 6, 10, 7, 10, 8, 5, 4, 10, 4, 8, 10, -1, -1, -1, -1},
00360         {6, 9, 5, 6, 11, 9, 11, 8, 9, -1, -1, -1, -1, -1, -1, -1},
00361         {3, 6, 11, 0, 6, 3, 0, 5, 6, 0, 9, 5, -1, -1, -1, -1},
00362         {0, 11, 8, 0, 5, 11, 0, 1, 5, 5, 6, 11, -1, -1, -1, -1},
00363         {6, 11, 3, 6, 3, 5, 5, 3, 1, -1, -1, -1, -1, -1, -1, -1},
00364         {1, 2, 10, 9, 5, 11, 9, 11, 8, 11, 5, 6, -1, -1, -1, -1},
00365         {0, 11, 3, 0, 6, 11, 0, 9, 6, 5, 6, 9, 1, 2, 10, -1},
00366         {11, 8, 5, 11, 5, 6, 8, 0, 5, 10, 5, 2, 0, 2, 5, -1},
00367         {6, 11, 3, 6, 3, 5, 2, 10, 3, 10, 5, 3, -1, -1, -1, -1},
00368         {5, 8, 9, 5, 2, 8, 5, 6, 2, 3, 8, 2, -1, -1, -1, -1},
00369         {9, 5, 6, 9, 6, 0, 0, 6, 2, -1, -1, -1, -1, -1, -1, -1},
00370         {1, 5, 8, 1, 8, 0, 5, 6, 8, 3, 8, 2, 6, 2, 8, -1},
00371         {1, 5, 6, 2, 1, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00372         {1, 3, 6, 1, 6, 10, 3, 8, 6, 5, 6, 9, 8, 9, 6, -1},
00373         {10, 1, 0, 10, 0, 6, 9, 5, 0, 5, 6, 0, -1, -1, -1, -1},
00374         {0, 3, 8, 5, 6, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00375         {10, 5, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00376         {11, 5, 10, 7, 5, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00377         {11, 5, 10, 11, 7, 5, 8, 3, 0, -1, -1, -1, -1, -1, -1, -1},
00378         {5, 11, 7, 5, 10, 11, 1, 9, 0, -1, -1, -1, -1, -1, -1, -1},
00379         {10, 7, 5, 10, 11, 7, 9, 8, 1, 8, 3, 1, -1, -1, -1, -1},
00380         {11, 1, 2, 11, 7, 1, 7, 5, 1, -1, -1, -1, -1, -1, -1, -1},
00381         {0, 8, 3, 1, 2, 7, 1, 7, 5, 7, 2, 11, -1, -1, -1, -1},
00382         {9, 7, 5, 9, 2, 7, 9, 0, 2, 2, 11, 7, -1, -1, -1, -1},
00383         {7, 5, 2, 7, 2, 11, 5, 9, 2, 3, 2, 8, 9, 8, 2, -1},
00384         {2, 5, 10, 2, 3, 5, 3, 7, 5, -1, -1, -1, -1, -1, -1, -1},
00385         {8, 2, 0, 8, 5, 2, 8, 7, 5, 10, 2, 5, -1, -1, -1, -1},
00386         {9, 0, 1, 5, 10, 3, 5, 3, 7, 3, 10, 2, -1, -1, -1, -1},
00387         {9, 8, 2, 9, 2, 1, 8, 7, 2, 10, 2, 5, 7, 5, 2, -1},
00388         {1, 3, 5, 3, 7, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00389         {0, 8, 7, 0, 7, 1, 1, 7, 5, -1, -1, -1, -1, -1, -1, -1},
00390         {9, 0, 3, 9, 3, 5, 5, 3, 7, -1, -1, -1, -1, -1, -1, -1},
00391         {9, 8, 7, 5, 9, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00392         {5, 8, 4, 5, 10, 8, 10, 11, 8, -1, -1, -1, -1, -1, -1, -1},
00393         {5, 0, 4, 5, 11, 0, 5, 10, 11, 11, 3, 0, -1, -1, -1, -1},
00394         {0, 1, 9, 8, 4, 10, 8, 10, 11, 10, 4, 5, -1, -1, -1, -1},
00395         {10, 11, 4, 10, 4, 5, 11, 3, 4, 9, 4, 1, 3, 1, 4, -1},
00396         {2, 5, 1, 2, 8, 5, 2, 11, 8, 4, 5, 8, -1, -1, -1, -1},
00397         {0, 4, 11, 0, 11, 3, 4, 5, 11, 2, 11, 1, 5, 1, 11, -1},
00398         {0, 2, 5, 0, 5, 9, 2, 11, 5, 4, 5, 8, 11, 8, 5, -1},
00399         {9, 4, 5, 2, 11, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00400         {2, 5, 10, 3, 5, 2, 3, 4, 5, 3, 8, 4, -1, -1, -1, -1},
00401         {5, 10, 2, 5, 2, 4, 4, 2, 0, -1, -1, -1, -1, -1, -1, -1},
00402         {3, 10, 2, 3, 5, 10, 3, 8, 5, 4, 5, 8, 0, 1, 9, -1},
00403         {5, 10, 2, 5, 2, 4, 1, 9, 2, 9, 4, 2, -1, -1, -1, -1},
00404         {8, 4, 5, 8, 5, 3, 3, 5, 1, -1, -1, -1, -1, -1, -1, -1},
00405         {0, 4, 5, 1, 0, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00406         {8, 4, 5, 8, 5, 3, 9, 0, 5, 0, 3, 5, -1, -1, -1, -1},
00407         {9, 4, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00408         {4, 11, 7, 4, 9, 11, 9, 10, 11, -1, -1, -1, -1, -1, -1, -1},
00409         {0, 8, 3, 4, 9, 7, 9, 11, 7, 9, 10, 11, -1, -1, -1, -1},
00410         {1, 10, 11, 1, 11, 4, 1, 4, 0, 7, 4, 11, -1, -1, -1, -1},
00411         {3, 1, 4, 3, 4, 8, 1, 10, 4, 7, 4, 11, 10, 11, 4, -1},
00412         {4, 11, 7, 9, 11, 4, 9, 2, 11, 9, 1, 2, -1, -1, -1, -1},
00413         {9, 7, 4, 9, 11, 7, 9, 1, 11, 2, 11, 1, 0, 8, 3, -1},
00414         {11, 7, 4, 11, 4, 2, 2, 4, 0, -1, -1, -1, -1, -1, -1, -1},
00415         {11, 7, 4, 11, 4, 2, 8, 3, 4, 3, 2, 4, -1, -1, -1, -1},
00416         {2, 9, 10, 2, 7, 9, 2, 3, 7, 7, 4, 9, -1, -1, -1, -1},
00417         {9, 10, 7, 9, 7, 4, 10, 2, 7, 8, 7, 0, 2, 0, 7, -1},
00418         {3, 7, 10, 3, 10, 2, 7, 4, 10, 1, 10, 0, 4, 0, 10, -1},
00419         {1, 10, 2, 8, 7, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00420         {4, 9, 1, 4, 1, 7, 7, 1, 3, -1, -1, -1, -1, -1, -1, -1},
00421         {4, 9, 1, 4, 1, 7, 0, 8, 1, 8, 7, 1, -1, -1, -1, -1},
00422         {4, 0, 3, 7, 4, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00423         {4, 8, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00424         {9, 10, 8, 10, 11, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00425         {3, 0, 9, 3, 9, 11, 11, 9, 10, -1, -1, -1, -1, -1, -1, -1},
00426         {0, 1, 10, 0, 10, 8, 8, 10, 11, -1, -1, -1, -1, -1, -1, -1},
00427         {3, 1, 10, 11, 3, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00428         {1, 2, 11, 1, 11, 9, 9, 11, 8, -1, -1, -1, -1, -1, -1, -1},
00429         {3, 0, 9, 3, 9, 11, 1, 2, 9, 2, 11, 9, -1, -1, -1, -1},
00430         {0, 2, 11, 8, 0, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00431         {3, 2, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00432         {2, 3, 8, 2, 8, 10, 10, 8, 9, -1, -1, -1, -1, -1, -1, -1},
00433         {9, 10, 2, 0, 9, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00434         {2, 3, 8, 2, 8, 10, 0, 1, 8, 1, 10, 8, -1, -1, -1, -1},
00435         {1, 10, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00436         {1, 3, 8, 9, 1, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00437         {0, 9, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00438         {0, 3, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
00439         {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}
00440 };
00441 
00442 ColorRGBGenerator::ColorRGBGenerator()
00443         : rgbmode(0), originx(0), originy(0), originz(0), inner(0.0), outer(0.0), minimum(0.0), maximum(0.0), needtorecolor(1), colormap(0), em_data(0)
00444 {
00445 }
00446 
00447 ColorRGBGenerator::ColorRGBGenerator(EMData* data)
00448         : rgbmode(0), minimum(0.0), maximum(0.0), needtorecolor(1), colormap(0)
00449 {
00450         set_data(data);
00451 }
00452 
00453 void ColorRGBGenerator::set_data(EMData* data)
00454 {
00455         em_data = data;
00456         originx = data->get_xsize()/2;
00457         originy = data->get_ysize()/2;
00458         originz = data->get_zsize()/2;
00459         inner = 0;
00460         outer = (float)originx;
00461 }
00462 
00463 void ColorRGBGenerator::set_cmap_data(EMData* data)
00464 {
00465         cmap = data;
00466         setMinMax(data->get_attr("minimum"), data->get_attr("maximum"));
00467 }
00468         
00469 void ColorRGBGenerator::generateRadialColorMap()
00470 {
00471 
00472         int size = int(sqrt((float)originx*originx + originy*originy + originz*originz));
00473         
00474         if(colormap) delete colormap;
00475         colormap = new float[size*3];
00476         
00477         for(int i = 0; i < size; i++){
00478                 float normrad = 4.189f*(i - inner)/(outer - inner);
00479                 if(normrad < 2.094){
00480                         if (normrad < 0.0) normrad = 0.0;
00481                         colormap[i*3] = 0.5f*(1 + cos(normrad)/cos(1.047f - normrad));
00482                         colormap[i*3 + 1] = 1.5f - colormap[i*3];
00483                         colormap[i*3 + 2] = 0.0;
00484                 }
00485                 if(normrad >= 2.094){
00486                         if (normrad > 4.189f) normrad = 4.189f;
00487                         normrad =- 2.094f;
00488                         colormap[i*3] = 0.0;
00489                         colormap[i*3 + 1] = 0.5f*(1 + cos(normrad)/cos(1.047f - normrad));
00490                         colormap[i*3 + 2] = 1.5f - colormap[i*3 + 1];
00491                 }
00492         }
00493 }
00494 
00495 float* ColorRGBGenerator::getRGBColor(int x, int y, int z)
00496 {
00497         // Get data using radial color info
00498         if (rgbmode == 1){
00499                 //calculate radius
00500 #ifdef _WIN32
00501                 float rad = (float)sqrt(pow(double(x-originx),2) + pow(double(y-originy),2) + pow(double(z-originz),2));
00502 #else
00503                 float rad = sqrt(float(pow(x-originx,2) + pow(y-originy,2) + pow(z-originz,2)));
00504 #endif  //_WIN32
00505         
00506                 //This indicates that new color info needs to be bound (to the GPU)
00507                 if(needtorecolor){ 
00508                         generateRadialColorMap();
00509                         needtorecolor = false;
00510                 }
00511                 
00512                 return &colormap[int(rad)*3];
00513         }
00514         // Get data using color map
00515         if (rgbmode == 2){
00516                 float value = cmap->get_value_at(x, y, z);
00517                 value = 4.189f*(value - minimum)/(maximum - minimum);
00518                 if(value < 2.094){
00519                         if (value < 0.0) value = 0.0;
00520                         rgb[0] = 0.5f*(1 + cos(value)/cos(1.047f - value));
00521                         rgb[1] = 1.5f - rgb[0];
00522                         rgb[2] = 0.0;
00523                 }
00524                 if(value >= 2.094){
00525                         if (value > 4.189f) value = 4.189f;
00526                         value =- 2.094f;
00527                         rgb[0] = 0.0;
00528                         rgb[1] = 0.5f*(1 + cos(value)/cos(1.047f - value));
00529                         rgb[2] = 1.5f - rgb[1];
00530                 }
00531                 
00532                 return &rgb[0];
00533         }
00534         
00535         return &colormap[0];    
00536 }
00537 
00538 MarchingCubes::MarchingCubes()
00539         : _isodl(0), needtobind(1)
00540 {
00541         rgbgenerator = ColorRGBGenerator();
00542 
00543 #ifdef _WIN32
00544         typedef void (APIENTRYP PFNGLGENBUFFERSPROC) (GLsizei n, GLuint *buffers);
00545         PFNGLGENBUFFERSPROC glGenBuffers;
00546         glGenBuffers = (PFNGLGENBUFFERSPROC) wglGetProcAddress("glGenBuffers");
00547 #endif  //_WIN32
00548 
00549         glGenBuffers(4, buffer);
00550 
00551 }
00552 
00553 MarchingCubes::MarchingCubes(EMData * em)
00554         : _isodl(0), needtobind(1)
00555 {
00556         rgbgenerator = ColorRGBGenerator();
00557 
00558 #ifdef _WIN32
00559         typedef void (APIENTRYP PFNGLGENBUFFERSPROC) (GLsizei n, GLuint *buffers);
00560         PFNGLGENBUFFERSPROC glGenBuffers;
00561         glGenBuffers = (PFNGLGENBUFFERSPROC) wglGetProcAddress("glGenBuffers");
00562 #endif  //_WIN32
00563 
00564         glGenBuffers(4, buffer);
00565         set_data(em);
00566 }
00567 
00568 void MarchingCubes::clear_min_max_vals()
00569 {
00570         for (vector<EMData*>::iterator it = minvals.begin(); it != minvals.end(); ++it)
00571         {
00572                 if ( (*it) != 0 ) delete *it;
00573         }
00574         minvals.clear();
00575 
00576         for (vector<EMData*>::iterator it = maxvals.begin(); it != maxvals.end(); ++it)
00577         {
00578                 if ( (*it) != 0 ) delete *it;
00579         }
00580         maxvals.clear();
00581 }
00582 
00583 bool MarchingCubes::calculate_min_max_vals()
00584 {
00585 
00586         if (_emdata == NULL ) throw NullPointerException("Error, cannot generate search tree if the overriding EMData object is NULL");
00587 
00588         clear_min_max_vals();
00589 
00590         int nx = _emdata->get_xsize();
00591         int ny = _emdata->get_ysize();
00592         int nz = _emdata->get_zsize();
00593 
00594         // Create the binary min max tree
00595         while ( nx > 1 || ny > 1 || nz > 1 )
00596         {
00597                 int size = minvals.size();
00598 
00599                 if ( size == 0 ){
00600                         Dict a;
00601                         // Setting search to 3 at the bottom level is most important.
00602                         // FIXME - d.woolford add comments heere
00603                         a["search"] = 3;
00604                         minvals.push_back(_emdata->process("math.minshrink",a));
00605                         maxvals.push_back(_emdata->process("math.maxshrink",a));
00606                 }else {
00607                         minvals.push_back(minvals[size-1]->process("math.minshrink"));
00608                         maxvals.push_back(maxvals[size-1]->process("math.maxshrink"));
00609                 }
00610 
00611                 nx = minvals[size]->get_xsize();
00612                 ny = minvals[size]->get_ysize();
00613                 nz = minvals[size]->get_zsize();
00614 #if MARCHING_CUBES_DEBUG
00615                 cout << "dims are " << nx << " " << ny << " " << nz << endl;
00616 #endif
00617         }
00618 
00619         drawing_level = -1;
00620 
00621         return true;
00622 }
00623 
00624 MarchingCubes::~MarchingCubes() {
00625         clear_min_max_vals();
00626 
00627 #ifdef _WIN32
00628         typedef void (APIENTRYP PFNGLDELETEBUFFERSPROC) (GLsizei n, const GLuint *buffers);
00629         PFNGLDELETEBUFFERSPROC glDeleteBuffers;
00630         glDeleteBuffers = (PFNGLDELETEBUFFERSPROC) wglGetProcAddress("glDeleteBuffers");
00631 #endif  //_WIN32
00632 
00633         glDeleteBuffers(4, buffer);
00634 }
00635 
00636 Dict MarchingCubes::get_isosurface()
00637 {
00638         calculate_surface();
00639         Dict d;
00640         d.put("points", (float*)pp.get_data());
00641         for (unsigned int i = 0; i < ff.elem(); ++i ) ff[i] /= 3;
00642         d.put("faces", (unsigned int*)ff.get_data());
00643         d.put("normals", (float*)nn.get_data());
00644         d.put("size", ff.elem());
00645         return d;
00646 }
00647 
00648 void MarchingCubes::surface_face_z()
00649 {
00650         float* f = pp.get_data();
00651         float* n = nn.get_data();
00652         for (unsigned int i = 0; i < pp.elem(); i += 3 ) {
00653                 if (f[i+2] == 0.5) continue;
00654                 Vec3f z(0,0,1.0);
00655                 Vec3f axis(-n[i+1],n[i],0);
00656                 axis.normalize();
00657         
00658                 Dict d;
00659                 d["type"] = "spin";
00660                 d["Omega"] =  90.f;
00661                 d["n1"] = axis[0];
00662                 d["n2"] = axis[1];
00663                 d["n3"] = 0;
00664                 Transform t(d);
00665                 Vec3f delta = t*z;
00666                 
00667                 f[i] += delta[0]*.25f;
00668                 f[i+1] += delta[1]*.25f;
00669                 f[i+2] = 0.5;
00670         }
00671         
00672         for (unsigned int i = 0; i < nn.elem(); i += 3 ) {
00673                 n[i] = 0;
00674                 n[i+1] = 0;
00675                 n[i+2] = 1;
00676         }
00677 }
00678 
00679 void MarchingCubes::set_data(EMData* data)
00680 {
00681         if ( data->get_zsize() == 1 ) throw ImageDimensionException("The z dimension of the image must be greater than 1");
00682         _emdata = data;
00683         calculate_min_max_vals();
00684         rgbgenerator.set_data(data);
00685 }
00686 
00687 void MarchingCubes::set_surface_value(const float value) {
00688 
00689         if(_surf_value == value) return;
00690 
00691         _surf_value = value;
00692 
00693 }
00694 
00695 void MarchingCubes::calculate_surface() {
00696 
00697         if ( _emdata == 0 ) throw NullPointerException("Error, attempt to generate isosurface, but the emdata image object has not been set");
00698         if ( minvals.size() == 0 || maxvals.size() == 0 ) throw NotExistingObjectException("Vector of EMData pointers", "Error, the min and max val search trees have not been created");
00699         
00700         point_map.clear();
00701         pp.clear();
00702         nn.clear();
00703         ff.clear();
00704         vv.clear();
00705 
00706 #if MARCHING_CUBES_DEBUG
00707         int time0 = clock();
00708 #endif
00709 
00710         float min = minvals[minvals.size()-1]->get_value_at(0,0,0);
00711         float max = maxvals[minvals.size()-1]->get_value_at(0,0,0);
00712         if ( min < _surf_value &&  max > _surf_value) draw_cube(0,0,0,minvals.size()-1);
00713         
00714 #if MARCHING_CUBES_DEBUG
00715         int time1 = clock();
00716         cout << "It took " << (time1-time0) << " " << (float)(time1-time0)/CLOCKS_PER_SEC << " to traverse the search tree and generate polygons" << endl;
00717         cout << "... using surface value " << _surf_value << endl;
00718 #endif
00719 }
00720 
00721 void MarchingCubes::draw_cube(const int x, const int y, const int z, const int cur_level ) {
00722 
00723         if ( cur_level == drawing_level )
00724         {
00725                 EMData* e;
00726                 if ( drawing_level == - 1 ) e = _emdata;
00727                 else e = minvals[drawing_level];
00728                 if ( x < (e->get_xsize()-1) && y < (e->get_ysize()-1) && z < (e->get_zsize()-1))
00729                         marching_cube(x,y,z, cur_level);
00730         }
00731         else
00732         {
00733                 EMData* e;
00734                 if ( cur_level > 0 ) {
00735                         int xsize = minvals[cur_level-1]->get_xsize();
00736                         int ysize = minvals[cur_level-1]->get_ysize();
00737                         int zsize = minvals[cur_level-1]->get_zsize();
00738                         e = minvals[cur_level-1];
00739                         for(int i=0; i<8; ++i ) {
00740                                 int xx = 2*x+a2fVertexOffset[i][0];
00741                                 if ( xx >= xsize ) continue;
00742                                 int yy = 2*y+a2fVertexOffset[i][1];
00743                                 if ( yy >= ysize ) continue;
00744                                 int zz = 2*z+a2fVertexOffset[i][2];
00745                                 if ( zz >= zsize ) continue;
00746 
00747                                 float min = minvals[cur_level-1]->get_value_at(xx,yy,zz);
00748                                 float max = maxvals[cur_level-1]->get_value_at(xx,yy,zz);
00749                                 if ( min < _surf_value &&  max > _surf_value)
00750                                         draw_cube(xx,yy,zz,cur_level-1);
00751                         }
00752                 }
00753                 else {
00754                         e = _emdata;
00755                         for(int i=0; i<8; ++i ) {
00756                                         draw_cube(2*x+a2fVertexOffset[i][0],2*y+a2fVertexOffset[i][1],2*z+a2fVertexOffset[i][2],cur_level-1);
00757                         }
00758                 }
00759 
00760                 if ( x == (minvals[cur_level]->get_xsize()-1) ) {
00761                         if ( e->get_xsize() > 2*x ){
00762                                 for(int i=0; i<4; ++i ) {
00763                                         draw_cube(2*x+a2fPosXOffset[i][0],2*y+a2fPosXOffset[i][1],2*z+a2fPosXOffset[i][2],cur_level-1);
00764                                 }
00765                         }
00766                         if ( y == (minvals[cur_level]->get_ysize()-1) ) {
00767                                 if ( e->get_ysize() > 2*y ) {
00768                                         for(int i=0; i<2; ++i ) {
00769                                                 draw_cube(2*x+a2fPosXPosYOffset[i][0],2*y+a2fPosXPosYOffset[i][1],2*z+a2fPosXPosYOffset[i][2],cur_level-1);
00770                                         }
00771                                 }
00772                                 if (  z == (minvals[cur_level]->get_zsize()-1) ){
00773                                         if ( e->get_zsize() > 2*z ) {
00774                                                 draw_cube(2*x+2,2*y+2,2*z+2,cur_level-1);
00775                                         }
00776                                 }
00777                         }
00778                         if ( z == (minvals[cur_level]->get_zsize()-1) ) {
00779                                 if ( e->get_zsize() > 2*z ) {
00780                                         for(int i=0; i<2; ++i ) {
00781                                                 draw_cube(2*x+a2fPosXPosZOffset[i][0],2*y+a2fPosXPosZOffset[i][1],2*z+a2fPosXPosZOffset[i][2],cur_level-1);
00782                                         }
00783                                 }
00784                         }
00785                 }
00786                 if ( y == (minvals[cur_level]->get_ysize()-1) ) {
00787                         if ( e->get_ysize() > 2*y ) {
00788                                 for(int i=0; i<4; ++i ) {
00789                                         draw_cube(2*x+a2fPosYOffset[i][0],2*y+a2fPosYOffset[i][1],2*z+a2fPosYOffset[i][2],cur_level-1);
00790                                 }
00791                         }
00792                         if ( z == (minvals[cur_level]->get_ysize()-1) ) {
00793                                 if ( e->get_zsize() > 2*z ) {
00794                                         for(int i=0; i<2; ++i ) {
00795                                                 draw_cube(2*x+a2fPosYPosZOffset[i][0],2*y+a2fPosYPosZOffset[i][1],2*z+a2fPosYPosZOffset[i][2],cur_level-1);
00796                                         }
00797                                 }
00798                         }
00799                 }
00800                 if ( z == (minvals[cur_level]->get_zsize()-1) ) {
00801                         if ( e->get_zsize() > 2*z ) {
00802                                 for(int i=0; i<4; ++i ) {
00803                                         draw_cube(2*x+a2fPosZOffset[i][0],2*y+a2fPosZOffset[i][1],2*z+a2fPosZOffset[i][2],cur_level-1);
00804                                 }
00805                         }
00806                 }
00807 
00808         }
00809 }
00810 
00811 void MarchingCubes::get_normal(Vector3 &normal, int fX, int fY, int fZ)
00812 {
00813     normal[0] = _emdata->get_value_at(fX-1, fY, fZ) - _emdata->get_value_at(fX+1, fY, fZ);
00814     normal[1] = _emdata->get_value_at(fX, fY-1, fZ) - _emdata->get_value_at(fX, fY+1, fZ);
00815     normal[2] = _emdata->get_value_at(fX, fY, fZ-1) - _emdata->get_value_at(fX, fY, fZ+1);
00816     normal.normalize();
00817 }
00818 
00819 float MarchingCubes::get_offset(float fValue1, float fValue2, float fValueDesired)
00820 {
00821         float fDelta = fValue2 - fValue1;
00822 
00823         if(fDelta == 0.0f)
00824         {
00825                 return 0.5f;
00826         }
00827         return (fValueDesired - fValue1)/fDelta;
00828 }
00829 
00830 int MarchingCubes::get_edge_num(int x, int y, int z, int edge) {
00831         // edge direction is right, down, back (x, y, z)
00832         unsigned int index = 0;
00833         index = (x << 22) | (y << 12) | (z << 2) | edge;
00834         return index;
00835 }
00836 
00837 void MarchingCubes::color_vertices()
00838 {
00839         cc.clear();
00840 #ifdef _WIN32
00841         int scaling = (int)pow(2.0,drawing_level + 1);          // Needed to account for sampling rate
00842 #else
00843         int scaling = pow(2,drawing_level + 1);         // Needed to account for sampling rate
00844 #endif  //_WIN32        
00845         //Color vertices. We don't need to rerun marching cubes on color vertices, so this method improves effciency
00846         for(unsigned int i = 0; i < vv.elem(); i+=3){
00847                 float* color = rgbgenerator.getRGBColor(scaling*vv[i], scaling*vv[i+1], scaling*vv[i+2]);
00848                 cc.push_back_3(color);
00849         }
00850         rgbgenerator.setNeedToRecolor(false);
00851 }
00852                 
00853 void MarchingCubes::marching_cube(int fX, int fY, int fZ, int cur_level)
00854 {
00855 //      extern int aiCubeEdgeFlags[256];
00856 //      extern int a2iTriangleConnectionTable[256][16];
00857 
00858         int iCorner, iVertex, iVertexTest, iEdge, iTriangle, iFlagIndex, iEdgeFlags;
00859         float fOffset;
00860         Vector3 sColor;
00861         float afCubeValue[8];
00862         float asEdgeVertex[12][3];
00863         int pointIndex[12];
00864 
00865         int fxScale = 1, fyScale = 1, fzScale = 1;
00866         if ( cur_level != -1 )
00867         {
00868                 fxScale = _emdata->get_xsize()/minvals[cur_level]->get_xsize();
00869                 fyScale = _emdata->get_ysize()/minvals[cur_level]->get_ysize();
00870                 fzScale = _emdata->get_zsize()/minvals[cur_level]->get_zsize();
00871                 for(iVertex = 0; iVertex < 8; iVertex++)
00872                 {
00873                         afCubeValue[iVertex] = _emdata->get_value_at( fxScale*(fX + a2fVertexOffset[iVertex][0]),
00874                                                 fyScale*(fY + a2fVertexOffset[iVertex][1]), fzScale*(fZ + a2fVertexOffset[iVertex][2]));
00875                 }
00876         }
00877         else
00878         {
00879                 //Make a local copy of the values at the cube's corners
00880                 for(iVertex = 0; iVertex < 8; iVertex++)
00881                 {
00882                         afCubeValue[iVertex] = _emdata->get_value_at( fX + a2fVertexOffset[iVertex][0],
00883                                         fY + a2fVertexOffset[iVertex][1], fZ + a2fVertexOffset[iVertex][2]);
00884                 }
00885         }
00886 
00887         //Find which vertices are inside of the surface and which are outside
00888         iFlagIndex = 0;
00889         for(iVertexTest = 0; iVertexTest < 8; iVertexTest++)
00890         {
00891                 if (_surf_value >= 0 ){
00892                         if(afCubeValue[iVertexTest] <= _surf_value)
00893                                 iFlagIndex |= 1<<iVertexTest;
00894                 }
00895                 else {
00896                         if(afCubeValue[iVertexTest] >= _surf_value)
00897                                 iFlagIndex |= 1<<iVertexTest;
00898                 }
00899         }
00900 
00901         //Find which edges are intersected by the surface
00902         iEdgeFlags = aiCubeEdgeFlags[iFlagIndex];
00903 
00904         //If the cube is entirely inside or outside of the surface, then there will be no intersections
00905         if(iEdgeFlags == 0) return;
00906 
00907         //Find the point of intersection of the surface with each edge
00908         //Then find the normal to the surface at those points
00909         for(iEdge = 0; iEdge < 12; iEdge++)
00910         {
00911                 //if there is an intersection on this edge
00912                 if(iEdgeFlags & (1<<iEdge))
00913                 {
00914                         fOffset = get_offset(afCubeValue[ a2iEdgeConnection[iEdge][0] ],
00915                                                                  afCubeValue[ a2iEdgeConnection[iEdge][1] ], _surf_value);
00916 
00917                         if ( cur_level == -1 ){
00918                                 asEdgeVertex[iEdge][0] = fX + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][0]  +  fOffset * a2fEdgeDirection[iEdge][0]) + 0.5f;
00919                                 asEdgeVertex[iEdge][1] = fY + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][1]  +  fOffset * a2fEdgeDirection[iEdge][1]) + 0.5f;
00920                                 asEdgeVertex[iEdge][2] = fZ + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][2]  +  fOffset * a2fEdgeDirection[iEdge][2]) + 0.5f;
00921                         } else {
00922                                 asEdgeVertex[iEdge][0] = fxScale*(fX + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][0]  +  fOffset * a2fEdgeDirection[iEdge][0])) + 0.5f;
00923                                 asEdgeVertex[iEdge][1] = fyScale*(fY + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][1]  +  fOffset * a2fEdgeDirection[iEdge][1])) + 0.5f;
00924                                 asEdgeVertex[iEdge][2] = fzScale*(fZ + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][2]  +  fOffset * a2fEdgeDirection[iEdge][2])) + 0.5f;
00925                         }
00926 
00927                         pointIndex[iEdge] = get_edge_num(fX+edgeLookUp[iEdge][0], fY+edgeLookUp[iEdge][1], fZ+edgeLookUp[iEdge][2], edgeLookUp[iEdge][3]);
00928                 }
00929         }
00930         
00931         //Save voxel coords for later color processing
00932         int vox[3] = {fX, fY, fZ};
00933         
00934         //Draw the triangles that were found.  There can be up to five per cube
00935         for(iTriangle = 0; iTriangle < 5; iTriangle++)
00936         {
00937                 if(a2iTriangleConnectionTable[iFlagIndex][3*iTriangle] < 0)
00938                         break;
00939 
00940                 float pts[3][3];
00941                 for(iCorner = 0; iCorner < 3; iCorner++)
00942                 {
00943                         iVertex = a2iTriangleConnectionTable[iFlagIndex][3*iTriangle+iCorner];
00944                         memcpy(&pts[iCorner][0],  &asEdgeVertex[iVertex][0], 3*sizeof(float));
00945                 }
00946 
00947                 
00948                 
00949                 float v1[3] = {pts[1][0]-pts[0][0],pts[1][1]-pts[0][1],pts[1][2]-pts[0][2]};
00950                 float v2[3] = {pts[2][0]-pts[1][0],pts[2][1]-pts[1][1],pts[2][2]-pts[1][2]};
00951                 
00952                 float n[3] = { v1[1]*v2[2] - v1[2]*v2[1], v1[2]*v2[0] - v1[0]*v2[2], v1[0]*v2[1] - v1[1]*v2[0] };
00953                         
00954                 
00955                 for(iCorner = 0; iCorner < 3; iCorner++)
00956                 {
00957                         // Without vertex normalization
00958 //                      iVertex = a2iTriangleConnectionTable[iFlagIndex][3*iTriangle+iCorner];
00959 //                      int ss = pp.elem();
00960 //                      pp.push_back_3(&pts[iCorner][0]);
00961 //                      nn.push_back_3(&n[0]);
00962 //                      ff.push_back(ss);
00963 
00964 //                      With vertex normalization
00965                         iVertex = a2iTriangleConnectionTable[iFlagIndex][3*iTriangle+iCorner];
00966                         map<int,int>::iterator it = point_map.find(pointIndex[iVertex]);
00967                         if ( it == point_map.end() ){
00968                                 vv.push_back_3(&vox[0]);
00969                                 int ss = pp.elem();
00970                                 pp.push_back_3(&pts[iCorner][0]);
00971                                 nn.push_back_3(&n[0]);
00972                                 ff.push_back(ss);
00973                                 point_map[pointIndex[iVertex]] = ss;
00974                         } else {
00975                                 int idx = it->second;
00976                                 ff.push_back(idx);
00977                                 nn[idx] += n[0];
00978                                 nn[idx+1] += n[1];
00979                                 nn[idx+2] += n[2];
00980                         }
00981                 }
00982         }
00983 }
00984 
00985 
00986 
00987 U3DWriter::U3DWriter() : DIFFUSE_COLOR_COUNT(1),
00988                 SPECULAR_COLOR_COUNT(1)
00989 {
00990 }
00991 U3DWriter::~U3DWriter() {}
00992 
00993 using std::ofstream;
00994 
00995 int U3DWriter::write(const string& filename) {
00996         // Must open the ofstream in binary mode
00997         ofstream of(filename.c_str(), ofstream::binary);
00998         write(of);
00999         of.close();
01000 
01001         return 1;
01002 }
01003 
01004 ostream& U3DWriter::write(ostream& os) {
01005         write_header(os);
01006         return os;
01007 }
01008 
01009 unsigned int U3DWriter::size_of_in_bytes(){
01010         // this is just the size in bytes of all of the entries in this object that are written to the binary
01011         // output. This is based only on the behavior of write_header and needs to be udpated
01012         unsigned size = 4+2+2+4+4+8+4+8; // 36 bytes
01013         return size;
01014 }
01015 ostream& U3DWriter::write_header(ostream& os)
01016 {
01017         // checks
01018         test_type_sizes();
01019 
01020         U32 block_type_file_header = 0x00443355; // This is the block type tag of a file header, as taken form the ECMA spec
01021         write( os, block_type_file_header);
01022 //
01023         I16 major_version = -1; // Compliance has not been tested for this encoder, so we must specify a value less than 0
01024         I16 minor_version =  0; // This is the version of this encoder, which we are calling '0'
01025         write( os, major_version);
01026         write( os, minor_version);
01027 
01028         U32 profile_identifier = 0x00000000; // Base profile - no optional features ares used
01029         write( os, profile_identifier);
01030 
01031         U32 declaration_size = size_of_in_bytes(); // This will have to be addressed at a later point, this is the size if the declaration block in bytes
01032         write( os, declaration_size);
01033 
01034         U64 file_size = size_of_in_bytes(); // This is the size of the file in bytes
01035         write( os, file_size);
01036 
01037         U32 character_encoding = 106; // UTF-8 MIB from the IANA
01038         write( os, character_encoding);
01039 
01040         F64 unit_scaling = 1.0; // This should eventually be the correct scaling of the objects
01041         write( os, unit_scaling);
01042 
01043 
01044         return os;
01045 }
01046 
01047 void U3DWriter::test_type_sizes()
01048 {
01049         bool error = false;
01050         if (sizeof(F64) != 8 ){
01051                 cout << "Error, size of double is not 64 bytes, it's " << sizeof(F64)*4 << endl;
01052                 error = true;
01053         }
01054         if (sizeof(F32) != 4 ){
01055                 cout << "Error, size of float is not 32 bytes, it's " << sizeof(F32)*4 << endl;
01056                 error = true;
01057         }
01058         if (sizeof(U64) != 8) {
01059                 cout << "Error, size of long unsigned int is not 64 bytes, it's " << sizeof(U64)*4 << endl;
01060                 error = true;
01061         }
01062         if (sizeof(U32) != 4) {
01063                 cout << "Error, size of unsigned int is not 32 bytes, it's " << sizeof(U32)*4 << endl;
01064                 error = true;
01065         }
01066         if (sizeof(I16) != 2) {
01067                 cout << "Error, size of short int is not 16 bytes, it's " << sizeof(I16)*4 << endl;
01068                 error = true;
01069         }
01070         if (sizeof(U16) != 2) {
01071                 cout << "Error, size of short unsigned int is not 16 bytes, it's " << sizeof(U16)*4 << endl;
01072                 error = true;
01073         }
01074         if (sizeof(U8) != 1) {
01075                 cout << "Error, size of unsigned char is not  bytes, it's " << sizeof(U8)*4 << endl;
01076                 error = true;
01077         }
01078 
01079         if (error) {
01080                 throw;
01081         }
01082 }
01083 
01084 ostream& U3DWriter::write_clod_mesh_generator_node(ostream& os)
01085 {
01086         /*
01087         CLOD Mesh Declaration
01088         */
01089         U32 block_type_clod_mesh_generator = 0xFFFFFF31; // This is the block type tag of the continuous level of detail mesh generator, as taken form the ECMA spec
01090         write( os, block_type_clod_mesh_generator);
01091 
01092         string name("testmesh"); // The unique name, we get to make this up ourselves. It could be an empty string, we'd still have to call write_string(os,"")
01093         write(os,name);
01094 
01095         U32 chain_index = 0x00000000; // the value of Chain Index shall be zero for this type - as in the SPEC
01096         write( os, chain_index);
01097 
01098         /*
01099         Max Mesh Description
01100         */
01101         U32 mesh_attributes = 0x00000000; // Faces in the mesh have a normal index at each corner 0x00000001 is used to exlude normals POTENTIALLY USEFUL
01102         write(os,mesh_attributes);
01103         U32 face_count = ff.elem(); // The number of faces TO FILL IN LATER
01104         write(os,face_count);
01105         U32 position_count = pp.elem(); // The number of positions in the position array TO FILL IN LATER
01106         write(os,position_count);
01107         U32 normal_count = nn.elem(); // The number of normals in the normal array TO FILL IN LATER
01108         write(os,normal_count);
01109         U32 diffuse_color_count = DIFFUSE_COLOR_COUNT; // The number of colors in the diffuse color array TO FILL IN LATER
01110         write(os,diffuse_color_count);
01111         U32 specular_color_count = SPECULAR_COLOR_COUNT; // The number of colors in the specular color array TO FILL IN LATER
01112         write(os,specular_color_count);
01113         U32 texture_coord_count = 0x00000000; // The number of texture coordinates in teh texture coordinate array POTENTIALLY USEFUL
01114         write(os,texture_coord_count);
01115         U32 shading_count = 1; // The number of shading descriptions used in the mesh. This must correspond with the shader list in the shading group (see directly below we are using only one shader
01116         write(os,shading_count);
01117 
01118         /*
01119         Shading  Description
01120         */
01121         U32 shading_attributes = 0x00000003; // 0 means use neither diffuse or specular colors, 1 means use per vertex diffuse, 2 means use per vertex specular, 3 means use both POTENTIALLY USEFUL
01122         write(os,shading_attributes);
01123         U32 texture_layout_count = 0x00000000; // The number of texture layers used by this shader list
01124         write(os,texture_layout_count);
01125         U32 texture_coord_dimensions = 0x00000002; // The number of dimensions in the texture coordinate vector. I chose default to be 2. No particular reason. POTENTIALLY USEFUL
01126         write(os,texture_coord_dimensions);
01127         U32 original_shading_id = 0; // Original shading index for this shader list. Informative parameter only. This is shader 0
01128         write(os,original_shading_id);
01129 
01130         /*
01131         CLOD Description - describes the range of resolutions available for the continuous level of detail mesh
01132         If there were more than one level of detail than these two numbers would be different
01133         */
01134         U32 minimum_resolution = position_count; // the number of positions in the base mesh
01135         write(os,minimum_resolution);
01136         U32 final_maximum_resolution = position_count; // the number of positions in the Max Mesh Description (by definition)
01137         write(os,final_maximum_resolution);
01138 
01139         /*
01140         Resource Description
01141         */
01142 
01143         /*
01144         Quality Factors
01145         for information only. Not used by the renderer. Helpful for conveying information to the user
01146         */
01147         U32 position_quality_factor = 0x00000000; // position quality factor. Descriptive information only
01148         write(os,position_quality_factor);
01149         U32 normal_quality_factor = 0x00000000; // normal quality factor. Descriptive information only
01150         write(os,normal_quality_factor);
01151         U32 texture_coord_quality_factor = 0x00000000; // texture coordinate quality factor. Descriptive information only
01152         write(os,texture_coord_quality_factor);
01153 
01154         /*
01155         Inverse Quantization
01156         used to reconstruct floating point values that have been quantized.
01157         */
01158         F32 postion_inverse_quant = 1.0; // inverse quantization factor used in the reconstruction of the position vectors
01159         write(os,postion_inverse_quant);
01160         F32 normal_inverse_quant = 1.0; // inverse quantization factor used in the reconstruction of the normal vectors
01161         write(os,normal_inverse_quant);
01162         F32 texture_coord_inverse_quant = 1.0; // inverse quantization factor used in the reconstruction of the texture coordinates
01163         write(os,texture_coord_inverse_quant);
01164         F32 diffuse_color_inverse_quant = 1.0; // inverse quantization factor used in the reconstruction of the diffuse colors
01165         write(os,diffuse_color_inverse_quant);
01166         F32 specular_color_inverse_quant = 1.0; // inverse quantization factor used in the reconstruction of the specular colors
01167         write(os,specular_color_inverse_quant);
01168 
01169         /*
01170         Resource Parameters
01171         parameters that help to define the conversion from the Author Mesh format to the Render Mesh format
01172         */
01173 
01174         F32 normal_crease_parameter = 1.0; // A dot product value in the range -1 to 1, used to decide whether or not normals at the same position will be merged. 1.0 means never, -1.0 means always. 0 means the two normals must be separated by an angle no greater than 90 degrees to be merged. Think in angles.
01175         write(os,normal_crease_parameter);
01176         F32 normal_update_parameter = 0.0; // A strange a parameter that I couldn't make sense of - I think it means it will change the actual file itself if it encounters what it deems 'normal errors'
01177         write(os,normal_update_parameter);
01178         F32 normal_tolerance_parameter = 0.0; // Normals which are closer together than this value are considered equivalent in the Render Mesh. This is useful for compression
01179         write(os,normal_tolerance_parameter);
01180 
01181         /*
01182         Skeleton description
01183         used in bones-based animation by the Animation Modifier
01184         */
01185         U32 bone_count = 0x00000000; // The number of bones associated with this mesh. We will always have 0, but this could change (unlikely).
01186         write(os,bone_count);
01187 
01188         // WARNING - if bone_count is ever greater than one, then more writing would have to occur here
01189 
01191         /*
01192         Base mesh is basically the minimum LOD mesh - it must exist
01193         */
01194         {
01195         U32 block_type_clod_base_mesh_continuation = 0xFFFFFF3B; // This is the block type tag of the CLOD Base Mesh Continuation, as taken form the ECMA spec
01196         write( os, block_type_clod_base_mesh_continuation);
01197 
01198         write(os,name); // We use the same name as above
01199 
01200         U32 chain_index = 0x00000000; // the value of Chain Index shall be zero for this type - as in the SPEC
01201         write( os, chain_index);
01202 
01203         /*
01204         Base Mesh Description
01205         */
01206         U32     base_face_count = ff.elem(); // The number of faces in the base mesh TO FILL IN LATER
01207         write( os, base_face_count);
01208         U32     base_position_count = pp.elem(); // The number of positions used by base mesh in the position array TO FILL IN LATER
01209         write( os, base_position_count);
01210         U32 base_normal_count = nn.elem(); // The number of normals used by the base mesh in the normal array TO FILL IN LATER
01211         write( os, base_normal_count);
01212         U32 base_diffuse_color_count = DIFFUSE_COLOR_COUNT; // The number of diffuse colors used by the base mesh in the diffuse color array TO FILL IN LATER
01213         write( os, base_diffuse_color_count);
01214         U32 base_specular_color_count = SPECULAR_COLOR_COUNT; // The number of specular colors used by the base mesh in the specular color array TO FILL IN LATER
01215         write( os, base_specular_color_count);
01216         U32 base_texture_coord_count = 0x00000000; // The number of texture coordinates used by the base mesh in texture coordinate array TO FILL IN LATER
01217         write( os, base_texture_coord_count);
01218 
01219         /*
01220         Base mesh data
01221         */
01222 
01223         // Write position data
01224         F32* data = pp.get_data();
01225         for(unsigned int i = 0; i < pp.elem(); ++i) {
01226                 write(os,data[i]);
01227         }
01228 
01229         // Write normal data
01230         data = nn.get_data();
01231         for(unsigned int i = 0; i < nn.elem(); ++i) {
01232                 write(os,data[i]);
01233         }
01234 
01235         // Write Diffuse color, this is in rgba format
01236         F32 diffuse_rgba[4] = {1.0,0.0,0.0,1.0};
01237         for (unsigned int i = 0; i < 4; ++i) {
01238                 write(os,diffuse_rgba[i]);
01239         }
01240 
01241         // Write Specular color, this is in rgba format
01242         F32 specular_rgba[4] = {1.0,0.0,0.0,1.0};
01243         for (unsigned int i = 0; i < 4; ++i) {
01244                 write(os,specular_rgba[i]);
01245         }
01246 
01247         // THERE ARE NO TEXTURE COORDINATES, BUT IF THERE WERE WE WOULD HAVE TO WRITE THEM HERE
01248         // i.e. for i in range(0,base_texture_coord_count) write texture coords
01249 
01250         // Write normal data
01251         U32* faces = ff.get_data();
01252         for(unsigned int i = 0; i < pp.elem(); i = i+3) {
01253                 U32 shading_id = 0; // We only have one shader defined. This could could changed in future
01254                 write(os,shading_id);
01255 
01256                 // write base corner info - there are always three corners
01257                 for (unsigned int j =0; j < 3; ++j){
01258                         U32 position_index = faces[i+j];
01259                         write(os,position_index); // Write the position index
01260 
01261                         U32 normal_index = position_index;
01262                         write(os,normal_index); // Write the normal index, it is exactly the same as the normal index!
01263 
01264                         U32 base_diffuse_color_index = 0; // only using one diffuse color
01265                         write(os,base_diffuse_color_index);
01266 
01267                         U32 base_specular_color_index = 0; // only using one specular color
01268                         write(os,base_specular_color_index);
01269 
01270                         // Would need to write texture data here if we were doing that
01271 
01272                 }
01273 
01274         }
01275 
01276         }
01277         return os;
01278 }
01279 
01280 template<>
01281 ostream& U3DWriter::write(ostream& os, const string& s )
01282 {
01283         // WARNING - I AM NOT SURE IF THIS APPROACH WILL PRESENT UTF8 PROBLEMS.
01284         // See character_encoding in the file header writing module above
01285         test_type_sizes();
01286 
01287         short unsigned int size = s.size(); // To write a string you must first place its
01288         write(os,size);
01289 
01290         // Write the characters
01291         for(unsigned int i = 0; i<size; ++i) {
01292                 write(os,s[i]);
01293         }
01294 
01295         return os;
01296 }
01297 /*
01298 #include <boost/python.hpp>
01299 using namespace boost::python;
01300 
01301 BOOST_PYTHON_MODULE(gorgon)
01302 {
01303     class_<MarchingCubes>("MarchingCubes", init<>())
01304         .def("drawMesh", &MarchingCubes::drawMesh)
01305                 .def("setSurfaceValue", &MarchingCubes::setSurfaceValue)
01306                 .def("getSurfaceValue", &MarchingCubes::getSurfaceValue)
01307                 .def("set_sample_density", &MarchingCubes::set_sample_density)
01308                 .def("getSampleDensity", &MarchingCubes::getSampleDensity)
01309                 .def("loadMRC", &MarchingCubes::loadMRC)
01310     ;
01311 }
01312 */
01313 #endif

Generated on Thu May 3 10:06:25 2012 for EMAN2 by  doxygen 1.4.7