File: BrainModelVolumeTopologyGraph.h

package info (click to toggle)
caret 5.6.4~dfsg.1-3
  • links: PTS, VCS
  • area: main
  • in suites: jessie, jessie-kfreebsd, wheezy
  • size: 31,904 kB
  • ctags: 28,901
  • sloc: cpp: 378,050; python: 6,718; ansic: 5,507; makefile: 333; sh: 46
file content (417 lines) | stat: -rw-r--r-- 15,090 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417

#ifndef __BRAIN_MODEL_VOLUME_TOPOLOGY_GRAPH_H__
#define __BRAIN_MODEL_VOLUME_TOPOLOGY_GRAPH_H__

/*LICENSE_START*/
/*
 *  Copyright 1995-2002 Washington University School of Medicine
 *
 *  http://brainmap.wustl.edu
 *
 *  This file is part of CARET.
 *
 *  CARET is free software; you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation; either version 2 of the License, or
 *  (at your option) any later version.
 *
 *  CARET is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU General Public License for more details.
 *
 *  You should have received a copy of the GNU General Public License
 *  along with CARET; if not, write to the Free Software
 *  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 *
 */
/*LICENSE_END*/

#include <ios>
#include <map>
#include <set>
#include <vector>

#include "BrainModelAlgorithm.h"
#include "VoxelIJK.h"

class VolumeFile;

/// create graph of segmentation volume topology
/// This algorithm is based off the paper "Automated Topology Correction
/// for Human Brain Segmentation" by Lin Chen and Gudrun Wagenknecht, 
/// MICCAI 2006, LNCS 4191, pp. 316-323, 2006.  Springer-Verlag.
class BrainModelVolumeTopologyGraph : public BrainModelAlgorithm {
   public:
      // search axis
      enum SEARCH_AXIS {
         // search along X axis
         SEARCH_AXIS_X,
         // search along Y axis
         SEARCH_AXIS_Y,
         // search along Z axis
         SEARCH_AXIS_Z
      };
      
      // voxel neighbor connectivity
      enum VOXEL_NEIGHBOR_CONNECTIVITY {
         // 6-connected
         VOXEL_NEIGHBOR_CONNECTIVITY_6,
         // 18-connected
         VOXEL_NEIGHBOR_CONNECTIVITY_18,
         // 26-connected
         VOXEL_NEIGHBOR_CONNECTIVITY_26
      };
      
      /// Edge in the graph
      class GraphEdge {
         public:
            // constructor
            GraphEdge(const int vertexNumberIn,
                      const int strengthIn) 
               { vertexNumber = vertexNumberIn; strength = strengthIn; }

            // destructor
            ~GraphEdge() { }
            
            // get the vertex number
            int getVertexNumber() const { return vertexNumber; }
            
            // get the strength
            int getStrength() const { return strength; }
            
         protected:
            // the vertex number
            int vertexNumber;
            
            // the strength
            int strength;
      };
      
      /// Vertex in the graph
      class GraphVertex {
         public:
            // constructor
            GraphVertex(const int sliceNumberIn) 
               { sliceNumber = sliceNumberIn; identifier = -1; }
            
            // destructor
            ~GraphVertex() { voxels.clear(); }
            
            /// add a voxel to the vertex
            void addVoxel(const VoxelIJK& v) { voxels.push_back(v); }
            
            /// add connected graph vertex index
            void addConnectedGraphVertex(const int vertexIndex, const int strength) 
               { connectedGraphEdges.push_back(GraphEdge(vertexIndex, strength)); }
               
            /// get number of connected graph edges
            int getNumberOfConnectedGraphEdges() const 
               { return connectedGraphEdges.size(); }
         
            /// get connected graph vertex index
            const GraphEdge* getConnectedGraphEdge(const int indx) const 
               { return &connectedGraphEdges[indx]; }
               
            /// get the number of voxels in the vertex
            int getNumberOfVoxels() const { return voxels.size(); }
            
            /// get a voxel
            const VoxelIJK* getVoxel(const int indx) const { return &voxels[indx]; }

            /// get the slice number
            int getSliceNumber() const { return sliceNumber; }
            
            /// get the identifier
            int getIdentifier() const { return identifier; }
            
            /// set the identifier
            void setIdentifier(const int id) { identifier = id; }
            
            /// get a descriptive name (slice_number and identifier)
            QString getDescriptiveName() const
               { return ("S" + QString::number(sliceNumber) 
                         //+ "I" + QString::number(identifier)
                         + "N" + QString::number(voxels.size())); }
                            
