This repository has been archived by the owner on Jan 19, 2020. It is now read-only.
/
vtk2itk.cxx
370 lines (335 loc) · 11.5 KB
/
vtk2itk.cxx
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
/*=========================================================================
Program: Insight Segmentation & Registration Toolkit
Module: vtk2itk.cxx
Language: C++
Date: $Date$
Version: $Revision$
Copyright (c) 2002 Insight Consortium. All rights reserved.
See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
This software is distributed WITHOUT ANY WARRANTY; without even
the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
PURPOSE. See the above copyright notices for more information.
=========================================================================*/
// Disable warning for long symbol names in this file only
#ifdef _MSC_VER
#pragma warning ( disable : 4786 )
#endif
#include <iostream>
#include "itkMesh.h"
#include "itkTriangleCell.h"
#include "itkQuadrilateralCell.h"
#include "vtkDataSetWriter.h"
#include "vtkDataSetMapper.h"
#include "vtkRenderer.h"
#include "vtkRenderWindow.h"
#include "vtkActor.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkDataSetReader.h"
#include "vtkUnstructuredGrid.h"
#include "vtkDataSet.h"
#include "vtkCellArray.h"
#ifndef vtkFloatingPointType
#define vtkFloatingPointType float
#endif
typedef itk::Mesh<vtkFloatingPointType, 3,
itk::DefaultStaticMeshTraits< vtkFloatingPointType, 3, 3, vtkFloatingPointType, vtkFloatingPointType > > floatMesh;
void Display(vtkUnstructuredGrid* itkgrid, vtkUnstructuredGrid* vtkgrid)
{
// Create the renderer and window stuff
vtkRenderer* ren1 = vtkRenderer::New();
vtkRenderer* ren2 = vtkRenderer::New();
vtkRenderWindow* renWin = vtkRenderWindow::New();
renWin->AddRenderer(ren1);
renWin->AddRenderer(ren2);
vtkRenderWindowInteractor* inter = vtkRenderWindowInteractor::New();
inter->SetRenderWindow(renWin);
vtkDataSetMapper* mapper = vtkDataSetMapper::New();
#if VTK_MAJOR_VERSION <= 5
mapper->SetInput(itkgrid);
#else
mapper->SetInputData(itkgrid);
#endif
vtkActor* actor = vtkActor::New();
actor->SetMapper(mapper);
vtkDataSetMapper* mapper2 = vtkDataSetMapper::New();
#if VTK_MAJOR_VERSION <= 5
mapper2->SetInput(vtkgrid);
#else
mapper2->SetInputData(vtkgrid);
#endif
vtkActor* actor2 = vtkActor::New();
actor2->SetMapper(mapper2);
ren1->SetViewport(0.0, 0.0, 0.5, 1.0);
ren2->SetViewport(0.5, 0.0, 1.0, 1.0);
ren2->AddActor(actor2);
// add the actor and start the render loop
ren1->AddActor(actor);
renWin->Render();
inter->Start();
mapper->Delete();
actor->Delete();
ren1->Delete();
mapper2->Delete();
actor2->Delete();
ren2->Delete();
renWin->Delete();
}
#ifndef VTK_HAS_ID_TYPE
// if you get a syntax error hear, then your VTK
// has vtkIdType already defined
typedef long vtkIdType;
#endif
floatMesh::Pointer MeshFromUnstructuredGrid(vtkUnstructuredGrid* grid)
{
// Create a new mesh
floatMesh::Pointer mesh(floatMesh::New());
// Get the points from vtk
vtkPoints* vtkpoints = grid->GetPoints();
int numPoints = vtkpoints->GetNumberOfPoints();
// Create a compatible point container for the mesh
// the mesh is created with a null points container
floatMesh::PointsContainer::Pointer points =
floatMesh::PointsContainer::New();
// Resize the point container to be able to fit the vtk points
points->Reserve(numPoints);
// Set the point container on the mesh
mesh->SetPoints(points);
for(int i =0; i < numPoints; i++)
{
vtkFloatingPointType* apoint = vtkpoints->GetPoint(i);
mesh->SetPoint(i, floatMesh::PointType(apoint));
}
vtkCellArray* vtkcells = grid->GetCells();
floatMesh::CellsContainerPointer cells = floatMesh::CellsContainer::New();
mesh->SetCells(cells);
// extract the cell id's from the vtkUnstructuredGrid
int numcells = vtkcells->GetNumberOfCells();
int* vtkCellTypes = new int[numcells];
int cellId =0;
for(; cellId < numcells; cellId++)
{
vtkCellTypes[cellId] = grid->GetCellType(cellId);
}
cells->Reserve(numcells);
vtkIdType npts;
vtkIdType* pts;
cellId = 0;
for(vtkcells->InitTraversal(); vtkcells->GetNextCell(npts, pts); cellId++)
{
floatMesh::CellAutoPointer c;
switch(vtkCellTypes[cellId])
{
case VTK_TRIANGLE:
{
typedef itk::CellInterface<vtkFloatingPointType, floatMesh::CellTraits> CellInterfaceType;
typedef itk::TriangleCell<CellInterfaceType> TriangleCellType;
TriangleCellType * t = new TriangleCellType;
TriangleCellType::PointIdentifier itkPts[3];
for (int ii = 0; ii < npts; ++ii)
{
itkPts[ii] = static_cast<TriangleCellType::PointIdentifier>(pts[ii]);
}
t->SetPointIds(itkPts);
c.TakeOwnership( t );
break;
}
case VTK_QUAD:
{
typedef itk::CellInterface<vtkFloatingPointType, floatMesh::CellTraits> CellInterfaceType;
typedef itk::QuadrilateralCell<CellInterfaceType> QuadrilateralCellType;
QuadrilateralCellType * t = new QuadrilateralCellType;
QuadrilateralCellType::PointIdentifier itkPts[3];
for (int ii = 0; ii < npts; ++ii)
{
itkPts[ii] = static_cast<QuadrilateralCellType::PointIdentifier>(pts[ii]);
}
t->SetPointIds(itkPts);
c.TakeOwnership( t );
break;
}
case VTK_EMPTY_CELL:
case VTK_VERTEX:
case VTK_POLY_VERTEX:
case VTK_LINE:
case VTK_POLY_LINE:
case VTK_TRIANGLE_STRIP:
case VTK_POLYGON:
case VTK_PIXEL:
case VTK_TETRA:
case VTK_VOXEL:
case VTK_HEXAHEDRON:
case VTK_WEDGE:
case VTK_PYRAMID:
case VTK_PARAMETRIC_CURVE:
case VTK_PARAMETRIC_SURFACE:
default:
std::cerr << "Warning unhandled cell type "
<< vtkCellTypes[cellId] << std::endl;
;
}
mesh->SetCell(cellId, c);
}
mesh->Print(std::cout);
return mesh;
}
class VistVTKCellsClass
{
vtkCellArray* m_Cells;
int* m_LastCell;
int* m_TypeArray;
public:
// typedef the itk cells we are interested in
typedef itk::CellInterface<
floatMesh::PixelType,
floatMesh::CellTraits > CellInterfaceType;
typedef itk::TriangleCell<CellInterfaceType> floatTriangleCell;
typedef itk::QuadrilateralCell<CellInterfaceType> floatQuadrilateralCell;
// Set the vtkCellArray that will be constructed
void SetCellArray(vtkCellArray* a)
{
m_Cells = a;
}
// Set the cell counter pointer
void SetCellCounter(int* i)
{
m_LastCell = i;
}
// Set the type array for storing the vtk cell types
void SetTypeArray(int* i)
{
m_TypeArray = i;
}
// Visit a triangle and create the VTK_TRIANGLE cell
void Visit(unsigned long, floatTriangleCell* t)
{
m_Cells->InsertNextCell(3, (vtkIdType*)t->PointIdsBegin());
m_TypeArray[*m_LastCell] = VTK_TRIANGLE;
(*m_LastCell)++;
}
// Visit a triangle and create the VTK_QUAD cell
void Visit(unsigned long, floatQuadrilateralCell* t)
{
m_Cells->InsertNextCell(4, (vtkIdType*)t->PointIdsBegin());
m_TypeArray[*m_LastCell] = VTK_QUAD;
(*m_LastCell)++;
}
};
typedef itk::CellInterfaceVisitorImplementation<
vtkFloatingPointType, floatMesh::CellTraits,
itk::TriangleCell< itk::CellInterface<floatMesh::PixelType, floatMesh::CellTraits > >,
VistVTKCellsClass> TriangleVisitor;
typedef itk::CellInterfaceVisitorImplementation<
vtkFloatingPointType, floatMesh::CellTraits,
itk::QuadrilateralCell< itk::CellInterface<floatMesh::PixelType, floatMesh::CellTraits > >,
VistVTKCellsClass> QuadrilateralVisitor;
vtkUnstructuredGrid* MeshToUnstructuredGrid(floatMesh* mesh)
{
// Get the number of points in the mesh
int numPoints = mesh->GetNumberOfPoints();
if(numPoints == 0)
{
mesh->Print(std::cerr);
std::cerr << "no points in Grid " << std::endl;
exit(-1);
}
// Create a vtkUnstructuredGrid
vtkUnstructuredGrid* vgrid = vtkUnstructuredGrid::New();
// Create the vtkPoints object and set the number of points
vtkPoints* vpoints = vtkPoints::New();
vpoints->SetNumberOfPoints(numPoints);
// iterate over all the points in the itk mesh filling in
// the vtkPoints object as we go
floatMesh::PointsContainer::Pointer points = mesh->GetPoints();
for(floatMesh::PointsContainer::Iterator i = points->Begin();
i != points->End(); ++i)
{
// Get the point index from the point container iterator
int idx = i->Index();
// Set the vtk point at the index with the the coord array from itk
// itk returns a const pointer, but vtk is not const correct, so
// we have to use a const cast to get rid of the const
vpoints->SetPoint(idx, const_cast<vtkFloatingPointType*>(i->Value().GetDataPointer()));
}
// Set the points on the vtk grid
vgrid->SetPoints(vpoints);
// Now create the cells using the MulitVisitor
// 1. Create a MultiVisitor
floatMesh::CellType::MultiVisitor::Pointer mv =
floatMesh::CellType::MultiVisitor::New();
// 2. Create a triangle and quadrilateral visitor
TriangleVisitor::Pointer tv = TriangleVisitor::New();
QuadrilateralVisitor::Pointer qv = QuadrilateralVisitor::New();
// 3. Set up the visitors
int vtkCellCount = 0; // running counter for current cell being inserted into vtk
int numCells = mesh->GetNumberOfCells();
int *types = new int[numCells]; // type array for vtk
// create vtk cells and estimate the size
vtkCellArray* cells = vtkCellArray::New();
cells->EstimateSize(numCells, 4);
// Set the TypeArray CellCount and CellArray for both visitors
tv->SetTypeArray(types);
tv->SetCellCounter(&vtkCellCount);
tv->SetCellArray(cells);
qv->SetTypeArray(types);
qv->SetCellCounter(&vtkCellCount);
qv->SetCellArray(cells);
// add the visitors to the multivisitor
mv->AddVisitor(tv);
mv->AddVisitor(qv);
// Now ask the mesh to accept the multivisitor which
// will Call Visit for each cell in the mesh that matches the
// cell types of the visitors added to the MultiVisitor
mesh->Accept(mv);
// Now set the cells on the vtk grid with the type array and cell array
vgrid->SetCells(types, cells);
// Clean up vtk objects (no vtkSmartPointer ... )
cells->Delete();
vpoints->Delete();
// return the vtkUnstructuredGrid
return vgrid;
}
int main(int ac, char** av)
{
const char* fname = "c:/blow.vtk";
if(ac > 1 )
{
fname = av[1];
}
vtkDataSetReader* reader = vtkDataSetReader::New();
reader->SetFileName(fname);
reader->SetScalarsName("thickness9");
reader->SetVectorsName("displacement9");
std::cerr << "Begin reading vtk data " << std::endl;
reader->Update();
std::cerr << "End reading vtk data " << std::endl;
vtkDataSet* data = reader->GetOutput();
if( !data )
{
std::cerr << "Data is null" << std::endl;
}
// Convert from vtk to itk
std::cerr << "Begin conversion to itkMesh " << std::endl;
floatMesh::Pointer m =
MeshFromUnstructuredGrid(reader->GetUnstructuredGridOutput());
std::cerr << "End conversion to itkMesh " << std::endl;
// Convert back from itk to vtk
std::cerr << "Begin conversion to vtkUnstructuredGrid " << std::endl;
vtkUnstructuredGrid* grid = MeshToUnstructuredGrid(m);
std::cerr << "End conversion to vtkUnstructuredGrid " << std::endl;
std::cerr << "Begin vtk dispaly of both grids " << std::endl;
Display(grid, reader->GetUnstructuredGridOutput());
vtkDataSetWriter* writer = vtkDataSetWriter::New();
#if VTK_MAJOR_VERSION <= 5
writer->SetInput(grid);
#else
writer->SetInputData(grid);
#endif
writer->SetFileName("./itkblow.vtk");
writer->Update();
#ifdef VTK_USE_ANSI_STDLIB
grid->Print(std::cout);
#endif
reader->Delete();
grid->Delete();
return 0;
}