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 
00542 if ((int(glGetString(GL_VERSION)[0])-48)>2){
00543         rgbgenerator = ColorRGBGenerator();
00544 
00545 #ifdef _WIN32
00546         typedef void (APIENTRYP PFNGLGENBUFFERSPROC) (GLsizei n, GLuint *buffers);
00547         PFNGLGENBUFFERSPROC glGenBuffers;
00548         glGenBuffers = (PFNGLGENBUFFERSPROC) wglGetProcAddress("glGenBuffers");
00549 #endif  //_WIN32
00550 
00551         glGenBuffers(4, buffer);
00552 }
00553 
00554 }
00555 
00556 MarchingCubes::MarchingCubes(EMData * em)
00557         : _isodl(0)
00558 {
00559 if ((int(glGetString(GL_VERSION)[0])-48)>2){
00560         rgbgenerator = ColorRGBGenerator();
00561 
00562 #ifdef _WIN32
00563         typedef void (APIENTRYP PFNGLGENBUFFERSPROC) (GLsizei n, GLuint *buffers);
00564         PFNGLGENBUFFERSPROC glGenBuffers;
00565         glGenBuffers = (PFNGLGENBUFFERSPROC) wglGetProcAddress("glGenBuffers");
00566 #endif  //_WIN32
00567 
00568         glGenBuffers(4, buffer);                
00569         set_data(em);}
00570 else{
00571         set_data(em);
00572         }
00573 }
00574 
00575 
00576 
00577 void MarchingCubes::clear_min_max_vals()
00578 {
00579         for (vector<EMData*>::iterator it = minvals.begin(); it != minvals.end(); ++it)
00580         {
00581                 if ( (*it) != 0 ) delete *it;
00582         }
00583         minvals.clear();
00584 
00585         for (vector<EMData*>::iterator it = maxvals.begin(); it != maxvals.end(); ++it)
00586         {
00587                 if ( (*it) != 0 ) delete *it;
00588         }
00589         maxvals.clear();
00590 }
00591 
00592 bool MarchingCubes::calculate_min_max_vals()
00593 {
00594 
00595         if (_emdata == NULL ) throw NullPointerException("Error, cannot generate search tree if the overriding EMData object is NULL");
00596 
00597         clear_min_max_vals();
00598 
00599         int nx = _emdata->get_xsize();
00600         int ny = _emdata->get_ysize();
00601         int nz = _emdata->get_zsize();
00602 
00603         // Create the binary min max tree
00604         while ( nx > 1 || ny > 1 || nz > 1 )
00605         {
00606                 int size = minvals.size();
00607 
00608                 if ( size == 0 ){
00609                         Dict a;
00610                         // Setting search to 3 at the bottom level is most important.
00611                         // FIXME - d.woolford add comments heere
00612                         a["search"] = 3;
00613                         minvals.push_back(_emdata->process("math.minshrink",a));
00614                         maxvals.push_back(_emdata->process("math.maxshrink",a));
00615                 }else {
00616                         minvals.push_back(minvals[size-1]->process("math.minshrink"));
00617                         maxvals.push_back(maxvals[size-1]->process("math.maxshrink"));
00618                 }
00619 
00620                 nx = minvals[size]->get_xsize();
00621                 ny = minvals[size]->get_ysize();
00622                 nz = minvals[size]->get_zsize();
00623 #if MARCHING_CUBES_DEBUG
00624                 cout << "dims are " << nx << " " << ny << " " << nz << endl;
00625 #endif
00626         }
00627 
00628         drawing_level = -1;
00629 
00630         return true;
00631 }
00632 
00633 MarchingCubes::~MarchingCubes() {
00634         clear_min_max_vals();
00635 
00636 if ((int(glGetString(GL_VERSION)[0])-48)>2){
00637 #ifdef _WIN32
00638         typedef void (APIENTRYP PFNGLDELETEBUFFERSPROC) (GLsizei n, const GLuint *buffers);
00639         PFNGLDELETEBUFFERSPROC glDeleteBuffers;
00640         glDeleteBuffers = (PFNGLDELETEBUFFERSPROC) wglGetProcAddress("glDeleteBuffers");
00641 #endif  //_WIN32
00642 
00643         glDeleteBuffers(4, buffer);
00644 }
00645 }
00646 
00647 Dict MarchingCubes::get_isosurface()
00648 {
00649         calculate_surface();
00650         Dict d;
00651         d.put("points", (float*)pp.get_data());
00652         for (unsigned int i = 0; i < ff.elem(); ++i ) ff[i] /= 3;
00653         d.put("faces", (unsigned int*)ff.get_data());
00654         d.put("normals", (float*)nn.get_data());
00655         d.put("size", ff.elem());
00656         return d;
00657 }
00658 
00659 void MarchingCubes::surface_face_z()
00660 {
00661         float* f = pp.get_data();
00662         float* n = nn.get_data();
00663         for (unsigned int i = 0; i < pp.elem(); i += 3 ) {
00664                 if (f[i+2] == 0.5) continue;
00665                 Vec3f z(0,0,1.0);
00666                 Vec3f axis(-n[i+1],n[i],0);
00667                 axis.normalize();
00668         
00669                 Dict d;
00670                 d["type"] = "spin";
00671                 d["Omega"] =  90.f;
00672                 d["n1"] = axis[0];
00673                 d["n2"] = axis[1];
00674                 d["n3"] = 0;
00675                 Transform t(d);
00676                 Vec3f delta = t*z;
00677                 
00678                 f[i] += delta[0]*.25f;
00679                 f[i+1] += delta[1]*.25f;
00680                 f[i+2] = 0.5;
00681         }
00682         
00683         for (unsigned int i = 0; i < nn.elem(); i += 3 ) {
00684                 n[i] = 0;
00685                 n[i+1] = 0;
00686                 n[i+2] = 1;
00687         }
00688 }
00689 
00690 void MarchingCubes::set_data(EMData* data)
00691 {
00692         if ( data->get_zsize() == 1 ) throw ImageDimensionException("The z dimension of the image must be greater than 1");
00693         _emdata = data;
00694         calculate_min_max_vals();
00695         rgbgenerator.set_data(data);
00696 }
00697 
00698 void MarchingCubes::set_surface_value(const float value) {
00699 
00700         if(_surf_value == value) return;
00701 
00702         _surf_value = value;
00703 
00704 }
00705 
00706 void MarchingCubes::calculate_surface() {
00707 
00708         if ( _emdata == 0 ) throw NullPointerException("Error, attempt to generate isosurface, but the emdata image object has not been set");
00709         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");
00710         
00711         point_map.clear();
00712         pp.clear();
00713         nn.clear();
00714         ff.clear();
00715         vv.clear();
00716 
00717 #if MARCHING_CUBES_DEBUG
00718         int time0 = clock();
00719 #endif
00720 
00721         float min = minvals[minvals.size()-1]->get_value_at(0,0,0);
00722         float max = maxvals[minvals.size()-1]->get_value_at(0,0,0);
00723         if ( min < _surf_value &&  max > _surf_value) draw_cube(0,0,0,minvals.size()-1);
00724         
00725 #if MARCHING_CUBES_DEBUG
00726         int time1 = clock();
00727         cout << "It took " << (time1-time0) << " " << (float)(time1-time0)/CLOCKS_PER_SEC << " to traverse the search tree and generate polygons" << endl;
00728         cout << "... using surface value " << _surf_value << endl;
00729 #endif
00730 }
00731 
00732 void MarchingCubes::draw_cube(const int x, const int y, const int z, const int cur_level ) {
00733 
00734         if ( cur_level == drawing_level )
00735         {
00736                 EMData* e;
00737                 if ( drawing_level == - 1 ) e = _emdata;
00738                 else e = minvals[drawing_level];
00739                 if ( x < (e->get_xsize()-1) && y < (e->get_ysize()-1) && z < (e->get_zsize()-1))
00740                         marching_cube(x,y,z, cur_level);
00741         }
00742         else
00743         {
00744                 EMData* e;
00745                 if ( cur_level > 0 ) {
00746                         int xsize = minvals[cur_level-1]->get_xsize();
00747                         int ysize = minvals[cur_level-1]->get_ysize();
00748                         int zsize = minvals[cur_level-1]->get_zsize();
00749                         e = minvals[cur_level-1];
00750                         for(int i=0; i<8; ++i ) {
00751                                 int xx = 2*x+a2fVertexOffset[i][0];
00752                                 if ( xx >= xsize ) continue;
00753                                 int yy = 2*y+a2fVertexOffset[i][1];
00754                                 if ( yy >= ysize ) continue;
00755                                 int zz = 2*z+a2fVertexOffset[i][2];
00756                                 if ( zz >= zsize ) continue;
00757 
00758                                 float min = minvals[cur_level-1]->get_value_at(xx,yy,zz);
00759                                 float max = maxvals[cur_level-1]->get_value_at(xx,yy,zz);
00760                                 if ( min < _surf_value &&  max > _surf_value)
00761                                         draw_cube(xx,yy,zz,cur_level-1);
00762                         }
00763                 }
00764                 else {
00765                         e = _emdata;
00766                         for(int i=0; i<8; ++i ) {
00767                                         draw_cube(2*x+a2fVertexOffset[i][0],2*y+a2fVertexOffset[i][1],2*z+a2fVertexOffset[i][2],cur_level-1);
00768                         }
00769                 }
00770 
00771                 if ( x == (minvals[cur_level]->get_xsize()-1) ) {
00772                         if ( e->get_xsize() > 2*x ){
00773                                 for(int i=0; i<4; ++i ) {
00774                                         draw_cube(2*x+a2fPosXOffset[i][0],2*y+a2fPosXOffset[i][1],2*z+a2fPosXOffset[i][2],cur_level-1);
00775                                 }
00776                         }
00777                         if ( y == (minvals[cur_level]->get_ysize()-1) ) {
00778                                 if ( e->get_ysize() > 2*y ) {
00779                                         for(int i=0; i<2; ++i ) {
00780                                                 draw_cube(2*x+a2fPosXPosYOffset[i][0],2*y+a2fPosXPosYOffset[i][1],2*z+a2fPosXPosYOffset[i][2],cur_level-1);
00781                                         }
00782                                 }
00783                                 if (  z == (minvals[cur_level]->get_zsize()-1) ){
00784                                         if ( e->get_zsize() > 2*z ) {
00785                                                 draw_cube(2*x+2,2*y+2,2*z+2,cur_level-1);
00786                                         }
00787                                 }
00788                         }
00789                         if ( z == (minvals[cur_level]->get_zsize()-1) ) {
00790                                 if ( e->get_zsize() > 2*z ) {
00791                                         for(int i=0; i<2; ++i ) {
00792                                                 draw_cube(2*x+a2fPosXPosZOffset[i][0],2*y+a2fPosXPosZOffset[i][1],2*z+a2fPosXPosZOffset[i][2],cur_level-1);
00793                                         }
00794                                 }
00795                         }
00796                 }
00797                 if ( y == (minvals[cur_level]->get_ysize()-1) ) {
00798                         if ( e->get_ysize() > 2*y ) {
00799                                 for(int i=0; i<4; ++i ) {
00800                                         draw_cube(2*x+a2fPosYOffset[i][0],2*y+a2fPosYOffset[i][1],2*z+a2fPosYOffset[i][2],cur_level-1);
00801                                 }
00802                         }
00803                         if ( z == (minvals[cur_level]->get_ysize()-1) ) {
00804                                 if ( e->get_zsize() > 2*z ) {
00805                                         for(int i=0; i<2; ++i ) {
00806                                                 draw_cube(2*x+a2fPosYPosZOffset[i][0],2*y+a2fPosYPosZOffset[i][1],2*z+a2fPosYPosZOffset[i][2],cur_level-1);
00807                                         }
00808                                 }
00809                         }
00810                 }
00811                 if ( z == (minvals[cur_level]->get_zsize()-1) ) {
00812                         if ( e->get_zsize() > 2*z ) {
00813                                 for(int i=0; i<4; ++i ) {
00814                                         draw_cube(2*x+a2fPosZOffset[i][0],2*y+a2fPosZOffset[i][1],2*z+a2fPosZOffset[i][2],cur_level-1);
00815                                 }
00816                         }
00817                 }
00818 
00819         }
00820 }
00821 
00822 void MarchingCubes::get_normal(Vector3 &normal, int fX, int fY, int fZ)
00823 {
00824     normal[0] = _emdata->get_value_at(fX-1, fY, fZ) - _emdata->get_value_at(fX+1, fY, fZ);
00825     normal[1] = _emdata->get_value_at(fX, fY-1, fZ) - _emdata->get_value_at(fX, fY+1, fZ);
00826     normal[2] = _emdata->get_value_at(fX, fY, fZ-1) - _emdata->get_value_at(fX, fY, fZ+1);
00827     normal.normalize();
00828 }
00829 
00830 float MarchingCubes::get_offset(float fValue1, float fValue2, float fValueDesired)
00831 {
00832         float fDelta = fValue2 - fValue1;
00833 
00834         if(fDelta == 0.0f)
00835         {
00836                 return 0.5f;
00837         }
00838         return (fValueDesired - fValue1)/fDelta;
00839 }
00840 
00841 int MarchingCubes::get_edge_num(int x, int y, int z, int edge) {
00842         // edge direction is right, down, back (x, y, z)
00843         unsigned int index = 0;
00844         index = (x << 22) | (y << 12) | (z << 2) | edge;
00845         return index;
00846 }
00847 
00848 void MarchingCubes::color_vertices()
00849 {
00850         cc.clear();
00851 #ifdef _WIN32
00852         int scaling = (int)pow(2.0,drawing_level + 1);          // Needed to account for sampling rate
00853 #else
00854         int scaling = pow(2,drawing_level + 1);         // Needed to account for sampling rate
00855 #endif  //_WIN32        
00856         //Color vertices. We don't need to rerun marching cubes on color vertices, so this method improves effciency
00857         for(unsigned int i = 0; i < vv.elem(); i+=3){
00858                 float* color = rgbgenerator.getRGBColor(scaling*vv[i], scaling*vv[i+1], scaling*vv[i+2]);
00859                 cc.push_back_3(color);
00860         }
00861         rgbgenerator.setNeedToRecolor(false);
00862 }
00863                 
00864 void MarchingCubes::marching_cube(int fX, int fY, int fZ, int cur_level)
00865 {
00866 //      extern int aiCubeEdgeFlags[256];
00867 //      extern int a2iTriangleConnectionTable[256][16];
00868 
00869         int iCorner, iVertex, iVertexTest, iEdge, iTriangle, iFlagIndex, iEdgeFlags;
00870         float fOffset;
00871         Vector3 sColor;
00872         float afCubeValue[8];
00873         float asEdgeVertex[12][3];
00874         int pointIndex[12];
00875 
00876         int fxScale = 1, fyScale = 1, fzScale = 1;
00877         if ( cur_level != -1 )
00878         {
00879                 fxScale = _emdata->get_xsize()/minvals[cur_level]->get_xsize();
00880                 fyScale = _emdata->get_ysize()/minvals[cur_level]->get_ysize();
00881                 fzScale = _emdata->get_zsize()/minvals[cur_level]->get_zsize();
00882                 for(iVertex = 0; iVertex < 8; iVertex++)
00883                 {
00884                         afCubeValue[iVertex] = _emdata->get_value_at( fxScale*(fX + a2fVertexOffset[iVertex][0]),
00885                                                 fyScale*(fY + a2fVertexOffset[iVertex][1]), fzScale*(fZ + a2fVertexOffset[iVertex][2]));
00886                 }
00887         }
00888         else
00889         {
00890                 //Make a local copy of the values at the cube's corners
00891                 for(iVertex = 0; iVertex < 8; iVertex++)
00892                 {
00893                         afCubeValue[iVertex] = _emdata->get_value_at( fX + a2fVertexOffset[iVertex][0],
00894                                         fY + a2fVertexOffset[iVertex][1], fZ + a2fVertexOffset[iVertex][2]);
00895                 }
00896         }
00897 
00898         //Find which vertices are inside of the surface and which are outside
00899         iFlagIndex = 0;
00900         for(iVertexTest = 0; iVertexTest < 8; iVertexTest++)
00901         {
00902                 if (_surf_value >= 0 ){
00903                         if(afCubeValue[iVertexTest] <= _surf_value)
00904                                 iFlagIndex |= 1<<iVertexTest;
00905                 }
00906                 else {
00907                         if(afCubeValue[iVertexTest] >= _surf_value)
00908                                 iFlagIndex |= 1<<iVertexTest;
00909                 }
00910         }
00911 
00912         //Find which edges are intersected by the surface
00913         iEdgeFlags = aiCubeEdgeFlags[iFlagIndex];
00914 
00915         //If the cube is entirely inside or outside of the surface, then there will be no intersections
00916         if(iEdgeFlags == 0) return;
00917 
00918         //Find the point of intersection of the surface with each edge
00919         //Then find the normal to the surface at those points
00920         for(iEdge = 0; iEdge < 12; iEdge++)
00921         {
00922                 //if there is an intersection on this edge
00923                 if(iEdgeFlags & (1<<iEdge))
00924                 {
00925                         fOffset = get_offset(afCubeValue[ a2iEdgeConnection[iEdge][0] ],
00926                                                                  afCubeValue[ a2iEdgeConnection[iEdge][1] ], _surf_value);
00927 
00928                         if ( cur_level == -1 ){
00929                                 asEdgeVertex[iEdge][0] = fX + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][0]  +  fOffset * a2fEdgeDirection[iEdge][0]) + 0.5f;
00930                                 asEdgeVertex[iEdge][1] = fY + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][1]  +  fOffset * a2fEdgeDirection[iEdge][1]) + 0.5f;
00931                                 asEdgeVertex[iEdge][2] = fZ + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][2]  +  fOffset * a2fEdgeDirection[iEdge][2]) + 0.5f;
00932                         } else {
00933                                 asEdgeVertex[iEdge][0] = fxScale*(fX + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][0]  +  fOffset * a2fEdgeDirection[iEdge][0])) + 0.5f;
00934                                 asEdgeVertex[iEdge][1] = fyScale*(fY + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][1]  +  fOffset * a2fEdgeDirection[iEdge][1])) + 0.5f;
00935                                 asEdgeVertex[iEdge][2] = fzScale*(fZ + (a2fVertexOffset[ a2iEdgeConnection[iEdge][0] ][2]  +  fOffset * a2fEdgeDirection[iEdge][2])) + 0.5f;
00936                         }
00937 
00938                         pointIndex[iEdge] = get_edge_num(fX+edgeLookUp[iEdge][0], fY+edgeLookUp[iEdge][1], fZ+edgeLookUp[iEdge][2], edgeLookUp[iEdge][3]);
00939                 }
00940         }
00941         
00942         //Save voxel coords for later color processing
00943         int vox[3] = {fX, fY, fZ};
00944         
00945         //Draw the triangles that were found.  There can be up to five per cube
00946         for(iTriangle = 0; iTriangle < 5; iTriangle++)
00947         {
00948                 if(a2iTriangleConnectionTable[iFlagIndex][3*iTriangle] < 0)
00949                         break;
00950 
00951                 float pts[3][3];
00952                 for(iCorner = 0; iCorner < 3; iCorner++)
00953                 {
00954                         iVertex = a2iTriangleConnectionTable[iFlagIndex][3*iTriangle+iCorner];
00955                         memcpy(&pts[iCorner][0],  &asEdgeVertex[iVertex][0], 3*sizeof(float));
00956                 }
00957 
00958                 
00959                 
00960                 float v1[3] = {pts[1][0]-pts[0][0],pts[1][1]-pts[0][1],pts[1][2]-pts[0][2]};
00961                 float v2[3] = {pts[2][0]-pts[1][0],pts[2][1]-pts[1][1],pts[2][2]-pts[1][2]};
00962                 
00963                 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] };
00964                         
00965                 
00966                 for(iCorner = 0; iCorner < 3; iCorner++)
00967                 {
00968                         // Without vertex normalization
00969 //                      iVertex = a2iTriangleConnectionTable[iFlagIndex][3*iTriangle+iCorner];
00970 //                      int ss = pp.elem();
00971 //                      pp.push_back_3(&pts[iCorner][0]);
00972 //                      nn.push_back_3(&n[0]);
00973 //                      ff.push_back(ss);
00974 
00975 //                      With vertex normalization
00976                         iVertex = a2iTriangleConnectionTable[iFlagIndex][3*iTriangle+iCorner];
00977                         map<int,int>::iterator it = point_map.find(pointIndex[iVertex]);
00978                         if ( it == point_map.end() ){
00979                                 vv.push_back_3(&vox[0]);
00980                                 int ss = pp.elem();
00981                                 pp.push_back_3(&pts[iCorner][0]);
00982                                 nn.push_back_3(&n[0]);
00983                                 ff.push_back(ss);
00984                                 point_map[pointIndex[iVertex]] = ss;
00985                         } else {
00986                                 int idx = it->second;
00987                                 ff.push_back(idx);
00988                                 nn[idx] += n[0];
00989                                 nn[idx+1] += n[1];
00990                                 nn[idx+2] += n[2];
00991                         }
00992                 }
00993         }
00994 }
00995 
00996 
00997 
00998 U3DWriter::U3DWriter() : DIFFUSE_COLOR_COUNT(1),
00999                 SPECULAR_COLOR_COUNT(1)
01000 {
01001 }
01002 U3DWriter::~U3DWriter() {}
01003 
01004 using std::ofstream;
01005 
01006 int U3DWriter::write(const string& filename) {
01007         // Must open the ofstream in binary mode
01008         ofstream of(filename.c_str(), ofstream::binary);
01009         write(of);
01010         of.close();
01011 
01012         return 1;
01013 }
01014 
01015 ostream& U3DWriter::write(ostream& os) {
01016         write_header(os);
01017         return os;
01018 }
01019 
01020 unsigned int U3DWriter::size_of_in_bytes(){
01021         // this is just the size in bytes of all of the entries in this object that are written to the binary
01022         // output. This is based only on the behavior of write_header and needs to be udpated
01023         unsigned size = 4+2+2+4+4+8+4+8; // 36 bytes
01024         return size;
01025 }
01026 ostream& U3DWriter::write_header(ostream& os)
01027 {
01028         // checks
01029         test_type_sizes();
01030 
01031         U32 block_type_file_header = 0x00443355; // This is the block type tag of a file header, as taken form the ECMA spec
01032         write( os, block_type_file_header);
01033 //
01034         I16 major_version = -1; // Compliance has not been tested for this encoder, so we must specify a value less than 0
01035         I16 minor_version =  0; // This is the version of this encoder, which we are calling '0'
01036         write( os, major_version);
01037         write( os, minor_version);
01038 
01039         U32 profile_identifier = 0x00000000; // Base profile - no optional features ares used
01040         write( os, profile_identifier);
01041 
01042         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
01043         write( os, declaration_size);
01044 
01045         U64 file_size = size_of_in_bytes(); // This is the size of the file in bytes
01046         write( os, file_size);
01047 
01048         U32 character_encoding = 106; // UTF-8 MIB from the IANA
01049         write( os, character_encoding);
01050 
01051         F64 unit_scaling = 1.0; // This should eventually be the correct scaling of the objects
01052         write( os, unit_scaling);
01053 
01054 
01055         return os;
01056 }
01057 
01058 void U3DWriter::test_type_sizes()
01059 {
01060         bool error = false;
01061         if (sizeof(F64) != 8 ){
01062                 cout << "Error, size of double is not 64 bytes, it's " << sizeof(F64)*4 << endl;
01063                 error = true;
01064         }
01065         if (sizeof(F32) != 4 ){
01066                 cout << "Error, size of float is not 32 bytes, it's " << sizeof(F32)*4 << endl;
01067                 error = true;
01068         }
01069         if (sizeof(U64) != 8) {
01070                 cout << "Error, size of long unsigned int is not 64 bytes, it's " << sizeof(U64)*4 << endl;
01071                 error = true;
01072         }
01073         if (sizeof(U32) != 4) {
01074                 cout << "Error, size of unsigned int is not 32 bytes, it's " << sizeof(U32)*4 << endl;
01075                 error = true;
01076         }
01077         if (sizeof(I16) != 2) {
01078                 cout << "Error, size of short int is not 16 bytes, it's " << sizeof(I16)*4 << endl;
01079                 error = true;
01080         }
01081         if (sizeof(U16) != 2) {
01082                 cout << "Error, size of short unsigned int is not 16 bytes, it's " << sizeof(U16)*4 << endl;
01083                 error = true;
01084         }
01085         if (sizeof(U8) != 1) {
01086                 cout << "Error, size of unsigned char is not  bytes, it's " << sizeof(U8)*4 << endl;
01087                 error = true;
01088         }
01089 
01090         if (error) {
01091                 throw;
01092         }
01093 }
01094 
01095 ostream& U3DWriter::write_clod_mesh_generator_node(ostream& os)
01096 {
01097         /*
01098         CLOD Mesh Declaration
01099         */
01100         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
01101         write( os, block_type_clod_mesh_generator);
01102 
01103         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,"")
01104         write(os,name);
01105 
01106         U32 chain_index = 0x00000000; // the value of Chain Index shall be zero for this type - as in the SPEC
01107         write( os, chain_index);
01108 
01109         /*
01110         Max Mesh Description
01111         */
01112         U32 mesh_attributes = 0x00000000; // Faces in the mesh have a normal index at each corner 0x00000001 is used to exlude normals POTENTIALLY USEFUL
01113         write(os,mesh_attributes);
01114         U32 face_count = ff.elem(); // The number of faces TO FILL IN LATER
01115         write(os,face_count);
01116         U32 position_count = pp.elem(); // The number of positions in the position array TO FILL IN LATER
01117         write(os,position_count);
01118         U32 normal_count = nn.elem(); // The number of normals in the normal array TO FILL IN LATER
01119         write(os,normal_count);
01120         U32 diffuse_color_count = DIFFUSE_COLOR_COUNT; // The number of colors in the diffuse color array TO FILL IN LATER
01121         write(os,diffuse_color_count);
01122         U32 specular_color_count = SPECULAR_COLOR_COUNT; // The number of colors in the specular color array TO FILL IN LATER
01123         write(os,specular_color_count);
01124         U32 texture_coord_count = 0x00000000; // The number of texture coordinates in teh texture coordinate array POTENTIALLY USEFUL
01125         write(os,texture_coord_count);
01126         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
01127         write(os,shading_count);
01128 
01129         /*
01130         Shading  Description
01131         */
01132         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
01133         write(os,shading_attributes);
01134         U32 texture_layout_count = 0x00000000; // The number of texture layers used by this shader list
01135         write(os,texture_layout_count);
01136         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
01137         write(os,texture_coord_dimensions);
01138         U32 original_shading_id = 0; // Original shading index for this shader list. Informative parameter only. This is shader 0
01139         write(os,original_shading_id);
01140 
01141         /*
01142         CLOD Description - describes the range of resolutions available for the continuous level of detail mesh
01143         If there were more than one level of detail than these two numbers would be different
01144         */
01145         U32 minimum_resolution = position_count; // the number of positions in the base mesh
01146         write(os,minimum_resolution);
01147         U32 final_maximum_resolution = position_count; // the number of positions in the Max Mesh Description (by definition)
01148         write(os,final_maximum_resolution);
01149 
01150         /*
01151         Resource Description
01152         */
01153 
01154         /*
01155         Quality Factors
01156         for information only. Not used by the renderer. Helpful for conveying information to the user
01157         */
01158         U32 position_quality_factor = 0x00000000; // position quality factor. Descriptive information only
01159         write(os,position_quality_factor);
01160         U32 normal_quality_factor = 0x00000000; // normal quality factor. Descriptive information only
01161         write(os,normal_quality_factor);
01162         U32 texture_coord_quality_factor = 0x00000000; // texture coordinate quality factor. Descriptive information only
01163         write(os,texture_coord_quality_factor);
01164 
01165         /*
01166         Inverse Quantization
01167         used to reconstruct floating point values that have been quantized.
01168         */
01169         F32 postion_inverse_quant = 1.0; // inverse quantization factor used in the reconstruction of the position vectors
01170         write(os,postion_inverse_quant);
01171         F32 normal_inverse_quant = 1.0; // inverse quantization factor used in the reconstruction of the normal vectors
01172         write(os,normal_inverse_quant);
01173         F32 texture_coord_inverse_quant = 1.0; // inverse quantization factor used in the reconstruction of the texture coordinates
01174         write(os,texture_coord_inverse_quant);
01175         F32 diffuse_color_inverse_quant = 1.0; // inverse quantization factor used in the reconstruction of the diffuse colors
01176         write(os,diffuse_color_inverse_quant);
01177         F32 specular_color_inverse_quant = 1.0; // inverse quantization factor used in the reconstruction of the specular colors
01178         write(os,specular_color_inverse_quant);
01179 
01180         /*
01181         Resource Parameters
01182         parameters that help to define the conversion from the Author Mesh format to the Render Mesh format
01183         */
01184 
01185         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.
01186         write(os,normal_crease_parameter);
01187         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'
01188         write(os,normal_update_parameter);
01189         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
01190         write(os,normal_tolerance_parameter);
01191 
01192         /*
01193         Skeleton description
01194         used in bones-based animation by the Animation Modifier
01195         */
01196         U32 bone_count = 0x00000000; // The number of bones associated with this mesh. We will always have 0, but this could change (unlikely).
01197         write(os,bone_count);
01198 
01199         // WARNING - if bone_count is ever greater than one, then more writing would have to occur here
01200 
01202         /*
01203         Base mesh is basically the minimum LOD mesh - it must exist
01204         */
01205         {
01206         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
01207         write( os, block_type_clod_base_mesh_continuation);
01208 
01209         write(os,name); // We use the same name as above
01210 
01211         U32 chain_index = 0x00000000; // the value of Chain Index shall be zero for this type - as in the SPEC
01212         write( os, chain_index);
01213 
01214         /*
01215         Base Mesh Description
01216         */
01217         U32     base_face_count = ff.elem(); // The number of faces in the base mesh TO FILL IN LATER
01218         write( os, base_face_count);
01219         U32     base_position_count = pp.elem(); // The number of positions used by base mesh in the position array TO FILL IN LATER
01220         write( os, base_position_count);
01221         U32 base_normal_count = nn.elem(); // The number of normals used by the base mesh in the normal array TO FILL IN LATER
01222         write( os, base_normal_count);
01223         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
01224         write( os, base_diffuse_color_count);
01225         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
01226         write( os, base_specular_color_count);
01227         U32 base_texture_coord_count = 0x00000000; // The number of texture coordinates used by the base mesh in texture coordinate array TO FILL IN LATER
01228         write( os, base_texture_coord_count);
01229 
01230         /*
01231         Base mesh data
01232         */
01233 
01234         // Write position data
01235         F32* data = pp.get_data();
01236         for(unsigned int i = 0; i < pp.elem(); ++i) {
01237                 write(os,data[i]);
01238         }
01239 
01240         // Write normal data
01241         data = nn.get_data();
01242         for(unsigned int i = 0; i < nn.elem(); ++i) {
01243                 write(os,data[i]);
01244         }
01245 
01246         // Write Diffuse color, this is in rgba format
01247         F32 diffuse_rgba[4] = {1.0,0.0,0.0,1.0};
01248         for (unsigned int i = 0; i < 4; ++i) {
01249                 write(os,diffuse_rgba[i]);
01250         }
01251 
01252         // Write Specular color, this is in rgba format
01253         F32 specular_rgba[4] = {1.0,0.0,0.0,1.0};
01254         for (unsigned int i = 0; i < 4; ++i) {
01255                 write(os,specular_rgba[i]);
01256         }
01257 
01258         // THERE ARE NO TEXTURE COORDINATES, BUT IF THERE WERE WE WOULD HAVE TO WRITE THEM HERE
01259         // i.e. for i in range(0,base_texture_coord_count) write texture coords
01260 
01261         // Write normal data
01262         U32* faces = ff.get_data();
01263         for(unsigned int i = 0; i < pp.elem(); i = i+3) {
01264                 U32 shading_id = 0; // We only have one shader defined. This could could changed in future
01265                 write(os,shading_id);
01266 
01267                 // write base corner info - there are always three corners
01268                 for (unsigned int j =0; j < 3; ++j){
01269                         U32 position_index = faces[i+j];
01270                         write(os,position_index); // Write the position index
01271 
01272                         U32 normal_index = position_index;
01273                         write(os,normal_index); // Write the normal index, it is exactly the same as the normal index!
01274 
01275                         U32 base_diffuse_color_index = 0; // only using one diffuse color
01276                         write(os,base_diffuse_color_index);
01277 
01278                         U32 base_specular_color_index = 0; // only using one specular color
01279                         write(os,base_specular_color_index);
01280 
01281                         // Would need to write texture data here if we were doing that
01282 
01283                 }
01284 
01285         }
01286 
01287         }
01288         return os;
01289 }
01290 
01291 template<>
01292 ostream& U3DWriter::write(ostream& os, const string& s )
01293 {
01294         // WARNING - I AM NOT SURE IF THIS APPROACH WILL PRESENT UTF8 PROBLEMS.
01295         // See character_encoding in the file header writing module above
01296         test_type_sizes();
01297 
01298         short unsigned int size = s.size(); // To write a string you must first place its
01299         write(os,size);
01300 
01301         // Write the characters
01302         for(unsigned int i = 0; i<size; ++i) {
01303                 write(os,s[i]);
01304         }
01305 
01306         return os;
01307 }
01308 /*
01309 #include <boost/python.hpp>
01310 using namespace boost::python;
01311 
01312 BOOST_PYTHON_MODULE(gorgon)
01313 {
01314     class_<MarchingCubes>("MarchingCubes", init<>())
01315         .def("drawMesh", &MarchingCubes::drawMesh)
01316                 .def("setSurfaceValue", &MarchingCubes::setSurfaceValue)
01317                 .def("getSurfaceValue", &MarchingCubes::getSurfaceValue)
01318                 .def("set_sample_density", &MarchingCubes::set_sample_density)
01319                 .def("getSampleDensity", &MarchingCubes::getSampleDensity)
01320                 .def("loadMRC", &MarchingCubes::loadMRC)
01321     ;
01322 }
01323 */
01324 #endif

Generated on Tue Jun 11 12:40:23 2013 for EMAN2 by  doxygen 1.4.7