            /// the slice number
            int sliceNumber;
            
            /// voxels in vertex
            std::vector<VoxelIJK> voxels;
            
            /// connected graph edges
            std::vector<GraphEdge> connectedGraphEdges;
            
            /// identifier
            int identifier;
      };
      
      /// cycle in the graph
      class GraphCycle {
         public:
            // constructor
            GraphCycle();
            
            // destructor
            ~GraphCycle();

            // clear the cycle
            void clear();
            
            /// is cycle empty
            bool empty() const { return cycle.empty(); }
            
            // number of graph vertices in cycle
            int getNumberOfGraphVerticesInCycle() const { return cycle.size(); }
            
            // get index of graph vertex in cycle
            int getGraphVertexIndex(const int indx) const { return cycle[indx]; }
            
            // set the cycle
            void set(const std::vector<int>& cycleVerticesIn,
                     const std::vector<int>& cycleSlicesIn);
            
            /// get the handle vertices
            std::vector<int> getHandleVertices() const { return handleVertices; }
            
            /// get the number of voxels that make up the handle
            int getHandleSizeInVoxels() const { return numVoxelsInHandle; }
            
            // set the vertices that form the handle
            void setHandleVertices(const std::vector<int>& handleVerticesIn,
                                   const int numVoxelsInHandleIn);
            
            // equality operator
            bool operator==(const GraphCycle& c) const;
            
            // comparison operator
            bool operator<(const GraphCycle& c) const;
            
            // get the cycle
            std::vector<int> getCycle() const;

         protected:
            /// the cycle
            std::vector<int> cycle;
            
            /// the cycle in sorted order (used for comparisons)
            std::vector<int> cycleSorted;
            
            /// the vertices that form the handle in the cycle
            std::vector<int> handleVertices;
            
            /// number of voxels in the handle
            int numVoxelsInHandle;
      };

      // constructor
      BrainModelVolumeTopologyGraph(BrainSet* bsIn,
                                    const VolumeFile* segmentationVolumeFileIn,
                                    const SEARCH_AXIS searchAxisIn,
                                    const VOXEL_NEIGHBOR_CONNECTIVITY voxelConnectivityIn);
                                    
      // destructor
      ~BrainModelVolumeTopologyGraph();
      
      // execute the algorithm
      void execute() throw (BrainModelAlgorithmException);
      
      /// get the number of vertices in the graph
      int getNumberOfGraphVertices() const { return graphVertices.size(); }
      
      /// get a vertex in the graph
      GraphVertex* getGraphVertex(const int indx) { return graphVertices[indx]; }
      
      /// get a vertex in the graph (const method)
      const GraphVertex* getGraphVertex(const int indx) const { return graphVertices[indx]; }
      
      /// get number of cycles in graph
      int getNumberOfGraphCycles() const { return graphCycles.size(); }
      
      /// get a cycle in the graph
      GraphCycle* getGraphCycle(const int indx) { return &graphCycles[indx]; }
      
      /// get a cycle in the graph (const method)
      const GraphCycle* getGraphCycle(const int indx) const { return &graphCycles[indx]; }
      
      /// get cycle with smallest vertex (vertex that contains fewest voxels)
      void getGraphCycleWithSmallestVertex(int &cycleIndexOut,
                                           int &vertexIndexOut,
                                           int &numberofVoxelsOut) const;
      
      /// get cycle with smallest handle (handle that contains fewest voxels)
      void getGraphCycleWithSmallestHandle(int &cycleIndexOut,
                                           std::vector<int>& vertexIndicesOut,
                                           int &numberOfVoxelsOut) const;
      
      // print the results
      void printResults() const;
      
      // get the search axis
      SEARCH_AXIS getSearchAxis() const { return searchAxis; }
      
      // create an paint volume file containing the handles
      void createHandlesPaintVolumeFile(VolumeFile& handlesPaintVolumeFile);
      
      // write the graph to a paint volume file
      void writeGraphToPaintVolumeFile(const QString& paintVolumeFileName) const throw (BrainModelAlgorithmException);
      
      // write graph to graphviz file
      void writeGraphVizDotFile(const QString& dotFileName) const throw (BrainModelAlgorithmException);
      
   protected:    
      /// slice neighbor connectivity
      enum SLICE_NEIGHBOR_CONNECTIVITY {
         /// 4-connected
         SLICE_NEIGHBOR_CONNECTIVITY_4,
         /// 8-connected
         SLICE_NEIGHBOR_CONNECTIVITY_8
      };
        
      /// vertex visited status
      enum VISITED_STATUS { 
         NOT_VISITED = 0,
         VISITED = 1
      };
      
