iMSTK
Interactive Medical Simulation Toolkit
imstkLocalMarchingCubes.cpp
1 /*
2 ** This file is part of the Interactive Medical Simulation Toolkit (iMSTK)
3 ** iMSTK is distributed under the Apache License, Version 2.0.
4 ** See accompanying NOTICE for details.
5 */
6 
7 #include "imstkLocalMarchingCubes.h"
8 #include "imstkGeometryUtilities.h"
9 #include "imstkLogger.h"
10 #include "imstkSurfaceMesh.h"
11 #include "imstkImageData.h"
12 #include "imstkVecDataArray.h"
13 
14 #include <vtkImageData.h>
15 #include <vtkFlyingEdges3D.h>
16 
17 namespace imstk
18 {
19 // v7_______e6_____________v6
20 // /| /|
21 // / | / |
22 // e7/ | e5/ |
23 // /___|______e4_________/ |
24 // v4| | |v5 |e10
25 // | | | |
26 // | |e11 |e9 |
27 // e8| | | |
28 // | |_________________|___| +
29 // | / v3 e2 | /v2 | +
30 // | / | / | /
31 // | /e3 | /e1 |/
32 // |/_____________________|/ |------ +
33 // v0 e0 v1
34 
38 static int edgeTable[256] =
39 {
40  0x0, 0x109, 0x203, 0x30a, 0x406, 0x50f, 0x605, 0x70c,
41  0x80c, 0x905, 0xa0f, 0xb06, 0xc0a, 0xd03, 0xe09, 0xf00,
42  0x190, 0x99, 0x393, 0x29a, 0x596, 0x49f, 0x795, 0x69c,
43  0x99c, 0x895, 0xb9f, 0xa96, 0xd9a, 0xc93, 0xf99, 0xe90,
44  0x230, 0x339, 0x33, 0x13a, 0x636, 0x73f, 0x435, 0x53c,
45  0xa3c, 0xb35, 0x83f, 0x936, 0xe3a, 0xf33, 0xc39, 0xd30,
46  0x3a0, 0x2a9, 0x1a3, 0xaa, 0x7a6, 0x6af, 0x5a5, 0x4ac,
47  0xbac, 0xaa5, 0x9af, 0x8a6, 0xfaa, 0xea3, 0xda9, 0xca0,
48  0x460, 0x569, 0x663, 0x76a, 0x66, 0x16f, 0x265, 0x36c,
49  0xc6c, 0xd65, 0xe6f, 0xf66, 0x86a, 0x963, 0xa69, 0xb60,
50  0x5f0, 0x4f9, 0x7f3, 0x6fa, 0x1f6, 0xff, 0x3f5, 0x2fc,
51  0xdfc, 0xcf5, 0xfff, 0xef6, 0x9fa, 0x8f3, 0xbf9, 0xaf0,
52  0x650, 0x759, 0x453, 0x55a, 0x256, 0x35f, 0x55, 0x15c,
53  0xe5c, 0xf55, 0xc5f, 0xd56, 0xa5a, 0xb53, 0x859, 0x950,
54  0x7c0, 0x6c9, 0x5c3, 0x4ca, 0x3c6, 0x2cf, 0x1c5, 0xcc,
55  0xfcc, 0xec5, 0xdcf, 0xcc6, 0xbca, 0xac3, 0x9c9, 0x8c0,
56  0x8c0, 0x9c9, 0xac3, 0xbca, 0xcc6, 0xdcf, 0xec5, 0xfcc,
57  0xcc, 0x1c5, 0x2cf, 0x3c6, 0x4ca, 0x5c3, 0x6c9, 0x7c0,
58  0x950, 0x859, 0xb53, 0xa5a, 0xd56, 0xc5f, 0xf55, 0xe5c,
59  0x15c, 0x55, 0x35f, 0x256, 0x55a, 0x453, 0x759, 0x650,
60  0xaf0, 0xbf9, 0x8f3, 0x9fa, 0xef6, 0xfff, 0xcf5, 0xdfc,
61  0x2fc, 0x3f5, 0xff, 0x1f6, 0x6fa, 0x7f3, 0x4f9, 0x5f0,
62  0xb60, 0xa69, 0x963, 0x86a, 0xf66, 0xe6f, 0xd65, 0xc6c,
63  0x36c, 0x265, 0x16f, 0x66, 0x76a, 0x663, 0x569, 0x460,
64  0xca0, 0xda9, 0xea3, 0xfaa, 0x8a6, 0x9af, 0xaa5, 0xbac,
65  0x4ac, 0x5a5, 0x6af, 0x7a6, 0xaa, 0x1a3, 0x2a9, 0x3a0,
66  0xd30, 0xc39, 0xf33, 0xe3a, 0x936, 0x83f, 0xb35, 0xa3c,
67  0x53c, 0x435, 0x73f, 0x636, 0x13a, 0x33, 0x339, 0x230,
68  0xe90, 0xf99, 0xc93, 0xd9a, 0xa96, 0xb9f, 0x895, 0x99c,
69  0x69c, 0x795, 0x49f, 0x596, 0x29a, 0x393, 0x99, 0x190,
70  0xf00, 0xe09, 0xd03, 0xc0a, 0xb06, 0xa0f, 0x905, 0x80c,
71  0x70c, 0x605, 0x50f, 0x406, 0x30a, 0x203, 0x109, 0x0
72 };
73 
77 static int triTable[256][16] =
78 {
79  { -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
80  { 0, 8, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
81  { 0, 1, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
82  { 1, 8, 3, 9, 8, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
83  { 1, 2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
84  { 0, 8, 3, 1, 2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
85  { 9, 2, 10, 0, 2, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
86  { 2, 8, 3, 2, 10, 8, 10, 9, 8, -1, -1, -1, -1, -1, -1, -1 },
87  { 3, 11, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
88  { 0, 11, 2, 8, 11, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
89  { 1, 9, 0, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
90  { 1, 11, 2, 1, 9, 11, 9, 8, 11, -1, -1, -1, -1, -1, -1, -1 },
91  { 3, 10, 1, 11, 10, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
92  { 0, 10, 1, 0, 8, 10, 8, 11, 10, -1, -1, -1, -1, -1, -1, -1 },
93  { 3, 9, 0, 3, 11, 9, 11, 10, 9, -1, -1, -1, -1, -1, -1, -1 },
94  { 9, 8, 10, 10, 8, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
95  { 4, 7, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
96  { 4, 3, 0, 7, 3, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
97  { 0, 1, 9, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
98  { 4, 1, 9, 4, 7, 1, 7, 3, 1, -1, -1, -1, -1, -1, -1, -1 },
99  { 1, 2, 10, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
100  { 3, 4, 7, 3, 0, 4, 1, 2, 10, -1, -1, -1, -1, -1, -1, -1 },
101  { 9, 2, 10, 9, 0, 2, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1 },
102  { 2, 10, 9, 2, 9, 7, 2, 7, 3, 7, 9, 4, -1, -1, -1, -1 },
103  { 8, 4, 7, 3, 11, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
104  { 11, 4, 7, 11, 2, 4, 2, 0, 4, -1, -1, -1, -1, -1, -1, -1 },
105  { 9, 0, 1, 8, 4, 7, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1 },
106  { 4, 7, 11, 9, 4, 11, 9, 11, 2, 9, 2, 1, -1, -1, -1, -1 },
107  { 3, 10, 1, 3, 11, 10, 7, 8, 4, -1, -1, -1, -1, -1, -1, -1 },
108  { 1, 11, 10, 1, 4, 11, 1, 0, 4, 7, 11, 4, -1, -1, -1, -1 },
109  { 4, 7, 8, 9, 0, 11, 9, 11, 10, 11, 0, 3, -1, -1, -1, -1 },
110  { 4, 7, 11, 4, 11, 9, 9, 11, 10, -1, -1, -1, -1, -1, -1, -1 },
111  { 9, 5, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
112  { 9, 5, 4, 0, 8, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
113  { 0, 5, 4, 1, 5, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
114  { 8, 5, 4, 8, 3, 5, 3, 1, 5, -1, -1, -1, -1, -1, -1, -1 },
115  { 1, 2, 10, 9, 5, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
116  { 3, 0, 8, 1, 2, 10, 4, 9, 5, -1, -1, -1, -1, -1, -1, -1 },
117  { 5, 2, 10, 5, 4, 2, 4, 0, 2, -1, -1, -1, -1, -1, -1, -1 },
118  { 2, 10, 5, 3, 2, 5, 3, 5, 4, 3, 4, 8, -1, -1, -1, -1 },
119  { 9, 5, 4, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
120  { 0, 11, 2, 0, 8, 11, 4, 9, 5, -1, -1, -1, -1, -1, -1, -1 },
121  { 0, 5, 4, 0, 1, 5, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1 },
122  { 2, 1, 5, 2, 5, 8, 2, 8, 11, 4, 8, 5, -1, -1, -1, -1 },
123  { 10, 3, 11, 10, 1, 3, 9, 5, 4, -1, -1, -1, -1, -1, -1, -1 },
124  { 4, 9, 5, 0, 8, 1, 8, 10, 1, 8, 11, 10, -1, -1, -1, -1 },
125  { 5, 4, 0, 5, 0, 11, 5, 11, 10, 11, 0, 3, -1, -1, -1, -1 },
126  { 5, 4, 8, 5, 8, 10, 10, 8, 11, -1, -1, -1, -1, -1, -1, -1 },
127  { 9, 7, 8, 5, 7, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
128  { 9, 3, 0, 9, 5, 3, 5, 7, 3, -1, -1, -1, -1, -1, -1, -1 },
129  { 0, 7, 8, 0, 1, 7, 1, 5, 7, -1, -1, -1, -1, -1, -1, -1 },
130  { 1, 5, 3, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
131  { 9, 7, 8, 9, 5, 7, 10, 1, 2, -1, -1, -1, -1, -1, -1, -1 },
132  { 10, 1, 2, 9, 5, 0, 5, 3, 0, 5, 7, 3, -1, -1, -1, -1 },
133  { 8, 0, 2, 8, 2, 5, 8, 5, 7, 10, 5, 2, -1, -1, -1, -1 },
134  { 2, 10, 5, 2, 5, 3, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1 },
135  { 7, 9, 5, 7, 8, 9, 3, 11, 2, -1, -1, -1, -1, -1, -1, -1 },
136  { 9, 5, 7, 9, 7, 2, 9, 2, 0, 2, 7, 11, -1, -1, -1, -1 },
137  { 2, 3, 11, 0, 1, 8, 1, 7, 8, 1, 5, 7, -1, -1, -1, -1 },
138  { 11, 2, 1, 11, 1, 7, 7, 1, 5, -1, -1, -1, -1, -1, -1, -1 },
139  { 9, 5, 8, 8, 5, 7, 10, 1, 3, 10, 3, 11, -1, -1, -1, -1 },
140  { 5, 7, 0, 5, 0, 9, 7, 11, 0, 1, 0, 10, 11, 10, 0, -1 },
141  { 11, 10, 0, 11, 0, 3, 10, 5, 0, 8, 0, 7, 5, 7, 0, -1 },
142  { 11, 10, 5, 7, 11, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
143  { 10, 6, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
144  { 0, 8, 3, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
145  { 9, 0, 1, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
146  { 1, 8, 3, 1, 9, 8, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1 },
147  { 1, 6, 5, 2, 6, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
148  { 1, 6, 5, 1, 2, 6, 3, 0, 8, -1, -1, -1, -1, -1, -1, -1 },
149  { 9, 6, 5, 9, 0, 6, 0, 2, 6, -1, -1, -1, -1, -1, -1, -1 },
150  { 5, 9, 8, 5, 8, 2, 5, 2, 6, 3, 2, 8, -1, -1, -1, -1 },
151  { 2, 3, 11, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
152  { 11, 0, 8, 11, 2, 0, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1 },
153  { 0, 1, 9, 2, 3, 11, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1 },
154  { 5, 10, 6, 1, 9, 2, 9, 11, 2, 9, 8, 11, -1, -1, -1, -1 },
155  { 6, 3, 11, 6, 5, 3, 5, 1, 3, -1, -1, -1, -1, -1, -1, -1 },
156  { 0, 8, 11, 0, 11, 5, 0, 5, 1, 5, 11, 6, -1, -1, -1, -1 },
157  { 3, 11, 6, 0, 3, 6, 0, 6, 5, 0, 5, 9, -1, -1, -1, -1 },
158  { 6, 5, 9, 6, 9, 11, 11, 9, 8, -1, -1, -1, -1, -1, -1, -1 },
159  { 5, 10, 6, 4, 7, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
160  { 4, 3, 0, 4, 7, 3, 6, 5, 10, -1, -1, -1, -1, -1, -1, -1 },
161  { 1, 9, 0, 5, 10, 6, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1 },
162  { 10, 6, 5, 1, 9, 7, 1, 7, 3, 7, 9, 4, -1, -1, -1, -1 },
163  { 6, 1, 2, 6, 5, 1, 4, 7, 8, -1, -1, -1, -1, -1, -1, -1 },
164  { 1, 2, 5, 5, 2, 6, 3, 0, 4, 3, 4, 7, -1, -1, -1, -1 },
165  { 8, 4, 7, 9, 0, 5, 0, 6, 5, 0, 2, 6, -1, -1, -1, -1 },
166  { 7, 3, 9, 7, 9, 4, 3, 2, 9, 5, 9, 6, 2, 6, 9, -1 },
167  { 3, 11, 2, 7, 8, 4, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1 },
168  { 5, 10, 6, 4, 7, 2, 4, 2, 0, 2, 7, 11, -1, -1, -1, -1 },
169  { 0, 1, 9, 4, 7, 8, 2, 3, 11, 5, 10, 6, -1, -1, -1, -1 },
170  { 9, 2, 1, 9, 11, 2, 9, 4, 11, 7, 11, 4, 5, 10, 6, -1 },
171  { 8, 4, 7, 3, 11, 5, 3, 5, 1, 5, 11, 6, -1, -1, -1, -1 },
172  { 5, 1, 11, 5, 11, 6, 1, 0, 11, 7, 11, 4, 0, 4, 11, -1 },
173  { 0, 5, 9, 0, 6, 5, 0, 3, 6, 11, 6, 3, 8, 4, 7, -1 },
174  { 6, 5, 9, 6, 9, 11, 4, 7, 9, 7, 11, 9, -1, -1, -1, -1 },
175  { 10, 4, 9, 6, 4, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
176  { 4, 10, 6, 4, 9, 10, 0, 8, 3, -1, -1, -1, -1, -1, -1, -1 },
177  { 10, 0, 1, 10, 6, 0, 6, 4, 0, -1, -1, -1, -1, -1, -1, -1 },
178  { 8, 3, 1, 8, 1, 6, 8, 6, 4, 6, 1, 10, -1, -1, -1, -1 },
179  { 1, 4, 9, 1, 2, 4, 2, 6, 4, -1, -1, -1, -1, -1, -1, -1 },
180  { 3, 0, 8, 1, 2, 9, 2, 4, 9, 2, 6, 4, -1, -1, -1, -1 },
181  { 0, 2, 4, 4, 2, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
182  { 8, 3, 2, 8, 2, 4, 4, 2, 6, -1, -1, -1, -1, -1, -1, -1 },
183  { 10, 4, 9, 10, 6, 4, 11, 2, 3, -1, -1, -1, -1, -1, -1, -1 },
184  { 0, 8, 2, 2, 8, 11, 4, 9, 10, 4, 10, 6, -1, -1, -1, -1 },
185  { 3, 11, 2, 0, 1, 6, 0, 6, 4, 6, 1, 10, -1, -1, -1, -1 },
186  { 6, 4, 1, 6, 1, 10, 4, 8, 1, 2, 1, 11, 8, 11, 1, -1 },
187  { 9, 6, 4, 9, 3, 6, 9, 1, 3, 11, 6, 3, -1, -1, -1, -1 },
188  { 8, 11, 1, 8, 1, 0, 11, 6, 1, 9, 1, 4, 6, 4, 1, -1 },
189  { 3, 11, 6, 3, 6, 0, 0, 6, 4, -1, -1, -1, -1, -1, -1, -1 },
190  { 6, 4, 8, 11, 6, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
191  { 7, 10, 6, 7, 8, 10, 8, 9, 10, -1, -1, -1, -1, -1, -1, -1 },
192  { 0, 7, 3, 0, 10, 7, 0, 9, 10, 6, 7, 10, -1, -1, -1, -1 },
193  { 10, 6, 7, 1, 10, 7, 1, 7, 8, 1, 8, 0, -1, -1, -1, -1 },
194  { 10, 6, 7, 10, 7, 1, 1, 7, 3, -1, -1, -1, -1, -1, -1, -1 },
195  { 1, 2, 6, 1, 6, 8, 1, 8, 9, 8, 6, 7, -1, -1, -1, -1 },
196  { 2, 6, 9, 2, 9, 1, 6, 7, 9, 0, 9, 3, 7, 3, 9, -1 },
197  { 7, 8, 0, 7, 0, 6, 6, 0, 2, -1, -1, -1, -1, -1, -1, -1 },
198  { 7, 3, 2, 6, 7, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
199  { 2, 3, 11, 10, 6, 8, 10, 8, 9, 8, 6, 7, -1, -1, -1, -1 },
200  { 2, 0, 7, 2, 7, 11, 0, 9, 7, 6, 7, 10, 9, 10, 7, -1 },
201  { 1, 8, 0, 1, 7, 8, 1, 10, 7, 6, 7, 10, 2, 3, 11, -1 },
202  { 11, 2, 1, 11, 1, 7, 10, 6, 1, 6, 7, 1, -1, -1, -1, -1 },
203  { 8, 9, 6, 8, 6, 7, 9, 1, 6, 11, 6, 3, 1, 3, 6, -1 },
204  { 0, 9, 1, 11, 6, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
205  { 7, 8, 0, 7, 0, 6, 3, 11, 0, 11, 6, 0, -1, -1, -1, -1 },
206  { 7, 11, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
207  { 7, 6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
208  { 3, 0, 8, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
209  { 0, 1, 9, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
210  { 8, 1, 9, 8, 3, 1, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1 },
211  { 10, 1, 2, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
212  { 1, 2, 10, 3, 0, 8, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1 },
213  { 2, 9, 0, 2, 10, 9, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1 },
214  { 6, 11, 7, 2, 10, 3, 10, 8, 3, 10, 9, 8, -1, -1, -1, -1 },
215  { 7, 2, 3, 6, 2, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
216  { 7, 0, 8, 7, 6, 0, 6, 2, 0, -1, -1, -1, -1, -1, -1, -1 },
217  { 2, 7, 6, 2, 3, 7, 0, 1, 9, -1, -1, -1, -1, -1, -1, -1 },
218  { 1, 6, 2, 1, 8, 6, 1, 9, 8, 8, 7, 6, -1, -1, -1, -1 },
219  { 10, 7, 6, 10, 1, 7, 1, 3, 7, -1, -1, -1, -1, -1, -1, -1 },
220  { 10, 7, 6, 1, 7, 10, 1, 8, 7, 1, 0, 8, -1, -1, -1, -1 },
221  { 0, 3, 7, 0, 7, 10, 0, 10, 9, 6, 10, 7, -1, -1, -1, -1 },
222  { 7, 6, 10, 7, 10, 8, 8, 10, 9, -1, -1, -1, -1, -1, -1, -1 },
223  { 6, 8, 4, 11, 8, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
224  { 3, 6, 11, 3, 0, 6, 0, 4, 6, -1, -1, -1, -1, -1, -1, -1 },
225  { 8, 6, 11, 8, 4, 6, 9, 0, 1, -1, -1, -1, -1, -1, -1, -1 },
226  { 9, 4, 6, 9, 6, 3, 9, 3, 1, 11, 3, 6, -1, -1, -1, -1 },
227  { 6, 8, 4, 6, 11, 8, 2, 10, 1, -1, -1, -1, -1, -1, -1, -1 },
228  { 1, 2, 10, 3, 0, 11, 0, 6, 11, 0, 4, 6, -1, -1, -1, -1 },
229  { 4, 11, 8, 4, 6, 11, 0, 2, 9, 2, 10, 9, -1, -1, -1, -1 },
230  { 10, 9, 3, 10, 3, 2, 9, 4, 3, 11, 3, 6, 4, 6, 3, -1 },
231  { 8, 2, 3, 8, 4, 2, 4, 6, 2, -1, -1, -1, -1, -1, -1, -1 },
232  { 0, 4, 2, 4, 6, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
233  { 1, 9, 0, 2, 3, 4, 2, 4, 6, 4, 3, 8, -1, -1, -1, -1 },
234  { 1, 9, 4, 1, 4, 2, 2, 4, 6, -1, -1, -1, -1, -1, -1, -1 },
235  { 8, 1, 3, 8, 6, 1, 8, 4, 6, 6, 10, 1, -1, -1, -1, -1 },
236  { 10, 1, 0, 10, 0, 6, 6, 0, 4, -1, -1, -1, -1, -1, -1, -1 },
237  { 4, 6, 3, 4, 3, 8, 6, 10, 3, 0, 3, 9, 10, 9, 3, -1 },
238  { 10, 9, 4, 6, 10, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
239  { 4, 9, 5, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
240  { 0, 8, 3, 4, 9, 5, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1 },
241  { 5, 0, 1, 5, 4, 0, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1 },
242  { 11, 7, 6, 8, 3, 4, 3, 5, 4, 3, 1, 5, -1, -1, -1, -1 },
243  { 9, 5, 4, 10, 1, 2, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1 },
244  { 6, 11, 7, 1, 2, 10, 0, 8, 3, 4, 9, 5, -1, -1, -1, -1 },
245  { 7, 6, 11, 5, 4, 10, 4, 2, 10, 4, 0, 2, -1, -1, -1, -1 },
246  { 3, 4, 8, 3, 5, 4, 3, 2, 5, 10, 5, 2, 11, 7, 6, -1 },
247  { 7, 2, 3, 7, 6, 2, 5, 4, 9, -1, -1, -1, -1, -1, -1, -1 },
248  { 9, 5, 4, 0, 8, 6, 0, 6, 2, 6, 8, 7, -1, -1, -1, -1 },
249  { 3, 6, 2, 3, 7, 6, 1, 5, 0, 5, 4, 0, -1, -1, -1, -1 },
250  { 6, 2, 8, 6, 8, 7, 2, 1, 8, 4, 8, 5, 1, 5, 8, -1 },
251  { 9, 5, 4, 10, 1, 6, 1, 7, 6, 1, 3, 7, -1, -1, -1, -1 },
252  { 1, 6, 10, 1, 7, 6, 1, 0, 7, 8, 7, 0, 9, 5, 4, -1 },
253  { 4, 0, 10, 4, 10, 5, 0, 3, 10, 6, 10, 7, 3, 7, 10, -1 },
254  { 7, 6, 10, 7, 10, 8, 5, 4, 10, 4, 8, 10, -1, -1, -1, -1 },
255  { 6, 9, 5, 6, 11, 9, 11, 8, 9, -1, -1, -1, -1, -1, -1, -1 },
256  { 3, 6, 11, 0, 6, 3, 0, 5, 6, 0, 9, 5, -1, -1, -1, -1 },
257  { 0, 11, 8, 0, 5, 11, 0, 1, 5, 5, 6, 11, -1, -1, -1, -1 },
258  { 6, 11, 3, 6, 3, 5, 5, 3, 1, -1, -1, -1, -1, -1, -1, -1 },
259  { 1, 2, 10, 9, 5, 11, 9, 11, 8, 11, 5, 6, -1, -1, -1, -1 },
260  { 0, 11, 3, 0, 6, 11, 0, 9, 6, 5, 6, 9, 1, 2, 10, -1 },
261  { 11, 8, 5, 11, 5, 6, 8, 0, 5, 10, 5, 2, 0, 2, 5, -1 },
262  { 6, 11, 3, 6, 3, 5, 2, 10, 3, 10, 5, 3, -1, -1, -1, -1 },
263  { 5, 8, 9, 5, 2, 8, 5, 6, 2, 3, 8, 2, -1, -1, -1, -1 },
264  { 9, 5, 6, 9, 6, 0, 0, 6, 2, -1, -1, -1, -1, -1, -1, -1 },
265  { 1, 5, 8, 1, 8, 0, 5, 6, 8, 3, 8, 2, 6, 2, 8, -1 },
266  { 1, 5, 6, 2, 1, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
267  { 1, 3, 6, 1, 6, 10, 3, 8, 6, 5, 6, 9, 8, 9, 6, -1 },
268  { 10, 1, 0, 10, 0, 6, 9, 5, 0, 5, 6, 0, -1, -1, -1, -1 },
269  { 0, 3, 8, 5, 6, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
270  { 10, 5, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
271  { 11, 5, 10, 7, 5, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
272  { 11, 5, 10, 11, 7, 5, 8, 3, 0, -1, -1, -1, -1, -1, -1, -1 },
273  { 5, 11, 7, 5, 10, 11, 1, 9, 0, -1, -1, -1, -1, -1, -1, -1 },
274  { 10, 7, 5, 10, 11, 7, 9, 8, 1, 8, 3, 1, -1, -1, -1, -1 },
275  { 11, 1, 2, 11, 7, 1, 7, 5, 1, -1, -1, -1, -1, -1, -1, -1 },
276  { 0, 8, 3, 1, 2, 7, 1, 7, 5, 7, 2, 11, -1, -1, -1, -1 },
277  { 9, 7, 5, 9, 2, 7, 9, 0, 2, 2, 11, 7, -1, -1, -1, -1 },
278  { 7, 5, 2, 7, 2, 11, 5, 9, 2, 3, 2, 8, 9, 8, 2, -1 },
279  { 2, 5, 10, 2, 3, 5, 3, 7, 5, -1, -1, -1, -1, -1, -1, -1 },
280  { 8, 2, 0, 8, 5, 2, 8, 7, 5, 10, 2, 5, -1, -1, -1, -1 },
281  { 9, 0, 1, 5, 10, 3, 5, 3, 7, 3, 10, 2, -1, -1, -1, -1 },
282  { 9, 8, 2, 9, 2, 1, 8, 7, 2, 10, 2, 5, 7, 5, 2, -1 },
283  { 1, 3, 5, 3, 7, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
284  { 0, 8, 7, 0, 7, 1, 1, 7, 5, -1, -1, -1, -1, -1, -1, -1 },
285  { 9, 0, 3, 9, 3, 5, 5, 3, 7, -1, -1, -1, -1, -1, -1, -1 },
286  { 9, 8, 7, 5, 9, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
287  { 5, 8, 4, 5, 10, 8, 10, 11, 8, -1, -1, -1, -1, -1, -1, -1 },
288  { 5, 0, 4, 5, 11, 0, 5, 10, 11, 11, 3, 0, -1, -1, -1, -1 },
289  { 0, 1, 9, 8, 4, 10, 8, 10, 11, 10, 4, 5, -1, -1, -1, -1 },
290  { 10, 11, 4, 10, 4, 5, 11, 3, 4, 9, 4, 1, 3, 1, 4, -1 },
291  { 2, 5, 1, 2, 8, 5, 2, 11, 8, 4, 5, 8, -1, -1, -1, -1 },
292  { 0, 4, 11, 0, 11, 3, 4, 5, 11, 2, 11, 1, 5, 1, 11, -1 },
293  { 0, 2, 5, 0, 5, 9, 2, 11, 5, 4, 5, 8, 11, 8, 5, -1 },
294  { 9, 4, 5, 2, 11, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
295  { 2, 5, 10, 3, 5, 2, 3, 4, 5, 3, 8, 4, -1, -1, -1, -1 },
296  { 5, 10, 2, 5, 2, 4, 4, 2, 0, -1, -1, -1, -1, -1, -1, -1 },
297  { 3, 10, 2, 3, 5, 10, 3, 8, 5, 4, 5, 8, 0, 1, 9, -1 },
298  { 5, 10, 2, 5, 2, 4, 1, 9, 2, 9, 4, 2, -1, -1, -1, -1 },
299  { 8, 4, 5, 8, 5, 3, 3, 5, 1, -1, -1, -1, -1, -1, -1, -1 },
300  { 0, 4, 5, 1, 0, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
301  { 8, 4, 5, 8, 5, 3, 9, 0, 5, 0, 3, 5, -1, -1, -1, -1 },
302  { 9, 4, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
303  { 4, 11, 7, 4, 9, 11, 9, 10, 11, -1, -1, -1, -1, -1, -1, -1 },
304  { 0, 8, 3, 4, 9, 7, 9, 11, 7, 9, 10, 11, -1, -1, -1, -1 },
305  { 1, 10, 11, 1, 11, 4, 1, 4, 0, 7, 4, 11, -1, -1, -1, -1 },
306  { 3, 1, 4, 3, 4, 8, 1, 10, 4, 7, 4, 11, 10, 11, 4, -1 },
307  { 4, 11, 7, 9, 11, 4, 9, 2, 11, 9, 1, 2, -1, -1, -1, -1 },
308  { 9, 7, 4, 9, 11, 7, 9, 1, 11, 2, 11, 1, 0, 8, 3, -1 },
309  { 11, 7, 4, 11, 4, 2, 2, 4, 0, -1, -1, -1, -1, -1, -1, -1 },
310  { 11, 7, 4, 11, 4, 2, 8, 3, 4, 3, 2, 4, -1, -1, -1, -1 },
311  { 2, 9, 10, 2, 7, 9, 2, 3, 7, 7, 4, 9, -1, -1, -1, -1 },
312  { 9, 10, 7, 9, 7, 4, 10, 2, 7, 8, 7, 0, 2, 0, 7, -1 },
313  { 3, 7, 10, 3, 10, 2, 7, 4, 10, 1, 10, 0, 4, 0, 10, -1 },
314  { 1, 10, 2, 8, 7, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
315  { 4, 9, 1, 4, 1, 7, 7, 1, 3, -1, -1, -1, -1, -1, -1, -1 },
316  { 4, 9, 1, 4, 1, 7, 0, 8, 1, 8, 7, 1, -1, -1, -1, -1 },
317  { 4, 0, 3, 7, 4, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
318  { 4, 8, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
319  { 9, 10, 8, 10, 11, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
320  { 3, 0, 9, 3, 9, 11, 11, 9, 10, -1, -1, -1, -1, -1, -1, -1 },
321  { 0, 1, 10, 0, 10, 8, 8, 10, 11, -1, -1, -1, -1, -1, -1, -1 },
322  { 3, 1, 10, 11, 3, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
323  { 1, 2, 11, 1, 11, 9, 9, 11, 8, -1, -1, -1, -1, -1, -1, -1 },
324  { 3, 0, 9, 3, 9, 11, 1, 2, 9, 2, 11, 9, -1, -1, -1, -1 },
325  { 0, 2, 11, 8, 0, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
326  { 3, 2, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
327  { 2, 3, 8, 2, 8, 10, 10, 8, 9, -1, -1, -1, -1, -1, -1, -1 },
328  { 9, 10, 2, 0, 9, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
329  { 2, 3, 8, 2, 8, 10, 0, 1, 8, 1, 10, 8, -1, -1, -1, -1 },
330  { 1, 10, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
331  { 1, 3, 8, 9, 1, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
332  { 0, 9, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
333  { 0, 3, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
334  { -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 }
335 };
336 
340 static double
341 lerp(double val1, double val2, double isovalue, double spacing)
342 {
343  return (isovalue - val1) * spacing / (val2 - val1);
344 }
345 
346 LocalMarchingCubes::LocalMarchingCubes()
347 {
348  setNumInputPorts(1);
349  setRequiredInputType<ImageData>(0);
350 
351  setNumOutputPorts(0);
352 }
353 
354 std::shared_ptr<SurfaceMesh>
355 LocalMarchingCubes::getOutputMesh(const int i) const
356 {
357  return std::dynamic_pointer_cast<SurfaceMesh>(getOutput(i));
358 }
359 
360 void
361 LocalMarchingCubes::setInputImage(std::shared_ptr<ImageData> inputImage)
362 {
363  setInput(inputImage, 0);
364 }
365 
366 void
367 LocalMarchingCubes::setModified(const Vec3i& coord)
368 {
369  std::shared_ptr<ImageData> image = std::dynamic_pointer_cast<ImageData>(getInput(0));
370  const int voxelIndex = static_cast<int>(image->getScalarIndex(coord));
371  m_modifiedVoxels[voxelIndex] = coord;
372  //m_modifiedVoxels.push_back(std::pair<int, Vec3i>(static_cast<int>(voxelIndex), coord));
373 }
374 
375 void
376 LocalMarchingCubes::setNumberOfChunks(const Vec3i& numChunks)
377 {
378  m_chunkCount = numChunks[0] * numChunks[1] * numChunks[2];
379  m_numChunks = numChunks;
380  setNumOutputPorts(m_chunkCount);
381  for (size_t i = 0; i < m_chunkCount; i++)
382  {
383  setOutput(std::make_shared<SurfaceMesh>(), i);
384  }
385  m_allModified = true;
386 }
387 
388 void
389 mcSubImage(std::shared_ptr<ImageData> imageData, std::shared_ptr<SurfaceMesh> outputSurf, const Vec3i& start, const Vec3i& end, const double isoValue)
390 {
391  const Vec3i& fullDims = imageData->getDimensions();
392  const Vec3d& spacing = imageData->getSpacing();
393  const Vec3d& origin = imageData->getOrigin();
394  const Vec3d shift = origin + spacing * 0.5;
395 
396  std::list<Vec3d> vertices;
397  std::list<Vec3i> indices;
398 
399  // For every voxel, identify case
400  double* imgPtr = static_cast<double*>(imageData->getScalars()->getVoidPointer());
401  //int blockIndex = 0;
402  // Iterate along the dual, assigning case numbers per paper
403  for (int z1 = start[2]; z1 < end[2]; z1++)
404  {
405  const int z2 = z1 + 1;
406  for (int y1 = start[1]; y1 < end[1]; y1++)
407  {
408  const int y2 = y1 + 1;
409  for (int x1 = start[0]; x1 < end[0]; x1++ /*, blockIndex++*/)
410  {
411  const int x2 = x1 + 1;
412 
413  const int i000 = static_cast<int>(imageData->getScalarIndex(Vec3i(x1, y1, z1)));
414  const int i001 = i000 + 1;
415 
416  const int i010 = i000 + fullDims[0];
417  const int i011 = i010 + 1;
418 
419  const int i100 = i000 + fullDims[0] * fullDims[1];
420  const int i101 = i100 + 1;
421  const int i110 = i100 + fullDims[0];
422  const int i111 = i110 + 1;
423 
424  const double vals[8] =
425  {
426  imgPtr[i000], //v0
427  imgPtr[i001], //v1
428  imgPtr[i101], //v2
429  imgPtr[i100], //v3
430  imgPtr[i010], //v4
431  imgPtr[i011], //v5
432  imgPtr[i111], //v6
433  imgPtr[i110] //v7
434  };
435 
436  // Assign a case
437  int mcCase = 0;
438  {
439  for (int i = 0; i < 8; i++)
440  {
441  if (vals[i] < isoValue)
442  {
443  mcCase |= (1 << i);
444  }
445  }
446  }
447 
448  // Generate the vertices
449  int cubeIndices[12] = { -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 }; // Map of local cube indices -> global vertex indices
450  {
451  const Vec3d min = Vec3d(x1, y1, z1).cwiseProduct(spacing);
452  const Vec3d max = Vec3d(x2, y2, z2).cwiseProduct(spacing);
453 
454  const int packedEdges = edgeTable[mcCase];
455 
456  if (packedEdges == 0)
457  {
458  continue;
459  }
460 
461  // If vertex from e1 should be added
462  if (packedEdges & 1)
463  {
464  const double x = lerp(vals[0], vals[1], isoValue, spacing[0]);
465  cubeIndices[0] = static_cast<int>(vertices.size());
466  vertices.push_back(Vec3d(min[0] + x, min[1], min[2]) + shift);
467  }
468  // If vertex from e2 should be added
469  if (packedEdges & 2)
470  {
471  const double z = lerp(vals[1], vals[2], isoValue, spacing[2]);
472  cubeIndices[1] = static_cast<int>(vertices.size());
473  vertices.push_back(Vec3d(max[0], min[1], min[2] + z) + shift);
474  }
475  // If vertex from e3 should be added
476  if (packedEdges & 4)
477  {
478  const double x = lerp(vals[3], vals[2], isoValue, spacing[0]);
479  cubeIndices[2] = static_cast<int>(vertices.size());
480  vertices.push_back(Vec3d(min[0] + x, min[1], max[2]) + shift);
481  }
482  // If vertex from e4 should be added
483  if (packedEdges & 8)
484  {
485  const double z = lerp(vals[0], vals[3], isoValue, spacing[2]);
486  cubeIndices[3] = static_cast<int>(vertices.size());
487  vertices.push_back(Vec3d(min[0], min[1], min[2] + z) + shift);
488  }
489  // If vertex from e5 should be added
490  if (packedEdges & 16)
491  {
492  const double x = lerp(vals[4], vals[5], isoValue, spacing[0]);
493  cubeIndices[4] = static_cast<int>(vertices.size());
494  vertices.push_back(Vec3d(min[0] + x, max[1], min[2]) + shift);
495  }
496  // If vertex from e6 should be added
497  if (packedEdges & 32)
498  {
499  const double z = lerp(vals[5], vals[6], isoValue, spacing[2]);
500  cubeIndices[5] = static_cast<int>(vertices.size());
501  vertices.push_back(Vec3d(max[0], max[1], min[2] + z) + shift);
502  }
503  // If vertex from e7 should be added
504  if (packedEdges & 64)
505  {
506  const double x = lerp(vals[7], vals[6], isoValue, spacing[0]);
507  cubeIndices[6] = static_cast<int>(vertices.size());
508  vertices.push_back(Vec3d(min[0] + x, max[1], max[2]) + shift);
509  }
510  // If vertex from e8 should be added
511  if (packedEdges & 128)
512  {
513  const double z = lerp(vals[4], vals[7], isoValue, spacing[2]);
514  cubeIndices[7] = static_cast<int>(vertices.size());
515  vertices.push_back(Vec3d(min[0], max[1], min[2] + z) + shift);
516  }
517  // If vertex from e9 should be added
518  if (packedEdges & 256)
519  {
520  const double y = lerp(vals[0], vals[4], isoValue, spacing[1]);
521  cubeIndices[8] = static_cast<int>(vertices.size());
522  vertices.push_back(Vec3d(min[0], min[1] + y, min[2]) + shift);
523  }
524  // If vertex from e10 should be added
525  if (packedEdges & 512)
526  {
527  const double y = lerp(vals[1], vals[5], isoValue, spacing[1]);
528  cubeIndices[9] = static_cast<int>(vertices.size());
529  vertices.push_back(Vec3d(max[0], min[1] + y, min[2]) + shift);
530  }
531  // If vertex from e11 should be added
532  if (packedEdges & 1024)
533  {
534  const double y = lerp(vals[2], vals[6], isoValue, spacing[1]);
535  cubeIndices[10] = static_cast<int>(vertices.size());
536  vertices.push_back(Vec3d(max[0], min[1] + y, max[2]) + shift);
537  }
538  // If vertex from e12 should be added
539  if (packedEdges & 2048)
540  {
541  const double y = lerp(vals[3], vals[7], isoValue, spacing[1]);
542  cubeIndices[11] = static_cast<int>(vertices.size());
543  vertices.push_back(Vec3d(min[0], min[1] + y, max[2]) + shift);
544  }
545  }
546 
547  // Generate the triangles
548  {
549  for (int i = 0; triTable[mcCase][i] != -1; i += 3)
550  {
551  // Local cube index
552  const int i1Local = triTable[mcCase][i];
553  const int i2Local = triTable[mcCase][i + 1];
554  const int i3Local = triTable[mcCase][i + 2];
555 
556  indices.push_back(Vec3i(static_cast<size_t>(cubeIndices[i1Local]), static_cast<size_t>(cubeIndices[i2Local]), static_cast<size_t>(cubeIndices[i3Local])));
557  }
558  }
559  }
560  }
561  }
562 
563  // Temp
564  std::shared_ptr<VecDataArray<double, 3>> vecVerticesPtr = std::make_shared<VecDataArray<double, 3>>(static_cast<int>(vertices.size()));
565  VecDataArray<double, 3>& vecVertices = *vecVerticesPtr;
566  //std::copy(vertices.begin(), vertices.end(), vecVertices->begin());
567  std::shared_ptr<VecDataArray<int, 3>> vecIndicesPtr = std::make_shared<VecDataArray<int, 3>>(static_cast<int>(indices.size()));
568  VecDataArray<int, 3>& vecIndices = *vecIndicesPtr;
569  //std::copy(indices.begin(), indices.end(), vecIndices->begin());
570 
571  int j = 0;
572  for (auto i : indices)
573  {
574  vecIndices[j++] = i;
575  }
576  j = 0;
577  for (auto i : vertices)
578  {
579  vecVertices[j++] = i;
580  }
581 
582  outputSurf->initialize(vecVerticesPtr, vecIndicesPtr);
583  outputSurf->computeVertexNormals();
584 }
585 
586 void
587 LocalMarchingCubes::requestUpdate()
588 {
589  // Foreword:
590  // Voxels are referred to as the elements that make up the image, each with an intensity
591  // Blocks are referred to as the elements inbetween the voxels for MC. Sort of like the "dual" of the graph.
592  std::shared_ptr<ImageData> imageData = std::dynamic_pointer_cast<ImageData>(getInput(0));
593 
594  if (imageData->getScalarType() != IMSTK_DOUBLE)
595  {
596  LOG(WARNING) << "Local marching cubes only supports doubles";
597  return;
598  }
599 
600  //const Vec3d& spacing = imageData->getSpacing();
601  //const Vec3d& origin = imageData->getOrigin();
602  //const Vec3d shift = origin + spacing * 0.5;
603  const Vec3i& dims = imageData->getDimensions();
604 
605  // The dimensions must be divisible by number of subdivisions
606  // Increasingly adjust until divisible
607  const Vec3i m_previousNumChunks = m_numChunks;
608  while ((dims[0] - 1) % m_numChunks[0] != 0)
609  {
610  m_numChunks[0]++;
611  }
612  while ((dims[1] - 1) % m_numChunks[1] != 0)
613  {
614  m_numChunks[1]++;
615  }
616  while ((dims[2] - 1) % m_numChunks[2] != 0)
617  {
618  m_numChunks[2]++;
619  }
620  if (m_previousNumChunks != m_numChunks)
621  {
622  LOG(WARNING) << "LocalMarchingCubes: Provided with chunk dimensions that aren't divisible with the image dimensions, increasing to next divisor";
623  LOG(WARNING) << "Requested # chunks: " << m_previousNumChunks[0] << ", " << m_previousNumChunks[1] << ", " << m_previousNumChunks[2];
624  LOG(WARNING) << "Not divisible with: " << dims[0] - 1 << ", " << dims[1] - 1 << ", " << dims[2] - 1;
625  if (m_numChunks == dims)
626  {
627  LOG(WARNING) << "LocalMarchingCubes: No divisor found";
628  return;
629  }
630 
631  setNumberOfChunks(m_numChunks);
632  }
633 
634  const Vec3i chunkDimensions = (dims - Vec3i(1, 1, 1)).cwiseQuotient(m_numChunks);
635 
636  if (m_allModified)
637  {
638  // For every chunk
639  size_t i = 0;
640  for (int z = 0; z < m_numChunks[2]; z++)
641  {
642  for (int y = 0; y < m_numChunks[1]; y++)
643  {
644  for (int x = 0; x < m_numChunks[0]; x++, i++)
645  {
646  const Vec3i chunkCoord = Vec3i(x, y, z);
647  std::shared_ptr<SurfaceMesh> outputSurf = std::dynamic_pointer_cast<SurfaceMesh>(getOutput(i));
648  const Vec3i coordStart = chunkCoord.cwiseProduct(chunkDimensions);
649 
650  mcSubImage(imageData, outputSurf, coordStart, coordStart + chunkDimensions, m_isoValue);
651  }
652  }
653  }
654  m_allModified = false;
655  }
656  else
657  {
658  const Vec3i dim1 = dims - Vec3i(1, 1, 1);
659 
660  // Set of modified blocks (all in one map to remove duplicates)
661  std::unordered_map<int, Vec3i> modifiedBlocks;
662  {
663  // We need to compute the set of blocks that are modified given the modified voxels
664  // then we can compute the set of chunks modified
665  for (auto modifiedVoxel : m_modifiedVoxels)
666  {
667  const Vec3i& voxelCoord = modifiedVoxel.second;
668 
669  // 2x2 iteration, index space for blocks is -(1, 1, 1)
670  const Vec3i minCoord = (voxelCoord - Vec3i(1, 1, 1)).cwiseMax(Vec3i(0, 0, 0));
671  const Vec3i maxCoord = (voxelCoord + Vec3i(1, 1, 1)).cwiseMin(dim1 - Vec3i(1, 1, 1));
672  for (int z = minCoord[2]; z < maxCoord[2]; z++)
673  {
674  for (int y = minCoord[1]; y < maxCoord[1]; y++)
675  {
676  for (int x = minCoord[0]; x < maxCoord[0]; x++)
677  {
678  modifiedBlocks[x + dim1[0] * (y + z * dim1[1])] = Vec3i(x, y, z);
679  }
680  }
681  }
682  }
683  }
684 
685  // Set of modified chunks (all in one map to remove duplicates)
686  m_modifiedChunks.clear();
687  {
688  for (auto modifiedBlock : modifiedBlocks)
689  {
690  //const int blockIndex = modifiedBlock.first;
691  const Vec3i& blockCoord = modifiedBlock.second;
692 
693  // Compute the chunk
694  const Vec3i chunkCoord = blockCoord.cast<double>().cwiseQuotient(dim1.cast<double>()).cwiseProduct(m_numChunks.cast<double>()).cast<int>();
695 
696  const int chunkId = chunkCoord[0] + (chunkCoord[1] + chunkCoord[2] * m_numChunks[1]) * m_numChunks[0];
697  m_modifiedChunks.insert(std::pair<int, Vec3i>(chunkId, chunkCoord));
698  }
699  }
700  // Updates all the modified chunks
701  for (auto i : m_modifiedChunks)
702  {
703  const Vec3i& chunkCoord = i.second;
704  std::shared_ptr<SurfaceMesh> outputSurf = std::dynamic_pointer_cast<SurfaceMesh>(getOutput(i.first));
705  const Vec3i coordStart = chunkCoord.cwiseProduct(chunkDimensions);
706  mcSubImage(imageData, outputSurf, coordStart, coordStart + chunkDimensions, m_isoValue);
707  outputSurf->postModified();
708  }
709  m_modifiedVoxels.clear();
710  }
711 }
712 } // namespace imstk
imstk::LocalMarchingCubes::setModified
void setModified(const Vec3i &coord)
Set a voxel that was modified in the image (the neighboring dual voxels will be updated on the next r...
Definition: imstkLocalMarchingCubes.cpp:367
imstk::GeometryAlgorithm::getInput
std::shared_ptr< Geometry > getInput(size_t port=0) const
Returns input geometry given port, returns nullptr if doesn&#39;t exist.
Definition: imstkGeometryAlgorithm.h:50
imstk
Compound Geometry.
Definition: OctreeDebugModel.cpp:10
imstk::LocalMarchingCubes::setNumberOfChunks
void setNumberOfChunks(const Vec3i &numChunks)
Set the number of chunks. one minus the dimensions of image must be divisible by numChunks ((dimensio...
Definition: imstkLocalMarchingCubes.cpp:376
imstk::GeometryAlgorithm::setNumOutputPorts
void setNumOutputPorts(const size_t numPorts)
Get/Set the amount of output ports.
Definition: imstkGeometryAlgorithm.cpp:49
imstk::GeometryAlgorithm::getOutput
std::shared_ptr< Geometry > getOutput(size_t port=0) const
Returns output geometry given port, returns nullptr if doesn&#39;t exist.
Definition: imstkGeometryAlgorithm.h:58
imstk::GeometryAlgorithm::setInput
void setInput(std::shared_ptr< Geometry > inputGeometry, size_t port=0)
Set the input at the port.
Definition: imstkGeometryAlgorithm.cpp:14
imstk::LocalMarchingCubes::requestUpdate
void requestUpdate() override
Users can implement this for the logic to be run.
Definition: imstkLocalMarchingCubes.cpp:587
imstk::SurfaceMesh
Represents a set of triangles & vertices via an array of Vec3d double vertices & Vec3i integer indice...
Definition: imstkSurfaceMesh.h:124
imstk::GeometryAlgorithm::setOutput
void setOutput(std::shared_ptr< Geometry > inputGeometry, const size_t port=0)
Set the output at the port.
Definition: imstkGeometryAlgorithm.cpp:25
imstk::GeometryAlgorithm::setNumInputPorts
void setNumInputPorts(const size_t numPorts)
Get/Set the amount of input ports.
Definition: imstkGeometryAlgorithm.cpp:35
imstk::ImageData
Class to represent 1, 2, or 3D image data (i.e. structured points)
Definition: imstkImageData.h:22
Generated by   doxygen 1.8.13