Chaste  Build::
CellPopulationActorGenerator.cpp
1 /*
2 
3 Copyright (c) 2005-2016, University of Oxford.
4  All rights reserved.
5 
6  University of Oxford means the Chancellor, Masters and Scholars of the
7  University of Oxford, having an administrative office at Wellington
8  Square, Oxford OX1 2JD, UK.
9 
10  This file is CellPopulation of Chaste.
11 
12  Redistribution and use in source and binary forms, with or without
13  modification, are permitted provided that the following conditions are met:
14  * Redistributions of source code must retain the above copyright notice,
15  this list of conditions and the following disclaimer.
16  * Redistributions in binary form must reproduce the above copyright notice,
17  this list of conditions and the following disclaimer in the documentation
18  and/or other materials provided with the distribution.
19  * Neither the name of the University of Oxford nor the names of its
20  contributors may be used to endorse or promote products derived from this
21  software without specific prior written permission.
22 
23  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
24  AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
25  IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A CellPopulationICULAR PURPOSE
26  ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
27  LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
28  CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE
29  GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
30  HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
31  LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
32  OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
33 
34  */
35 
36 #include <boost/filesystem.hpp>
37 #include <boost/lexical_cast.hpp>
38 #define _BACKWARD_BACKWARD_WARNING_H 1 //Cut out the vtk deprecated warning
39 #include <vtkPoints.h>
40 #include <vtkPolyData.h>
41 #include <vtkPolyDataMapper.h>
42 #include <vtkActor.h>
43 #include <vtkProperty.h>
44 #include <vtkUnsignedCharArray.h>
45 #if VTK_MAJOR_VERSION > 5
46  #include <vtkNamedColors.h>
47 #endif
48 #include <vtkSphereSource.h>
49 #include <vtkGlyph3D.h>
50 #include <vtkGlyph2D.h>
51 #include <vtkCubeAxesActor2D.h>
52 #include <vtkImageData.h>
53 #include <vtkGeometryFilter.h>
54 #include <vtkTubeFilter.h>
55 #include <vtkExtractEdges.h>
56 #include <vtkUnstructuredGrid.h>
57 #include <vtkCell.h>
58 #include <vtkPolygon.h>
59 #include <vtkIdList.h>
60 #include <vtkFeatureEdges.h>
61 #include "UblasIncludes.hpp"
62 #include "UblasVectorInclude.hpp"
63 #include "Exception.hpp"
64 #include "BaseUnits.hpp"
65 #include "VesselNetworkWriter.hpp"
66 
67 #include "CellPopulationActorGenerator.hpp"
68 
69 
70 template<unsigned DIM>
72  : AbstractActorGenerator<DIM>(),
73  mpCellPopulation()
74 {
75 
76 }
77 
78 template<unsigned DIM>
80 {
81 
82 }
83 
84 template<unsigned DIM>
85 void CellPopulationActorGenerator<DIM>::AddActor(vtkSmartPointer<vtkRenderer> pRenderer)
86 {
88  {
89  // Check the cell population type
90  if(boost::dynamic_pointer_cast<MeshBasedCellPopulation<DIM> >(mpCellPopulation))
91  {
93  }
94  else
95  {
96  // Fall back to a centre based default
97  vtkSmartPointer<vtkPoints> p_points = vtkSmartPointer<vtkPoints>::New();
98  vtkSmartPointer<vtkPolyData> p_polydata = vtkSmartPointer<vtkPolyData>::New();
99 
100  for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = mpCellPopulation->Begin();
101  cell_iter != mpCellPopulation->End(); ++cell_iter)
102  {
103  c_vector<double, DIM> centre = mpCellPopulation->GetLocationOfCellCentre(*cell_iter);
104  if(DIM==3)
105  {
106  p_points->InsertNextPoint(centre[0], centre[1], centre[2]);
107  }
108  else
109  {
110  p_points->InsertNextPoint(centre[0], centre[1], 0.0);
111  }
112  }
113  p_polydata->SetPoints(p_points);
114 
115  vtkSmartPointer<vtkSphereSource> p_spheres = vtkSmartPointer<vtkSphereSource>::New();
116  p_spheres->SetRadius(0.5);
117  p_spheres->SetPhiResolution(16);
118  p_spheres->SetThetaResolution(16);
119 
120  vtkSmartPointer<vtkGlyph3D> p_glyph = vtkSmartPointer<vtkGlyph3D>::New();
121  #if VTK_MAJOR_VERSION <= 5
122  p_glyph->SetInput(p_polydata);
123  p_glyph->SetSource(p_spheres->GetOutput());
124  #else
125  p_glyph->SetInputData(p_polydata);
126  p_glyph->SetSourceConnection(p_spheres->GetOutputPort());
127  #endif
128  p_glyph->ClampingOff();
129  p_glyph->SetScaleModeToScaleByScalar();
130  p_glyph->SetScaleFactor(1.0);
131  p_glyph->Update();
132 
133  vtkSmartPointer<vtkPolyDataMapper> p_mapper = vtkSmartPointer<vtkPolyDataMapper>::New();
134  #if VTK_MAJOR_VERSION <= 5
135  p_mapper->SetInput(p_glyph->GetOutput());
136  #else
137  p_mapper->SetInputData(p_glyph->GetOutput());
138  #endif
139  p_mapper->ScalarVisibilityOff();
140 
141  vtkSmartPointer<vtkActor> p_actor = vtkSmartPointer<vtkActor>::New();
142  p_actor->SetMapper(p_mapper);
143  p_actor->GetProperty()->SetColor(0,1,0);
144  pRenderer->AddActor(p_actor);
145  }
146  }
147 }
148 
149 template<unsigned DIM>
150 void CellPopulationActorGenerator<DIM>::AddMeshBasedCellPopulationActor(vtkSmartPointer<vtkRenderer> pRenderer)
151 {
152  boost::shared_ptr<MeshBasedCellPopulation<DIM> > p_cell_population = boost::dynamic_pointer_cast<MeshBasedCellPopulation<DIM> >(mpCellPopulation);
153 
154  if(!p_cell_population)
155  {
156  EXCEPTION("Could not cast mesh to MeshBased type.");
157  }
158 
159  // Add the voronoi mesh
160  p_cell_population->CreateVoronoiTessellation();
161  vtkSmartPointer<vtkUnstructuredGrid> p_voronoi_grid = vtkSmartPointer<vtkUnstructuredGrid>::New();
162 
163  if(p_cell_population->GetVoronoiTessellation() != NULL)
164  {
165  vtkSmartPointer<vtkPoints> p_points = vtkSmartPointer<vtkPoints>::New();
166  p_points->GetData()->SetName("Vertex positions");
167  for (unsigned node_num=0; node_num<p_cell_population->GetVoronoiTessellation()->GetNumNodes(); node_num++)
168  {
169  c_vector<double, DIM> position = p_cell_population->GetVoronoiTessellation()->GetNode(node_num)->rGetLocation();
170  if (DIM==2)
171  {
172  p_points->InsertPoint(node_num, position[0], position[1], 0.0);
173  }
174  else
175  {
176  p_points->InsertPoint(node_num, position[0], position[1], position[2]);
177  }
178  }
179  p_voronoi_grid->SetPoints(p_points);
180 
181  for (typename VertexMesh<DIM,DIM>::VertexElementIterator iter = p_cell_population->GetVoronoiTessellation()->GetElementIteratorBegin();
182  iter != p_cell_population->GetVoronoiTessellation()->GetElementIteratorEnd(); ++iter)
183  {
184  vtkSmartPointer<vtkCell> p_cell;
185  if (DIM == 2)
186  {
187  p_cell = vtkSmartPointer<vtkPolygon>::New();
188  }
189  else
190  {
191  p_cell = vtkSmartPointer<vtkConvexPointSet>::New();
192  }
193  vtkSmartPointer<vtkIdList> p_cell_id_list = p_cell->GetPointIds();
194  p_cell_id_list->SetNumberOfIds(iter->GetNumNodes());
195  for (unsigned j=0; j<iter->GetNumNodes(); ++j)
196  {
197  p_cell_id_list->SetId(j, iter->GetNodeGlobalIndex(j));
198  }
199  p_voronoi_grid->InsertNextCell(p_cell->GetCellType(), p_cell_id_list);
200  }
201  }
202 
203  vtkSmartPointer<vtkGeometryFilter> p_geom_filter = vtkSmartPointer<vtkGeometryFilter>::New();
204  #if VTK_MAJOR_VERSION <= 5
205  p_geom_filter->SetInput(p_voronoi_grid);
206  #else
207  p_geom_filter->SetInputData(p_voronoi_grid);
208  #endif
209 
210  vtkSmartPointer<vtkPolyDataMapper> p_mapper = vtkSmartPointer<vtkPolyDataMapper>::New();
211  p_mapper->SetInputConnection(p_geom_filter->GetOutputPort());
212  p_mapper->ScalarVisibilityOff();
213 
214  vtkSmartPointer<vtkActor> p_actor = vtkSmartPointer<vtkActor>::New();
215  p_actor->SetMapper(p_mapper);
216  p_actor->GetProperty()->SetColor(1,1,0);
217  p_actor->GetProperty()->SetOpacity(0.8);
218  pRenderer->AddActor(p_actor);
219 
220  vtkSmartPointer<vtkFeatureEdges> p_voronoi_extract_edges = vtkSmartPointer<vtkFeatureEdges>::New();
221  p_voronoi_extract_edges->SetInputConnection(p_geom_filter->GetOutputPort());
222  p_voronoi_extract_edges->SetFeatureEdges(false);
223  p_voronoi_extract_edges->SetBoundaryEdges(true);
224  p_voronoi_extract_edges->SetManifoldEdges(true);
225  p_voronoi_extract_edges->SetNonManifoldEdges(false);
226 
227  vtkSmartPointer<vtkTubeFilter> p_voronoi_tubes = vtkSmartPointer<vtkTubeFilter>::New();
228  p_voronoi_tubes->SetInputConnection(p_voronoi_extract_edges->GetOutputPort());
229  p_voronoi_tubes->SetRadius(0.006);
230  p_voronoi_tubes->SetNumberOfSides(12);
231 
232  vtkSmartPointer<vtkPolyDataMapper> p_voronoi_tube_mapper = vtkSmartPointer<vtkPolyDataMapper>::New();
233  p_voronoi_tube_mapper->SetInputConnection(p_voronoi_tubes->GetOutputPort());
234  p_voronoi_tube_mapper->ScalarVisibilityOff();
235 
236  vtkSmartPointer<vtkActor> p_voronoi_tube_actor = vtkSmartPointer<vtkActor>::New();
237  p_voronoi_tube_actor->SetMapper(p_voronoi_tube_mapper);
238  p_voronoi_tube_actor->GetProperty()->SetColor(0,0,0);
239  pRenderer->AddActor(p_voronoi_tube_actor);
240 
241 // // Do the mutable mesh
242 // //Make the local mesh into a VtkMesh
243 // vtkSmartPointer<vtkUnstructuredGrid> p_mutable_grid = vtkSmartPointer<vtkUnstructuredGrid>::New();
244 // vtkSmartPointer<vtkPoints> p_points = vtkSmartPointer<vtkPoints>::New();
245 // p_points->GetData()->SetName("Vertex positions");
246 //
247 // for (typename AbstractMesh<DIM,DIM>::NodeIterator node_iter = p_cell_population->rGetMesh().GetNodeIteratorBegin();
248 // node_iter != p_cell_population->rGetMesh().GetNodeIteratorEnd();
249 // ++node_iter)
250 // {
251 // c_vector<double, DIM> current_item = node_iter->rGetLocation();
252 // if (DIM == 3)
253 // {
254 // p_points->InsertNextPoint(current_item[0], current_item[1], current_item[2]);
255 // }
256 // else if (DIM == 2)
257 // {
258 // p_points->InsertNextPoint(current_item[0], current_item[1], 0.0);
259 // }
260 // else // (DIM == 1)
261 // {
262 // p_points->InsertNextPoint(current_item[0], 0.0, 0.0);
263 // }
264 // }
265 //
266 // p_mutable_grid->SetPoints(p_points);
267 //
268 // for (typename AbstractTetrahedralMesh<DIM,DIM>::ElementIterator elem_iter = p_cell_population->rGetMesh().GetElementIteratorBegin();
269 // elem_iter != p_cell_population->rGetMesh().GetElementIteratorEnd();
270 // ++elem_iter)
271 // {
272 //
273 // vtkSmartPointer<vtkCell> p_cell;
274 // ///\todo This ought to look exactly like the other MakeVtkMesh
275 // if (DIM == 3)
276 // {
277 // p_cell = vtkSmartPointer<vtkTetra>::New();
278 // }
279 // else if (DIM == 2)
280 // {
281 // p_cell = vtkSmartPointer<vtkTriangle>::New();
282 // }
283 // else //(DIM == 1)
284 // {
285 // p_cell = vtkSmartPointer<vtkLine>::New();
286 // }
287 // vtkSmartPointer<vtkIdList> p_cell_id_list = p_cell->GetPointIds();
288 // for (unsigned j = 0; j < DIM+1; ++j)
289 // {
290 // unsigned global_node_index = elem_iter->GetNodeGlobalIndex(j);
291 // p_cell_id_list->SetId(j, global_node_index);
292 // }
293 // p_mutable_grid->InsertNextCell(p_cell->GetCellType(), p_cell_id_list);
294 // }
295 //
296 // vtkSmartPointer<vtkGeometryFilter> p_mutable_geom_filter = vtkSmartPointer<vtkGeometryFilter>::New();
297 // #if VTK_MAJOR_VERSION <= 5
298 // p_mutable_geom_filter->SetInput(p_mutable_grid);
299 // #else
300 // p_mutable_geom_filter->SetInputData(p_mutable_grid);
301 // #endif
302 //
303 // vtkSmartPointer<vtkFeatureEdges> p_extract_edges = vtkSmartPointer<vtkFeatureEdges>::New();
304 // p_extract_edges->SetInputConnection(p_mutable_geom_filter->GetOutputPort());
305 // p_extract_edges->SetFeatureEdges(false);
306 // p_extract_edges->SetBoundaryEdges(true);
307 // p_extract_edges->SetManifoldEdges(true);
308 // p_extract_edges->SetNonManifoldEdges(false);
309 //
310 // vtkSmartPointer<vtkTubeFilter> p_mutable_tubes = vtkSmartPointer<vtkTubeFilter>::New();
311 // p_mutable_tubes->SetInputConnection(p_extract_edges->GetOutputPort());
312 // p_mutable_tubes->SetRadius(0.02);
313 // p_mutable_tubes->SetNumberOfSides(12);
314 //
315 // vtkSmartPointer<vtkPolyDataMapper> p_mutable_mapper = vtkSmartPointer<vtkPolyDataMapper>::New();
316 // p_mutable_mapper->SetInputConnection(p_mutable_tubes->GetOutputPort());
317 //
318 // vtkSmartPointer<vtkActor> p_mutable_actor = vtkSmartPointer<vtkActor>::New();
319 // p_mutable_actor->SetMapper(p_mutable_mapper);
320 // p_mutable_actor->GetProperty()->SetColor(1,1,1);
321 // mpRenderer->AddActor(p_mutable_actor);
322 
323  // Render the centres
324  vtkSmartPointer<vtkPoints> p_centres_points = vtkSmartPointer<vtkPoints>::New();
325  vtkSmartPointer<vtkPolyData> p_polydata = vtkSmartPointer<vtkPolyData>::New();
326 
327  for (typename AbstractCellPopulation<DIM>::Iterator cell_iter = mpCellPopulation->Begin();
328  cell_iter != mpCellPopulation->End(); ++cell_iter)
329  {
330  c_vector<double, DIM> centre = mpCellPopulation->GetLocationOfCellCentre(*cell_iter);
331  if(DIM==3)
332  {
333  p_centres_points->InsertNextPoint(centre[0], centre[1], centre[2]);
334  }
335  else
336  {
337  p_centres_points->InsertNextPoint(centre[0], centre[1], 0.0);
338  }
339  }
340  p_polydata->SetPoints(p_centres_points);
341 
342  vtkSmartPointer<vtkSphereSource> p_spheres = vtkSmartPointer<vtkSphereSource>::New();
343  p_spheres->SetRadius(0.1);
344  p_spheres->SetPhiResolution(16);
345  p_spheres->SetThetaResolution(16);
346 
347  vtkSmartPointer<vtkGlyph3D> p_glyph = vtkSmartPointer<vtkGlyph3D>::New();
348  #if VTK_MAJOR_VERSION <= 5
349  p_glyph->SetInput(p_polydata);
350  p_glyph->SetSource(p_spheres->GetOutput());
351  #else
352  p_glyph->SetInputData(p_polydata);
353  p_glyph->SetSourceConnection(p_spheres->GetOutputPort());
354  #endif
355  p_glyph->ClampingOff();
356  p_glyph->SetScaleModeToScaleByScalar();
357  p_glyph->SetScaleFactor(1.0);
358  p_glyph->Update();
359 
360  vtkSmartPointer<vtkPolyDataMapper> p_centres_mapper = vtkSmartPointer<vtkPolyDataMapper>::New();
361  #if VTK_MAJOR_VERSION <= 5
362  p_centres_mapper->SetInput(p_glyph->GetOutput());
363  #else
364  p_centres_mapper->SetInputData(p_glyph->GetOutput());
365  #endif
366  p_centres_mapper->ScalarVisibilityOff();
367 
368  vtkSmartPointer<vtkActor> p_centres_actor = vtkSmartPointer<vtkActor>::New();
369  p_centres_actor->SetMapper(p_centres_mapper);
370  p_centres_actor->GetProperty()->SetColor(0.6,0.0,0.0);
371  pRenderer->AddActor(p_centres_actor);
372 
373 }
374 
375 template<unsigned DIM>
376 void CellPopulationActorGenerator<DIM>::SetCellPopulation(boost::shared_ptr<AbstractCellPopulation<DIM> > pCellPopulation)
377 {
378  this->mpCellPopulation = pCellPopulation;
379 }
380 
381 template class CellPopulationActorGenerator<2>;
382 template class CellPopulationActorGenerator<3>;
void SetCellPopulation(boost::shared_ptr< AbstractCellPopulation< DIM > > pCellPopulation)
void AddMeshBasedCellPopulationActor(vtkSmartPointer< vtkRenderer > pRenderer)
boost::shared_ptr< AbstractCellPopulation< DIM > > mpCellPopulation
void AddActor(vtkSmartPointer< vtkRenderer > pRenderer)