      /// A slice of a volume
      class VolumeSlice {
         public:
            // constructor
            VolumeSlice(const int dimIin, const int dimJin);
         
            // destructor
            ~VolumeSlice();
            
            // get indices valid
            bool getIJValid(const int i, const int j) const;
            
            // get dimension I
            int getDimI() const { return dimI; }
            
            // get dimension J
            int getDimJ() const { return dimJ; }
            
            // set a voxel in the slice
            void setVoxel(const SEARCH_AXIS searchAxis, const VoxelIJK& v, const int value);
            
            // get a voxel in the slice
            float getVoxel(const SEARCH_AXIS searchAxis, const VoxelIJK& v) const;
            
            // set a voxel in the slice
            void setVoxel(const int horizIndex, const int vertIndex, const int value);
            
            // get a voxel in the slice
            float getVoxel(const int horizIndex, const int vertIndex) const;
            
            // set all voxels in the slice
            void setAllVoxels(const int value);
         
         protected:
            // get the index
            int getIndex(const SEARCH_AXIS searchAxis, const VoxelIJK& v) const;
            
            // get the index
            int getIndex(const int horizIndex, const int vertIndex) const;
            
            /// the voxels
            int* voxels;
            
            /// dimension I
            int dimI;
            
            /// dimension J
            int dimJ;
      };
      
      // create the vertices in the graph
      void createGraphVertices() throw (BrainModelAlgorithmException);
      
      // create the edges (connections between vertices) in the graph
      void createGraphEdges() throw (BrainModelAlgorithmException);
      
      // search graph for cycles
      void searchGraphForCycles();
      
      // determine the handles
      void determineHandles();
      
      // perform breadth first search
      void breadthFirstSearchForCycles(const int startVertexIndex,
                                       const int searchForVertexIndex,
                                       GraphCycle& cycleOut);
      
      // get voxel slice neighbors
      void getVoxelSliceNeighbors(const VoxelIJK& v,
                             const VolumeSlice& slice,
                             const int searchForValue,
                             std::vector<VoxelIJK>& neighborsOut) const;
      
      // add a slice neighboring voxel if valid and not searched
      void addSliceNeighbor(const VolumeSlice& slice,
                       const int i,
                       const int j,
                       const int k, 
                       const int searchForValue,
                       std::vector<VoxelIJK>& neighborsOut) const;
      
      // get graph vertex connected neighbors
      void getGraphVertexConnectedNeighbors(const VoxelIJK& v,
                                 const bool adjoiningSlicesOnlyFlag,
                                 std::map<int,int>& graphVertexNeighborsInOut) const;
      
      // get graph vertex connected neighbors in next slice only
      void getGraphVertexConnectedNeighborsInNextSlice(const VoxelIJK& v,
                                 std::map<int,int>& graphVertexNeighborsInOut) const;
      
      // add a graph vertex neighbor
      void addGraphVertexNeighbor(const int i,
                                  const int j,
                                  const int k,
                                  std::set<int>& neighborsOut) const;
                             
      // print vertices in the graph
      void printGraphVertices() const;
      
      // print cycles in the graph
      void printGraphCycles() const;
      
      /// adjust IJK for axis being processed
      void ijkForSlice(int& i, int& j, int& k) const;
      
      /// get IJK from loop
      void ijkFromLoop(const int horiz, const int vert, const int slice,
                       int& i, int& j, int& k) const;
      
      /// Get voxel connects to voxels in another vertex.
      bool getVoxelConnectedToGraphVertex(const VoxelIJK& v1,
                                          const int otherGraphVertexIndex) const;
      
      /// input segmentation volume
      const VolumeFile* inputSegmentationVolumeFile;
      
      /// segmentation volume on which topology graph is built
      VolumeFile* segmentationVolumeFile;
      
      /// vertices in the graph
      std::vector<GraphVertex*> graphVertices;
      
      /// cycles in the graph
      std::vector<GraphCycle> graphCycles;
      
      /// the search axis
      SEARCH_AXIS searchAxis;

      /// the voxel connectivity
      VOXEL_NEIGHBOR_CONNECTIVITY volumeConnectivity;
      
      /// the slice connectivity
      SLICE_NEIGHBOR_CONNECTIVITY sliceConnectivity;
      
      // volume where the voxel values are the index of the graph 
      // vertex to which the voxel belongs
      VolumeFile* voxelGraphVertexIndexVolumeFile;
};

#endif // __BRAIN_MODEL_VOLUME_TOPOLOGY_GRAPH_H__