Developer Documentation
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
MarchingCubes.cpp
1 /*
2 Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
3 All rights reserved.
4 
5 Redistribution and use in source and binary forms, with or without modification,
6 are permitted provided that the following conditions are met:
7 
8 Redistributions of source code must retain the above copyright notice, this list of
9 conditions and the following disclaimer. Redistributions in binary form must reproduce
10 the above copyright notice, this list of conditions and the following disclaimer
11 in the documentation and/or other materials provided with the distribution.
12 
13 Neither the name of the Johns Hopkins University nor the names of its contributors
14 may be used to endorse or promote products derived from this software without specific
15 prior written permission.
16 
17 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
18 EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
19 OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
20 SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
21 INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
22 TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
23 BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
24 CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
25 ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
26 DAMAGE.
27 */
28 #include <math.h>
29 #include "MarchingCubes.h"
30 
32 // Square //
34 int Square::AntipodalCornerIndex(int idx){
35  int x,y;
36  FactorCornerIndex(idx,x,y);
37  return CornerIndex( (x+1)%2 , (y+1)%2 );
38 }
39 int Square::CornerIndex(int x,int y){return (y<<1)|x;}
40 void Square::FactorCornerIndex(int idx,int& x,int& y){
41  x=(idx>>0)%2;
42  y=(idx>>1)%2;
43 }
44 int Square::EdgeIndex(int orientation,int i){
45  switch(orientation){
46  case 0: // x
47  if(!i) {return 0;} // (0,0) -> (1,0)
48  else {return 2;} // (0,1) -> (1,1)
49  case 1: // y
50  if(!i) {return 3;} // (0,0) -> (0,1)
51  else {return 1;} // (1,0) -> (1,1)
52  };
53  return -1;
54 }
55 void Square::FactorEdgeIndex(int idx,int& orientation,int& i){
56  switch(idx){
57  case 0: case 2:
58  orientation=0;
59  i=idx/2;
60  return;
61  case 1: case 3:
62  orientation=1;
63  i=((idx/2)+1)%2;
64  return;
65  };
66 }
67 void Square::EdgeCorners(int idx,int& c1,int& c2){
68  int orientation,i;
69  FactorEdgeIndex(idx,orientation,i);
70  switch(orientation){
71  case 0:
72  c1=CornerIndex(0,i);
73  c2=CornerIndex(1,i);
74  break;
75  case 1:
76  c1=CornerIndex(i,0);
77  c2=CornerIndex(i,1);
78  break;
79  };
80 }
81 int Square::ReflectEdgeIndex(int idx,int edgeIndex){
82  int orientation=edgeIndex%2;
83  int o,i;
84  FactorEdgeIndex(idx,o,i);
85  if(o!=orientation){return idx;}
86  else{return EdgeIndex(o,(i+1)%2);}
87 }
88 int Square::ReflectCornerIndex(int idx,int edgeIndex){
89  int orientation=edgeIndex%2;
90  int x,y;
91  FactorCornerIndex(idx,x,y);
92  switch(orientation){
93  case 0: return CornerIndex((x+1)%2,y);
94  case 1: return CornerIndex(x,(y+1)%2);
95  };
96  return -1;
97 }
98 
99 
100 
102 // Cube //
104 int Cube::CornerIndex(int x,int y,int z){return (z<<2)|(y<<1)|x;}
105 void Cube::FactorCornerIndex(int idx,int& x,int& y,int& z){
106  x=(idx>>0)%2;
107  y=(idx>>1)%2;
108  z=(idx>>2)%2;
109 }
110 int Cube::EdgeIndex(int orientation,int i,int j){return (i | (j<<1))|(orientation<<2);}
111 void Cube::FactorEdgeIndex( int idx , int& orientation , int& i , int &j )
112 {
113  orientation=idx>>2;
114  i = (idx&1);
115  j = (idx&2)>>1;
116 }
117 int Cube::FaceIndex(int x,int y,int z){
118  if (x<0) {return 0;}
119  else if (x>0) {return 1;}
120  else if (y<0) {return 2;}
121  else if (y>0) {return 3;}
122  else if (z<0) {return 4;}
123  else if (z>0) {return 5;}
124  else {return -1;}
125 }
126 int Cube::FaceIndex(int dir,int offSet){return (dir<<1)|offSet;}
127 
128 void Cube::FactorFaceIndex(int idx,int& x,int& y,int& z){
129  x=y=z=0;
130  switch(idx){
131  case 0: x=-1; break;
132  case 1: x= 1; break;
133  case 2: y=-1; break;
134  case 3: y= 1; break;
135  case 4: z=-1; break;
136  case 5: z= 1; break;
137  };
138 }
139 void Cube::FactorFaceIndex(int idx,int& dir,int& offSet){
140  dir = idx>>1;
141  offSet=idx &1;
142 }
143 
144 int Cube::FaceAdjacentToEdges(int eIndex1,int eIndex2){
145  int f1,f2,g1,g2;
146  FacesAdjacentToEdge(eIndex1,f1,f2);
147  FacesAdjacentToEdge(eIndex2,g1,g2);
148  if(f1==g1 || f1==g2){return f1;}
149  if(f2==g1 || f2==g2){return f2;}
150  return -1;
151 }
152 
153 void Cube::FacesAdjacentToEdge(int eIndex,int& f1Index,int& f2Index){
154  int orientation,i1,i2;
155  FactorEdgeIndex(eIndex,orientation,i1,i2);
156  i1<<=1;
157  i2<<=1;
158  i1--;
159  i2--;
160  switch(orientation){
161  case 0:
162  f1Index=FaceIndex( 0,i1, 0);
163  f2Index=FaceIndex( 0, 0,i2);
164  break;
165  case 1:
166  f1Index=FaceIndex(i1, 0, 0);
167  f2Index=FaceIndex( 0, 0,i2);
168  break;
169  case 2:
170  f1Index=FaceIndex(i1, 0, 0);
171  f2Index=FaceIndex( 0,i2, 0);
172  break;
173  };
174 }
175 void Cube::EdgeCorners(int idx,int& c1,int& c2){
176  int orientation,i1,i2;
177  FactorEdgeIndex(idx,orientation,i1,i2);
178  switch(orientation){
179  case 0:
180  c1=CornerIndex(0,i1,i2);
181  c2=CornerIndex(1,i1,i2);
182  break;
183  case 1:
184  c1=CornerIndex(i1,0,i2);
185  c2=CornerIndex(i1,1,i2);
186  break;
187  case 2:
188  c1=CornerIndex(i1,i2,0);
189  c2=CornerIndex(i1,i2,1);
190  break;
191  };
192 }
193 void Cube::FaceCorners(int idx,int& c1,int& c2,int& c3,int& c4){
194  int i=idx%2;
195  switch(idx/2){
196  case 0:
197  c1=CornerIndex(i,0,0);
198  c2=CornerIndex(i,1,0);
199  c3=CornerIndex(i,0,1);
200  c4=CornerIndex(i,1,1);
201  return;
202  case 1:
203  c1=CornerIndex(0,i,0);
204  c2=CornerIndex(1,i,0);
205  c3=CornerIndex(0,i,1);
206  c4=CornerIndex(1,i,1);
207  return;
208  case 2:
209  c1=CornerIndex(0,0,i);
210  c2=CornerIndex(1,0,i);
211  c3=CornerIndex(0,1,i);
212  c4=CornerIndex(1,1,i);
213  return;
214  }
215 }
216 int Cube::AntipodalCornerIndex(int idx){
217  int x,y,z;
218  FactorCornerIndex(idx,x,y,z);
219  return CornerIndex((x+1)%2,(y+1)%2,(z+1)%2);
220 }
221 int Cube::FaceReflectFaceIndex(int idx,int faceIndex){
222  if(idx/2!=faceIndex/2){return idx;}
223  else{
224  if(idx%2) {return idx-1;}
225  else {return idx+1;}
226  }
227 }
228 int Cube::FaceReflectEdgeIndex(int idx,int faceIndex){
229  int orientation=faceIndex/2;
230  int o,i,j;
231  FactorEdgeIndex(idx,o,i,j);
232  if(o==orientation){return idx;}
233  switch(orientation){
234  case 0: return EdgeIndex(o,(i+1)%2,j);
235  case 1:
236  switch(o){
237  case 0: return EdgeIndex(o,(i+1)%2,j);
238  case 2: return EdgeIndex(o,i,(j+1)%2);
239  };
240  case 2: return EdgeIndex(o,i,(j+1)%2);
241  };
242  return -1;
243 }
244 int Cube::FaceReflectCornerIndex(int idx,int faceIndex){
245  int orientation=faceIndex/2;
246  int x,y,z;
247  FactorCornerIndex(idx,x,y,z);
248  switch(orientation){
249  case 0: return CornerIndex((x+1)%2,y,z);
250  case 1: return CornerIndex(x,(y+1)%2,z);
251  case 2: return CornerIndex(x,y,(z+1)%2);
252  };
253  return -1;
254 }
255 int Cube::EdgeReflectCornerIndex(int idx,int edgeIndex){
256  int orientation,x,y,z;
257  FactorEdgeIndex(edgeIndex,orientation,x,y);
258  FactorCornerIndex(idx,x,y,z);
259  switch(orientation){
260  case 0: return CornerIndex( x ,(y+1)%2,(z+1)%2);
261  case 1: return CornerIndex((x+1)%2, y ,(z+1)%2);
262  case 2: return CornerIndex((x+1)%2,(y+1)%2, z );
263  };
264  return -1;
265 }
266 int Cube::EdgeReflectEdgeIndex(int edgeIndex){
267  int o,i1,i2;
268  FactorEdgeIndex(edgeIndex,o,i1,i2);
269  return Cube::EdgeIndex(o,(i1+1)%2,(i2+1)%2);
270 }
271 
272 
274 // MarchingSquares //
276 /*
277 0} // (0,0) -> (1,0)
278 1} // (1,0) -> (1,1)
279 2} // (0,1) -> (1,1)
280 3} // (0,0) -> (0,1)
281 */
282 const int MarchingSquares::edgeMask[1<<Square::CORNERS]={
283  0, // 0 -> -> ->
284  9, // 1 -> 0 -> (0,0) -> 0,3 -> 9
285  3, // 2 -> 1 -> (1,0) -> 0,1 -> 3
286  10, // 3 -> 0,1 -> (0,0) (1,0) -> 1,3 -> 10
287  12, // 4 -> 2 -> (0,1) -> 2,3 -> 12
288  5, // 5 -> 0,2 -> (0,0) (0,1) -> 0,2 -> 5
289  15, // 6 -> 1,2 -> (1,0) (0,1) -> 0,1,2,3 -> 15
290  6, // 7 -> 0,1,2 -> (0,0) (1,0) (0,1) -> 1,2 -> 6
291  6, // 8 -> 3 -> (1,1) -> 1,2 -> 6
292  15, // 9 -> 0,3 -> (0,0) (1,1) -> 0,1,2,3 -> 15
293  5, // 10 -> 1,3 -> (1,0) (1,1) -> 0,2 -> 5
294  12, // 11 -> 0,1,3 -> (0,0) (1,0) (1,1) -> 2,3 -> 12
295  10, // 12 -> 2,3 -> (0,1) (1,1) -> 1,3 -> 10
296  3, // 13 -> 0,2,3 -> (0,0) (0,1) (1,1) -> 0,1 -> 3
297  9, // 14 -> 1,2,3 -> (1,0) (0,1) (1,1) -> 0,3 -> 9
298  0, // 15 -> 0,1,2,3 -> (0,0) (1,0) (0,1) (1,1) ->
299 };
300 const int MarchingSquares::edges[1<<Square::CORNERS][MAX_EDGES*2+1] = {
301  { -1, -1, -1, -1, -1}, //
302  { 3, 0, -1, -1, -1}, // (0,0)
303  { 0, 1, -1, -1, -1}, // (1,0)
304  { 3, 1, -1, -1, -1}, // (0,0) (1,0)
305  { 2, 3, -1, -1, -1}, // (0,1)
306  { 2, 0, -1, -1, -1}, // (0,0) (0,1)
307  { 0, 1, 2, 3, -1}, // (1,0) (0,1)
308  { 1, 2, -1, -1, -1}, // (0,0) (1,0) (0,1)
309  { 2, 1, -1, -1, -1}, // (1,1)
310  { 3, 0, 1, 2, -1}, // (0,0) (1,1)
311  { 0, 2, -1, -1, -1}, // (1,0) (1,1)
312  { 3, 2, -1, -1, -1}, // (0,0) (1,0) (1,1)
313  { 1, 3, -1, -1, -1}, // (0,1) (1,1)
314  { 1, 0, -1, -1, -1}, // (0,0) (0,1) (1,1)
315  { 0, 3, -1, -1, -1}, // (1,0) (0,1) (1,1)
316  { -1, -1, -1, -1, -1}, // (0,0) (1,0) (0,1) (1,1)
317 };
318 
319 double MarchingSquares::vertexList[Square::EDGES][2];
320 int MarchingSquares::GetIndex(const double v[Square::CORNERS],double iso){
321  int idx=0;
322  for(int i=0;i<int(Square::CORNERS);i++){if(v[i]<iso){idx|=(1<<i);}}
323  return idx;
324 }
325 
326 int MarchingSquares::IsAmbiguous(const double v[Square::CORNERS],double isoValue){
327  int idx=GetIndex(v,isoValue);
328  return (idx==5) || (idx==10);
329 }
330 int MarchingSquares::AddEdges(const double v[Square::CORNERS],double iso,Edge* isoEdges){
331  int idx,nEdges=0;
332  Edge e;
333 
334  idx=GetIndex(v,iso);
335 
336  /* Cube is entirely in/out of the surface */
337  if (!edgeMask[idx]) return 0;
338 
339  /* Find the vertices where the surface intersects the cube */
340  int i,j,ii=1;
341  for(i=0;i<12;i++){
342  if(edgeMask[idx] & ii){SetVertex(i,v,iso);}
343  ii<<=1;
344  }
345  /* Create the triangle */
346  for (i=0;edges[idx][i]!=-1;i+=2) {
347  for(j=0;j<2;j++){
348  e.p[0][j]=vertexList[edges[idx][i+0]][j];
349  e.p[1][j]=vertexList[edges[idx][i+1]][j];
350  }
351  isoEdges[nEdges++]=e;
352  }
353  return nEdges;
354 }
355 
356 int MarchingSquares::AddEdgeIndices(const double v[Square::CORNERS],double iso,int* isoIndices){
357  int idx,nEdges=0;
358 
359  idx=GetIndex(v,iso);
360 
361  /* Cube is entirely in/out of the surface */
362  if (!edgeMask[idx]) return 0;
363 
364  /* Create the triangle */
365  for(int i=0;edges[idx][i]!=-1;i+=2){
366  for(int j=0;j<2;j++){isoIndices[i+j]=edges[idx][i+j];}
367  nEdges++;
368  }
369  return nEdges;
370 }
371 void MarchingSquares::SetVertex(int e,const double values[Square::CORNERS],double iso){
372  int o,i,c1,c2;
373  Square::FactorEdgeIndex(e,o,i);
374  Square::EdgeCorners(e,c1,c2);
375  switch(o){
376  case 0:
377  vertexList[e][0]=Interpolate(values[c1]-iso,values[c2]-iso);
378  vertexList[e][1]=i;
379  break;
380  case 1:
381  vertexList[e][1]=Interpolate(values[c1]-iso,values[c2]-iso);
382  vertexList[e][0]=i;
383  break;
384  }
385 }
386 double MarchingSquares::Interpolate(double v1,double v2){return v1/(v1-v2);}
387 
388 
390 // MarchingCubes //
392 const int MarchingCubes::edgeMask[1<<Cube::CORNERS]={
393  0, 273, 545, 816, 2082, 2355, 2563, 2834,
394  1042, 1283, 1587, 1826, 3120, 3361, 3601, 3840,
395  324, 85, 869, 628, 2406, 2167, 2887, 2646,
396  1366, 1095, 1911, 1638, 3444, 3173, 3925, 3652,
397  644, 917, 165, 436, 2726, 2999, 2183, 2454,
398  1686, 1927, 1207, 1446, 3764, 4005, 3221, 3460,
399  960, 721, 481, 240, 3042, 2803, 2499, 2258,
400  2002, 1731, 1523, 1250, 4080, 3809, 3537, 3264,
401  2184, 2457, 2729, 3000, 170, 443, 651, 922,
402  3226, 3467, 3771, 4010, 1208, 1449, 1689, 1928,
403  2508, 2269, 3053, 2812, 494, 255, 975, 734,
404  3550, 3279, 4095, 3822, 1532, 1261, 2013, 1740,
405  2572, 2845, 2093, 2364, 558, 831, 15, 286,
406  3614, 3855, 3135, 3374, 1596, 1837, 1053, 1292,
407  2888, 2649, 2409, 2168, 874, 635, 331, 90,
408  3930, 3659, 3451, 3178, 1912, 1641, 1369, 1096,
409  1096, 1369, 1641, 1912, 3178, 3451, 3659, 3930,
410  90, 331, 635, 874, 2168, 2409, 2649, 2888,
411  1292, 1053, 1837, 1596, 3374, 3135, 3855, 3614,
412  286, 15, 831, 558, 2364, 2093, 2845, 2572,
413  1740, 2013, 1261, 1532, 3822, 4095, 3279, 3550,
414  734, 975, 255, 494, 2812, 3053, 2269, 2508,
415  1928, 1689, 1449, 1208, 4010, 3771, 3467, 3226,
416  922, 651, 443, 170, 3000, 2729, 2457, 2184,
417  3264, 3537, 3809, 4080, 1250, 1523, 1731, 2002,
418  2258, 2499, 2803, 3042, 240, 481, 721, 960,
419  3460, 3221, 4005, 3764, 1446, 1207, 1927, 1686,
420  2454, 2183, 2999, 2726, 436, 165, 917, 644,
421  3652, 3925, 3173, 3444, 1638, 1911, 1095, 1366,
422  2646, 2887, 2167, 2406, 628, 869, 85, 324,
423  3840, 3601, 3361, 3120, 1826, 1587, 1283, 1042,
424  2834, 2563, 2355, 2082, 816, 545, 273, 0
425 };
426 const int MarchingCubes::triangles[1<<Cube::CORNERS][MAX_TRIANGLES*3+1] = {
427  { -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
428  { 0, 4, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
429  { 5, 0, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
430  { 8, 9, 5, 8, 5, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
431  { 1, 5, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
432  { 0, 4, 8, 1, 5, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
433  { 9, 11, 1, 9, 1, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
434  { 8, 9, 11, 8, 11, 1, 8, 1, 4, -1, -1, -1, -1, -1, -1, -1},
435  { 4, 1, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
436  { 10, 8, 0, 10, 0, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
437  { 5, 0, 9, 4, 1, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
438  { 10, 8, 9, 10, 9, 5, 10, 5, 1, -1, -1, -1, -1, -1, -1, -1},
439  { 11, 10, 4, 11, 4, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
440  { 11, 10, 8, 11, 8, 0, 11, 0, 5, -1, -1, -1, -1, -1, -1, -1},
441  { 9, 11, 10, 9, 10, 4, 9, 4, 0, -1, -1, -1, -1, -1, -1, -1},
442  { 8, 9, 11, 8, 11, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
443  { 8, 6, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
444  { 6, 2, 0, 4, 6, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
445  { 6, 2, 8, 5, 0, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
446  { 5, 4, 6, 9, 5, 6, 2, 9, 6, -1, -1, -1, -1, -1, -1, -1},
447  { 1, 5, 11, 8, 6, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
448  { 1, 5, 11, 6, 2, 0, 4, 6, 0, -1, -1, -1, -1, -1, -1, -1},
449  { 6, 2, 8, 9, 11, 1, 9, 1, 0, -1, -1, -1, -1, -1, -1, -1},
450  { 9, 11, 2, 2, 11, 1, 2, 1, 6, 6, 1, 4, -1, -1, -1, -1},
451  { 1, 10, 4, 2, 8, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
452  { 2, 0, 1, 6, 2, 1, 10, 6, 1, -1, -1, -1, -1, -1, -1, -1},
453  { 5, 0, 9, 4, 1, 10, 8, 6, 2, -1, -1, -1, -1, -1, -1, -1},
454  { 5, 2, 9, 5, 6, 2, 5, 1, 6, 1, 10, 6, -1, -1, -1, -1},
455  { 2, 8, 6, 4, 5, 11, 4, 11, 10, -1, -1, -1, -1, -1, -1, -1},
456  { 5, 2, 0, 6, 2, 5, 11, 6, 5, 10, 6, 11, -1, -1, -1, -1},
457  { 9, 11, 10, 9, 10, 4, 9, 4, 0, 8, 6, 2, -1, -1, -1, -1},
458  { 9, 11, 2, 2, 11, 6, 10, 6, 11, -1, -1, -1, -1, -1, -1, -1},
459  { 9, 2, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
460  { 7, 9, 2, 4, 8, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
461  { 0, 2, 7, 0, 7, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
462  { 7, 5, 4, 2, 7, 4, 8, 2, 4, -1, -1, -1, -1, -1, -1, -1},
463  { 7, 9, 2, 5, 11, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
464  { 1, 5, 11, 0, 4, 8, 9, 2, 7, -1, -1, -1, -1, -1, -1, -1},
465  { 1, 0, 2, 1, 2, 7, 1, 7, 11, -1, -1, -1, -1, -1, -1, -1},
466  { 1, 7, 11, 1, 2, 7, 1, 4, 2, 4, 8, 2, -1, -1, -1, -1},
467  { 4, 1, 10, 9, 2, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
468  { 7, 9, 2, 0, 1, 10, 0, 10, 8, -1, -1, -1, -1, -1, -1, -1},
469  { 4, 1, 10, 2, 7, 5, 0, 2, 5, -1, -1, -1, -1, -1, -1, -1},
470  { 2, 10, 8, 1, 10, 2, 7, 1, 2, 5, 1, 7, -1, -1, -1, -1},
471  { 7, 9, 2, 10, 4, 5, 11, 10, 5, -1, -1, -1, -1, -1, -1, -1},
472  { 11, 10, 8, 11, 8, 0, 11, 0, 5, 9, 2, 7, -1, -1, -1, -1},
473  { 11, 10, 7, 7, 10, 4, 7, 4, 2, 2, 4, 0, -1, -1, -1, -1},
474  { 11, 10, 7, 7, 10, 2, 8, 2, 10, -1, -1, -1, -1, -1, -1, -1},
475  { 7, 9, 8, 6, 7, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
476  { 4, 6, 7, 0, 4, 7, 9, 0, 7, -1, -1, -1, -1, -1, -1, -1},
477  { 6, 7, 5, 8, 6, 5, 0, 8, 5, -1, -1, -1, -1, -1, -1, -1},
478  { 4, 6, 7, 5, 4, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
479  { 5, 11, 1, 8, 6, 7, 9, 8, 7, -1, -1, -1, -1, -1, -1, -1},
480  { 4, 6, 7, 0, 4, 7, 9, 0, 7, 11, 1, 5, -1, -1, -1, -1},
481  { 8, 1, 0, 11, 1, 8, 6, 11, 8, 7, 11, 6, -1, -1, -1, -1},
482  { 11, 6, 7, 1, 6, 11, 6, 1, 4, -1, -1, -1, -1, -1, -1, -1},
483  { 1, 10, 4, 6, 7, 9, 6, 9, 8, -1, -1, -1, -1, -1, -1, -1},
484  { 0, 1, 9, 9, 1, 10, 9, 10, 7, 7, 10, 6, -1, -1, -1, -1},
485  { 6, 7, 5, 8, 6, 5, 0, 8, 5, 1, 10, 4, -1, -1, -1, -1},
486  { 1, 7, 5, 10, 7, 1, 7, 10, 6, -1, -1, -1, -1, -1, -1, -1},
487  { 11, 10, 4, 11, 4, 5, 7, 9, 8, 6, 7, 8, -1, -1, -1, -1},
488  { 0, 6, 9, 9, 6, 7, 6, 0, 5, 5, 11, 10, 5, 10, 6, -1},
489  { 8, 7, 0, 6, 7, 8, 4, 0, 7, 11, 10, 4, 7, 11, 4, -1},
490  { 11, 10, 6, 11, 6, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
491  { 11, 7, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
492  { 0, 4, 8, 11, 7, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
493  { 9, 5, 0, 11, 7, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
494  { 11, 7, 3, 4, 8, 9, 5, 4, 9, -1, -1, -1, -1, -1, -1, -1},
495  { 3, 1, 5, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
496  { 0, 4, 8, 7, 3, 1, 5, 7, 1, -1, -1, -1, -1, -1, -1, -1},
497  { 3, 1, 0, 3, 0, 9, 3, 9, 7, -1, -1, -1, -1, -1, -1, -1},
498  { 7, 8, 9, 4, 8, 7, 3, 4, 7, 1, 4, 3, -1, -1, -1, -1},
499  { 1, 10, 4, 3, 11, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
500  { 3, 11, 7, 8, 0, 1, 10, 8, 1, -1, -1, -1, -1, -1, -1, -1},
501  { 4, 1, 10, 5, 0, 9, 11, 7, 3, -1, -1, -1, -1, -1, -1, -1},
502  { 10, 8, 9, 10, 9, 5, 10, 5, 1, 11, 7, 3, -1, -1, -1, -1},
503  { 4, 5, 7, 4, 7, 3, 4, 3, 10, -1, -1, -1, -1, -1, -1, -1},
504  { 10, 8, 3, 3, 8, 0, 3, 0, 7, 7, 0, 5, -1, -1, -1, -1},
505  { 4, 3, 10, 4, 7, 3, 4, 0, 7, 0, 9, 7, -1, -1, -1, -1},
506  { 10, 8, 3, 3, 8, 7, 9, 7, 8, -1, -1, -1, -1, -1, -1, -1},
507  { 11, 7, 3, 8, 6, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
508  { 11, 7, 3, 2, 0, 4, 2, 4, 6, -1, -1, -1, -1, -1, -1, -1},
509  { 11, 7, 3, 8, 6, 2, 5, 0, 9, -1, -1, -1, -1, -1, -1, -1},
510  { 5, 4, 6, 9, 5, 6, 2, 9, 6, 3, 11, 7, -1, -1, -1, -1},
511  { 8, 6, 2, 3, 1, 5, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1},
512  { 3, 1, 5, 3, 5, 7, 6, 2, 0, 4, 6, 0, -1, -1, -1, -1},
513  { 3, 1, 0, 3, 0, 9, 3, 9, 7, 2, 8, 6, -1, -1, -1, -1},
514  { 9, 4, 2, 2, 4, 6, 4, 9, 7, 7, 3, 1, 7, 1, 4, -1},
515  { 8, 6, 2, 11, 7, 3, 4, 1, 10, -1, -1, -1, -1, -1, -1, -1},
516  { 2, 0, 1, 6, 2, 1, 10, 6, 1, 11, 7, 3, -1, -1, -1, -1},
517  { 5, 0, 9, 4, 1, 10, 8, 6, 2, 11, 7, 3, -1, -1, -1, -1},
518  { 11, 7, 3, 5, 2, 9, 5, 6, 2, 5, 1, 6, 1, 10, 6, -1},
519  { 4, 5, 7, 4, 7, 3, 4, 3, 10, 6, 2, 8, -1, -1, -1, -1},
520  { 10, 5, 3, 3, 5, 7, 5, 10, 6, 6, 2, 0, 6, 0, 5, -1},
521  { 8, 6, 2, 4, 3, 10, 4, 7, 3, 4, 0, 7, 0, 9, 7, -1},
522  { 9, 7, 10, 10, 7, 3, 10, 6, 9, 6, 2, 9, -1, -1, -1, -1},
523  { 3, 11, 9, 2, 3, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
524  { 4, 8, 0, 2, 3, 11, 2, 11, 9, -1, -1, -1, -1, -1, -1, -1},
525  { 0, 2, 3, 0, 3, 11, 0, 11, 5, -1, -1, -1, -1, -1, -1, -1},
526  { 2, 3, 8, 8, 3, 11, 8, 11, 4, 4, 11, 5, -1, -1, -1, -1},
527  { 2, 3, 1, 2, 1, 5, 2, 5, 9, -1, -1, -1, -1, -1, -1, -1},
528  { 2, 3, 1, 2, 1, 5, 2, 5, 9, 0, 4, 8, -1, -1, -1, -1},
529  { 0, 2, 3, 0, 3, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
530  { 2, 3, 8, 8, 3, 4, 1, 4, 3, -1, -1, -1, -1, -1, -1, -1},
531  { 1, 10, 4, 9, 2, 3, 11, 9, 3, -1, -1, -1, -1, -1, -1, -1},
532  { 10, 8, 0, 10, 0, 1, 3, 11, 9, 2, 3, 9, -1, -1, -1, -1},
533  { 0, 2, 3, 0, 3, 11, 0, 11, 5, 1, 10, 4, -1, -1, -1, -1},
534  { 5, 2, 11, 11, 2, 3, 2, 5, 1, 1, 10, 8, 1, 8, 2, -1},
535  { 10, 2, 3, 9, 2, 10, 4, 9, 10, 5, 9, 4, -1, -1, -1, -1},
536  { 5, 10, 0, 0, 10, 8, 10, 5, 9, 9, 2, 3, 9, 3, 10, -1},
537  { 0, 2, 4, 4, 2, 10, 3, 10, 2, -1, -1, -1, -1, -1, -1, -1},
538  { 10, 8, 2, 10, 2, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
539  { 11, 9, 8, 3, 11, 8, 6, 3, 8, -1, -1, -1, -1, -1, -1, -1},
540  { 0, 11, 9, 3, 11, 0, 4, 3, 0, 6, 3, 4, -1, -1, -1, -1},
541  { 11, 5, 3, 5, 0, 3, 0, 6, 3, 0, 8, 6, -1, -1, -1, -1},
542  { 3, 4, 6, 11, 4, 3, 4, 11, 5, -1, -1, -1, -1, -1, -1, -1},
543  { 3, 1, 6, 6, 1, 5, 6, 5, 8, 8, 5, 9, -1, -1, -1, -1},
544  { 0, 6, 9, 4, 6, 0, 5, 9, 6, 3, 1, 5, 6, 3, 5, -1},
545  { 3, 1, 6, 6, 1, 8, 0, 8, 1, -1, -1, -1, -1, -1, -1, -1},
546  { 3, 1, 4, 3, 4, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
547  { 11, 9, 8, 3, 11, 8, 6, 3, 8, 4, 1, 10, -1, -1, -1, -1},
548  { 3, 9, 6, 11, 9, 3, 10, 6, 9, 0, 1, 10, 9, 0, 10, -1},
549  { 4, 1, 10, 11, 5, 3, 5, 0, 3, 0, 6, 3, 0, 8, 6, -1},
550  { 5, 10, 6, 1, 10, 5, 6, 11, 5, 6, 3, 11, -1, -1, -1, -1},
551  { 10, 5, 3, 4, 5, 10, 6, 3, 5, 9, 8, 6, 5, 9, 6, -1},
552  { 6, 3, 10, 9, 0, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
553  { 3, 10, 0, 0, 10, 4, 0, 8, 3, 8, 6, 3, -1, -1, -1, -1},
554  { 6, 3, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
555  { 10, 3, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
556  { 3, 6, 10, 0, 4, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
557  { 5, 0, 9, 10, 3, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
558  { 3, 6, 10, 8, 9, 5, 8, 5, 4, -1, -1, -1, -1, -1, -1, -1},
559  { 11, 1, 5, 10, 3, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
560  { 0, 4, 8, 1, 5, 11, 10, 3, 6, -1, -1, -1, -1, -1, -1, -1},
561  { 10, 3, 6, 0, 9, 11, 1, 0, 11, -1, -1, -1, -1, -1, -1, -1},
562  { 8, 9, 11, 8, 11, 1, 8, 1, 4, 10, 3, 6, -1, -1, -1, -1},
563  { 4, 1, 3, 6, 4, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
564  { 0, 1, 3, 8, 0, 3, 6, 8, 3, -1, -1, -1, -1, -1, -1, -1},
565  { 5, 0, 9, 3, 6, 4, 1, 3, 4, -1, -1, -1, -1, -1, -1, -1},
566  { 8, 9, 6, 6, 9, 5, 6, 5, 3, 3, 5, 1, -1, -1, -1, -1},
567  { 6, 4, 5, 6, 5, 11, 6, 11, 3, -1, -1, -1, -1, -1, -1, -1},
568  { 0, 6, 8, 0, 3, 6, 0, 5, 3, 5, 11, 3, -1, -1, -1, -1},
569  { 3, 9, 11, 0, 9, 3, 6, 0, 3, 4, 0, 6, -1, -1, -1, -1},
570  { 8, 9, 6, 6, 9, 3, 11, 3, 9, -1, -1, -1, -1, -1, -1, -1},
571  { 2, 8, 10, 3, 2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
572  { 3, 2, 0, 10, 3, 0, 4, 10, 0, -1, -1, -1, -1, -1, -1, -1},
573  { 5, 0, 9, 8, 10, 3, 8, 3, 2, -1, -1, -1, -1, -1, -1, -1},
574  { 9, 3, 2, 10, 3, 9, 5, 10, 9, 4, 10, 5, -1, -1, -1, -1},
575  { 11, 1, 5, 2, 8, 10, 3, 2, 10, -1, -1, -1, -1, -1, -1, -1},
576  { 3, 2, 0, 10, 3, 0, 4, 10, 0, 5, 11, 1, -1, -1, -1, -1},
577  { 9, 11, 1, 9, 1, 0, 2, 8, 10, 3, 2, 10, -1, -1, -1, -1},
578  { 10, 2, 4, 3, 2, 10, 1, 4, 2, 9, 11, 1, 2, 9, 1, -1},
579  { 1, 3, 2, 4, 1, 2, 8, 4, 2, -1, -1, -1, -1, -1, -1, -1},
580  { 0, 1, 3, 2, 0, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
581  { 1, 3, 2, 4, 1, 2, 8, 4, 2, 9, 5, 0, -1, -1, -1, -1},
582  { 9, 3, 2, 5, 3, 9, 3, 5, 1, -1, -1, -1, -1, -1, -1, -1},
583  { 3, 2, 11, 11, 2, 8, 11, 8, 5, 5, 8, 4, -1, -1, -1, -1},
584  { 5, 2, 0, 11, 2, 5, 2, 11, 3, -1, -1, -1, -1, -1, -1, -1},
585  { 4, 3, 8, 8, 3, 2, 3, 4, 0, 0, 9, 11, 0, 11, 3, -1},
586  { 9, 11, 3, 9, 3, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
587  { 10, 3, 6, 9, 2, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
588  { 9, 2, 7, 10, 3, 6, 0, 4, 8, -1, -1, -1, -1, -1, -1, -1},
589  { 10, 3, 6, 7, 5, 0, 7, 0, 2, -1, -1, -1, -1, -1, -1, -1},
590  { 7, 5, 4, 2, 7, 4, 8, 2, 4, 10, 3, 6, -1, -1, -1, -1},
591  { 10, 3, 6, 9, 2, 7, 1, 5, 11, -1, -1, -1, -1, -1, -1, -1},
592  { 10, 3, 6, 9, 2, 7, 1, 5, 11, 0, 4, 8, -1, -1, -1, -1},
593  { 1, 0, 2, 1, 2, 7, 1, 7, 11, 3, 6, 10, -1, -1, -1, -1},
594  { 10, 3, 6, 1, 7, 11, 1, 2, 7, 1, 4, 2, 4, 8, 2, -1},
595  { 9, 2, 7, 6, 4, 1, 6, 1, 3, -1, -1, -1, -1, -1, -1, -1},
596  { 0, 1, 3, 8, 0, 3, 6, 8, 3, 7, 9, 2, -1, -1, -1, -1},
597  { 0, 2, 7, 0, 7, 5, 4, 1, 3, 6, 4, 3, -1, -1, -1, -1},
598  { 2, 5, 8, 7, 5, 2, 6, 8, 5, 1, 3, 6, 5, 1, 6, -1},
599  { 6, 4, 5, 6, 5, 11, 6, 11, 3, 7, 9, 2, -1, -1, -1, -1},
600  { 9, 2, 7, 0, 6, 8, 0, 3, 6, 0, 5, 3, 5, 11, 3, -1},
601  { 3, 4, 11, 6, 4, 3, 7, 11, 4, 0, 2, 7, 4, 0, 7, -1},
602  { 11, 3, 8, 8, 3, 6, 8, 2, 11, 2, 7, 11, -1, -1, -1, -1},
603  { 9, 8, 10, 7, 9, 10, 3, 7, 10, -1, -1, -1, -1, -1, -1, -1},
604  { 9, 0, 7, 0, 4, 7, 4, 3, 7, 4, 10, 3, -1, -1, -1, -1},
605  { 8, 10, 0, 0, 10, 3, 0, 3, 5, 5, 3, 7, -1, -1, -1, -1},
606  { 10, 5, 4, 3, 5, 10, 5, 3, 7, -1, -1, -1, -1, -1, -1, -1},
607  { 9, 8, 10, 7, 9, 10, 3, 7, 10, 1, 5, 11, -1, -1, -1, -1},
608  { 1, 5, 11, 9, 0, 7, 0, 4, 7, 4, 3, 7, 4, 10, 3, -1},
609  { 11, 0, 7, 1, 0, 11, 3, 7, 0, 8, 10, 3, 0, 8, 3, -1},
610  { 7, 1, 4, 11, 1, 7, 4, 3, 7, 4, 10, 3, -1, -1, -1, -1},
611  { 4, 9, 8, 7, 9, 4, 1, 7, 4, 3, 7, 1, -1, -1, -1, -1},
612  { 7, 1, 3, 9, 1, 7, 1, 9, 0, -1, -1, -1, -1, -1, -1, -1},
613  { 8, 7, 0, 0, 7, 5, 7, 8, 4, 4, 1, 3, 4, 3, 7, -1},
614  { 5, 1, 3, 7, 5, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
615  { 3, 4, 11, 11, 4, 5, 4, 3, 7, 7, 9, 8, 7, 8, 4, -1},
616  { 3, 9, 0, 7, 9, 3, 0, 11, 3, 0, 5, 11, -1, -1, -1, -1},
617  { 3, 7, 11, 8, 4, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
618  { 3, 7, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
619  { 6, 10, 11, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
620  { 0, 4, 8, 10, 11, 7, 10, 7, 6, -1, -1, -1, -1, -1, -1, -1},
621  { 9, 5, 0, 6, 10, 11, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1},
622  { 8, 9, 5, 8, 5, 4, 6, 10, 11, 7, 6, 11, -1, -1, -1, -1},
623  { 5, 7, 6, 5, 6, 10, 5, 10, 1, -1, -1, -1, -1, -1, -1, -1},
624  { 5, 7, 6, 5, 6, 10, 5, 10, 1, 4, 8, 0, -1, -1, -1, -1},
625  { 1, 0, 10, 10, 0, 9, 10, 9, 6, 6, 9, 7, -1, -1, -1, -1},
626  { 1, 7, 10, 10, 7, 6, 7, 1, 4, 4, 8, 9, 4, 9, 7, -1},
627  { 7, 6, 4, 7, 4, 1, 7, 1, 11, -1, -1, -1, -1, -1, -1, -1},
628  { 11, 0, 1, 8, 0, 11, 7, 8, 11, 6, 8, 7, -1, -1, -1, -1},
629  { 7, 6, 4, 7, 4, 1, 7, 1, 11, 5, 0, 9, -1, -1, -1, -1},
630  { 11, 6, 1, 7, 6, 11, 5, 1, 6, 8, 9, 5, 6, 8, 5, -1},
631  { 4, 5, 7, 4, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
632  { 5, 7, 0, 0, 7, 8, 6, 8, 7, -1, -1, -1, -1, -1, -1, -1},
633  { 7, 6, 9, 9, 6, 0, 4, 0, 6, -1, -1, -1, -1, -1, -1, -1},
634  { 8, 9, 7, 8, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
635  { 8, 10, 11, 2, 8, 11, 7, 2, 11, -1, -1, -1, -1, -1, -1, -1},
636  { 10, 11, 4, 4, 11, 7, 4, 7, 0, 0, 7, 2, -1, -1, -1, -1},
637  { 8, 10, 11, 2, 8, 11, 7, 2, 11, 5, 0, 9, -1, -1, -1, -1},
638  { 9, 4, 2, 5, 4, 9, 7, 2, 4, 10, 11, 7, 4, 10, 7, -1},
639  { 1, 8, 10, 2, 8, 1, 5, 2, 1, 7, 2, 5, -1, -1, -1, -1},
640  { 1, 7, 10, 5, 7, 1, 4, 10, 7, 2, 0, 4, 7, 2, 4, -1},
641  { 7, 1, 9, 9, 1, 0, 1, 7, 2, 2, 8, 10, 2, 10, 1, -1},
642  { 7, 2, 9, 10, 1, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
643  { 8, 4, 2, 4, 1, 2, 1, 7, 2, 1, 11, 7, -1, -1, -1, -1},
644  { 11, 0, 1, 7, 0, 11, 0, 7, 2, -1, -1, -1, -1, -1, -1, -1},
645  { 5, 0, 9, 8, 4, 2, 4, 1, 2, 1, 7, 2, 1, 11, 7, -1},
646  { 2, 5, 1, 9, 5, 2, 1, 7, 2, 1, 11, 7, -1, -1, -1, -1},
647  { 4, 5, 8, 8, 5, 2, 7, 2, 5, -1, -1, -1, -1, -1, -1, -1},
648  { 7, 2, 0, 5, 7, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
649  { 7, 2, 4, 4, 2, 8, 4, 0, 7, 0, 9, 7, -1, -1, -1, -1},
650  { 7, 2, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
651  { 10, 11, 9, 6, 10, 9, 2, 6, 9, -1, -1, -1, -1, -1, -1, -1},
652  { 10, 11, 9, 6, 10, 9, 2, 6, 9, 0, 4, 8, -1, -1, -1, -1},
653  { 5, 10, 11, 6, 10, 5, 0, 6, 5, 2, 6, 0, -1, -1, -1, -1},
654  { 2, 5, 8, 8, 5, 4, 5, 2, 6, 6, 10, 11, 6, 11, 5, -1},
655  { 10, 1, 6, 1, 5, 6, 5, 2, 6, 5, 9, 2, -1, -1, -1, -1},
656  { 0, 4, 8, 10, 1, 6, 1, 5, 6, 5, 2, 6, 5, 9, 2, -1},
657  { 1, 0, 10, 10, 0, 6, 2, 6, 0, -1, -1, -1, -1, -1, -1, -1},
658  { 2, 6, 1, 1, 6, 10, 1, 4, 2, 4, 8, 2, -1, -1, -1, -1},
659  { 11, 9, 1, 1, 9, 2, 1, 2, 4, 4, 2, 6, -1, -1, -1, -1},
660  { 8, 1, 6, 0, 1, 8, 2, 6, 1, 11, 9, 2, 1, 11, 2, -1},
661  { 11, 6, 1, 1, 6, 4, 6, 11, 5, 5, 0, 2, 5, 2, 6, -1},
662  { 2, 6, 8, 11, 5, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
663  { 6, 4, 2, 2, 4, 9, 5, 9, 4, -1, -1, -1, -1, -1, -1, -1},
664  { 5, 9, 6, 6, 9, 2, 6, 8, 5, 8, 0, 5, -1, -1, -1, -1},
665  { 0, 2, 6, 0, 6, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
666  { 2, 6, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
667  { 8, 10, 11, 9, 8, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
668  { 0, 11, 9, 4, 11, 0, 11, 4, 10, -1, -1, -1, -1, -1, -1, -1},
669  { 5, 10, 11, 0, 10, 5, 10, 0, 8, -1, -1, -1, -1, -1, -1, -1},
670  { 4, 10, 11, 5, 4, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
671  { 1, 8, 10, 5, 8, 1, 8, 5, 9, -1, -1, -1, -1, -1, -1, -1},
672  { 9, 4, 10, 0, 4, 9, 10, 5, 9, 10, 1, 5, -1, -1, -1, -1},
673  { 0, 8, 10, 1, 0, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
674  { 10, 1, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
675  { 4, 9, 8, 1, 9, 4, 9, 1, 11, -1, -1, -1, -1, -1, -1, -1},
676  { 1, 11, 9, 0, 1, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
677  { 11, 0, 8, 5, 0, 11, 8, 1, 11, 8, 4, 1, -1, -1, -1, -1},
678  { 11, 5, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
679  { 5, 9, 8, 4, 5, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
680  { 9, 0, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
681  { 8, 4, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1},
682  { -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}
683 };
684 const int MarchingCubes::cornerMap[Cube::CORNERS]={0,1,3,2,4,5,7,6};
685 double MarchingCubes::vertexList[Cube::EDGES][3];
686 
687 int MarchingCubes::GetIndex(const double v[Cube::CORNERS],double iso){
688  int idx=0;
689  if (v[Cube::CornerIndex(0,0,0)] < iso) idx |= 1;
690  if (v[Cube::CornerIndex(1,0,0)] < iso) idx |= 2;
691  if (v[Cube::CornerIndex(1,1,0)] < iso) idx |= 4;
692  if (v[Cube::CornerIndex(0,1,0)] < iso) idx |= 8;
693  if (v[Cube::CornerIndex(0,0,1)] < iso) idx |= 16;
694  if (v[Cube::CornerIndex(1,0,1)] < iso) idx |= 32;
695  if (v[Cube::CornerIndex(1,1,1)] < iso) idx |= 64;
696  if (v[Cube::CornerIndex(0,1,1)] < iso) idx |= 128;
697  return idx;
698 }
699 int MarchingCubes::GetFaceIndex(const double values[Cube::CORNERS],double iso,int faceIndex){
700  int i,j,x,y,z,idx=0;
701  double v[2][2];
702  Cube::FactorFaceIndex(faceIndex,x,y,z);
703  if (x<0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(0,i,j)];}}}
704  else if (x>0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(1,i,j)];}}}
705  else if (y<0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(i,0,j)];}}}
706  else if (y>0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(i,1,j)];}}}
707  else if (z<0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(i,j,0)];}}}
708  else if (z>0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(i,j,1)];}}}
709  if (v[0][0] < iso) idx |= 1;
710  if (v[1][0] < iso) idx |= 2;
711  if (v[1][1] < iso) idx |= 4;
712  if (v[0][1] < iso) idx |= 8;
713  return idx;
714 }
715 int MarchingCubes::IsAmbiguous(const double v[Cube::CORNERS],double isoValue,int faceIndex){
716  int idx=GetFaceIndex(v,isoValue,faceIndex);
717  return (idx==5) || (idx==10);
718 }
719 int MarchingCubes::HasRoots(const double v[Cube::CORNERS],double isoValue,int faceIndex){
720  int idx=GetFaceIndex(v,isoValue,faceIndex);
721  return (idx!=0) && (idx !=15);
722 }
723 int MarchingCubes::HasRoots(const double v[Cube::CORNERS],double isoValue){
724  int idx=GetIndex(v,isoValue);
725  if(idx==0 || idx==255){return 0;}
726  else{return 1;}
727 }
728 int MarchingCubes::HasRoots(int mcIndex){
729  if(mcIndex==0 || mcIndex==255){return 0;}
730  else{return 1;}
731 }
732 int MarchingCubes::AddTriangles( const double v[Cube::CORNERS] , double iso , Triangle* isoTriangles )
733 {
734  int idx,ntriang=0;
735  Triangle tri;
736 
737  idx=GetIndex(v,iso);
738 
739  /* Cube is entirely in/out of the surface */
740  if (!edgeMask[idx]) return 0;
741 
742  /* Find the vertices where the surface intersects the cube */
743  int i,j,ii=1;
744  for(i=0;i<12;i++){
745  if(edgeMask[idx] & ii){SetVertex(i,v,iso);}
746  ii<<=1;
747  }
748  /* Create the triangle */
749  for( i=0 ; triangles[idx][i]!=-1 ; i+=3 )
750  {
751  for(j=0;j<3;j++){
752  tri.p[0][j]=vertexList[triangles[idx][i+0]][j];
753  tri.p[1][j]=vertexList[triangles[idx][i+1]][j];
754  tri.p[2][j]=vertexList[triangles[idx][i+2]][j];
755  }
756  isoTriangles[ntriang++]=tri;
757  }
758  return ntriang;
759 }
760 
761 int MarchingCubes::AddTriangleIndices(const double v[Cube::CORNERS],double iso,int* isoIndices){
762  int idx,ntriang=0;
763 
764  idx=GetIndex(v,iso);
765 
766  /* Cube is entirely in/out of the surface */
767  if (!edgeMask[idx]) return 0;
768 
769  /* Create the triangle */
770  for(int i=0;triangles[idx][i]!=-1;i+=3){
771  for(int j=0;j<3;j++){isoIndices[i+j]=triangles[idx][i+j];}
772  ntriang++;
773  }
774  return ntriang;
775 }
776 
777 void MarchingCubes::SetVertex( int e , const double values[Cube::CORNERS] , double iso )
778 {
779  double t;
780  int o , i1 , i2;
781  Cube::FactorEdgeIndex( e , o , i1 , i2 );
782  switch( o )
783  {
784  case 0:
785  t = Interpolate( values[ Cube::CornerIndex( 0 , i1 , i2 ) ] - iso , values[ Cube::CornerIndex( 1 , i1 , i2 ) ] - iso );
786  vertexList[e][0] = t , vertexList[e][1] = i1 , vertexList[e][2] = i2;
787  break;
788  case 1:
789  t = Interpolate( values[ Cube::CornerIndex( i1 , 0 , i2 ) ] - iso , values[ Cube::CornerIndex( i1 , 1 , i2 ) ] - iso );
790  vertexList[e][0] = i1 , vertexList[e][1] = t , vertexList[e][2] = i2;
791  break;
792  case 2:
793  t = Interpolate( values[ Cube::CornerIndex( i1 , i2 , 0 ) ] - iso , values[ Cube::CornerIndex( i1 , i2 , 1 ) ] - iso );
794  vertexList[e][0] = i1 , vertexList[e][1] = i2 , vertexList[e][2] = t;
795  break;
796  }
797 }
798 double MarchingCubes::Interpolate( double v1 , double v2 ) { return v1/(v1-v2); }
799 
800 
802 int MarchingCubes::GetIndex(const float v[Cube::CORNERS],float iso){
803  int idx=0;
804  if (v[Cube::CornerIndex(0,0,0)] < iso) idx |= 1;
805  if (v[Cube::CornerIndex(1,0,0)] < iso) idx |= 2;
806  if (v[Cube::CornerIndex(1,1,0)] < iso) idx |= 4;
807  if (v[Cube::CornerIndex(0,1,0)] < iso) idx |= 8;
808  if (v[Cube::CornerIndex(0,0,1)] < iso) idx |= 16;
809  if (v[Cube::CornerIndex(1,0,1)] < iso) idx |= 32;
810  if (v[Cube::CornerIndex(1,1,1)] < iso) idx |= 64;
811  if (v[Cube::CornerIndex(0,1,1)] < iso) idx |= 128;
812  return idx;
813 }
814 int MarchingCubes::GetFaceIndex(const float values[Cube::CORNERS],float iso,int faceIndex){
815  int i,j,x,y,z,idx=0;
816  double v[2][2];
817  Cube::FactorFaceIndex(faceIndex,x,y,z);
818  if (x<0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(0,i,j)];}}}
819  else if (x>0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(1,i,j)];}}}
820  else if (y<0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(i,0,j)];}}}
821  else if (y>0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(i,1,j)];}}}
822  else if (z<0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(i,j,0)];}}}
823  else if (z>0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=values[Cube::CornerIndex(i,j,1)];}}}
824  if (v[0][0] < iso) idx |= 1;
825  if (v[1][0] < iso) idx |= 2;
826  if (v[1][1] < iso) idx |= 4;
827  if (v[0][1] < iso) idx |= 8;
828  return idx;
829 }
830 int MarchingCubes::GetFaceIndex(int mcIndex,int faceIndex){
831  int i,j,x,y,z,idx=0;
832  int v[2][2];
833  Cube::FactorFaceIndex(faceIndex,x,y,z);
834  if (x<0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=mcIndex&(1<<MarchingCubes::cornerMap[Cube::CornerIndex(0,i,j)]);}}}
835  else if (x>0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=mcIndex&(1<<MarchingCubes::cornerMap[Cube::CornerIndex(1,i,j)]);}}}
836  else if (y<0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=mcIndex&(1<<MarchingCubes::cornerMap[Cube::CornerIndex(i,0,j)]);}}}
837  else if (y>0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=mcIndex&(1<<MarchingCubes::cornerMap[Cube::CornerIndex(i,1,j)]);}}}
838  else if (z<0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=mcIndex&(1<<MarchingCubes::cornerMap[Cube::CornerIndex(i,j,1)]);}}}
839  else if (z>0){for(i=0;i<2;i++){for(j=0;j<2;j++){v[i][j]=mcIndex&(1<<MarchingCubes::cornerMap[Cube::CornerIndex(i,j,1)]);}}}
840  if (v[0][0]) idx |= 1;
841  if (v[1][0]) idx |= 2;
842  if (v[1][1]) idx |= 4;
843  if (v[0][1]) idx |= 8;
844  return idx;
845 }
846 int MarchingCubes::IsAmbiguous(const float v[Cube::CORNERS],float isoValue,int faceIndex){
847  int idx=GetFaceIndex(v,isoValue,faceIndex);
848  return (idx==5) || (idx==10);
849 }
850 int MarchingCubes::IsAmbiguous(int mcIndex,int faceIndex){
851  int idx=GetFaceIndex(mcIndex,faceIndex);
852  return (idx==5) || (idx==10);
853 }
854 int MarchingCubes::HasRoots(const float v[Cube::CORNERS],float isoValue){
855  int idx=GetIndex(v,isoValue);
856  if(idx==0 || idx==255){return 0;}
857  else{return 1;}
858 }
859 int MarchingCubes::HasRoots(const float v[Cube::CORNERS],float isoValue,int faceIndex){
860  int idx=GetFaceIndex(v,isoValue,faceIndex);
861  return (idx!=0) && (idx!=15);
862 }
863 int MarchingCubes::HasFaceRoots(int mcIndex,int faceIndex){
864  int idx=GetFaceIndex(mcIndex,faceIndex);
865  return (idx!=0) && (idx!=15);
866 }
867 int MarchingCubes::HasEdgeRoots(int mcIndex,int edgeIndex){
868  int c1,c2;
869  Cube::EdgeCorners(edgeIndex,c1,c2);
870  if( ( (mcIndex&(1<<MarchingCubes::cornerMap[c1])) && (mcIndex&(1<<MarchingCubes::cornerMap[c2]))) ||
871  (!(mcIndex&(1<<MarchingCubes::cornerMap[c1])) && !(mcIndex&(1<<MarchingCubes::cornerMap[c2])))){return 0;}
872  else{return 1;}
873 }
874 int MarchingCubes::AddTriangles(const float v[Cube::CORNERS],float iso,Triangle* isoTriangles){
875  int idx,ntriang=0;
876  Triangle tri;
877 
878  idx=GetIndex(v,iso);
879 
880  /* Cube is entirely in/out of the surface */
881  if (!edgeMask[idx]) return 0;
882 
883  /* Find the vertices where the surface intersects the cube */
884  int i,j,ii=1;
885  for( i=0 ; i<12 ; i++ )
886  {
887  if(edgeMask[idx] & ii) SetVertex( i , v , iso );
888  ii<<=1;
889  }
890  /* Create the triangle */
891  for (i=0;triangles[idx][i]!=-1;i+=3) {
892  for(j=0;j<3;j++){
893  tri.p[0][j]=vertexList[triangles[idx][i+0]][j];
894  tri.p[1][j]=vertexList[triangles[idx][i+1]][j];
895  tri.p[2][j]=vertexList[triangles[idx][i+2]][j];
896  }
897  isoTriangles[ntriang++]=tri;
898  }
899  return ntriang;
900 }
901 
902 int MarchingCubes::AddTriangleIndices( const float v[Cube::CORNERS] , float iso , int* isoIndices )
903 {
904  int idx,ntriang=0;
905  idx=GetIndex(v,iso);
906  /* Cube is entirely in/out of the surface */
907  if (!edgeMask[idx]) return 0;
908  /* Create the triangle */
909  for(int i=0;triangles[idx][i]!=-1;i+=3){
910  for(int j=0;j<3;j++){isoIndices[i+j]=triangles[idx][i+j];}
911  ntriang++;
912  }
913  return ntriang;
914 }
915 int MarchingCubes::AddTriangleIndices( int idx , int* isoIndices )
916 {
917  int ntriang=0;
918 
919  /* Cube is entirely in/out of the surface */
920  if (!edgeMask[idx]) return 0;
921 
922  /* Create the triangle */
923  for(int i=0;triangles[idx][i]!=-1;i+=3){
924  for(int j=0;j<3;j++){isoIndices[i+j]=triangles[idx][i+j];}
925  ntriang++;
926  }
927  return ntriang;
928 }
929 
930 void MarchingCubes::SetVertex( int e , const float values[Cube::CORNERS] , float iso )
931 {
932  double t;
933  int o , i1 , i2;
934  Cube::FactorEdgeIndex( e , o , i1 , i2 );
935  switch( o )
936  {
937  case 0:
938  t = Interpolate( values[ Cube::CornerIndex( 0 , i1 , i2 ) ] - iso , values[ Cube::CornerIndex( 1 , i1 , i2 ) ] - iso );
939  vertexList[e][0] = t , vertexList[e][1] = i1 , vertexList[e][2] = i2;
940  break;
941  case 1:
942  t = Interpolate( values[ Cube::CornerIndex( i1 , 0 , i2 ) ] - iso , values[ Cube::CornerIndex( i1 , 1 , i2 ) ] - iso );
943  vertexList[e][0] = i1 , vertexList[e][1] = t , vertexList[e][2] = i2;
944  break;
945  case 2:
946  t = Interpolate( values[ Cube::CornerIndex( i1 , i2 , 0 ) ] - iso , values[ Cube::CornerIndex( i1 , i2 , 1 ) ] - iso );
947  vertexList[e][0] = i1 , vertexList[e][1] = i2 , vertexList[e][2] = t;
948  break;
949  }
950 }
951 float MarchingCubes::Interpolate( float v1 , float v2 ){ return v1/(v1-v2); }
Definition: Geometry.h:184