MOAB: Mesh Oriented datABase
(version 5.2.1)
|
00001 #ifndef __smoab_detail_LoadGeometry_h 00002 #define __smoab_detail_LoadGeometry_h 00003 00004 #include "SimpleMoab.h" 00005 #include "MixedCellConnectivity.h" 00006 00007 #include <vtkCellArray.h> 00008 #include <vtkIdTypeArray.h> 00009 #include <vtkNew.h> 00010 #include <vtkPointData.h> 00011 #include <vtkPoints.h> 00012 #include <vtkPointSet.h> 00013 #include <vtkPolyData.h> 00014 #include <vtkUnsignedCharArray.h> 00015 #include <vtkUnstructuredGrid.h> 00016 namespace smoab { namespace detail { 00017 00018 namespace { 00019 00020 //empty definition, only allow specializations to compile 00021 template<typename VTKGridType> struct setCells; 00022 00023 template<> struct setCells<vtkUnstructuredGrid> 00024 { 00025 typedef vtkUnstructuredGrid VTKGridType; 00026 void operator()(VTKGridType* grid, int, 00027 vtkUnsignedCharArray* types, 00028 vtkIdTypeArray* locations, 00029 vtkCellArray* cells) const 00030 { 00031 grid->SetCells(types,locations,cells,NULL,NULL); 00032 } 00033 }; 00034 00035 template<> struct setCells<vtkPolyData> 00036 { 00037 typedef vtkPolyData VTKGridType; 00038 void operator()(VTKGridType* grid, int dim, 00039 vtkUnsignedCharArray*, 00040 vtkIdTypeArray*, 00041 vtkCellArray* cells) const 00042 { 00043 if( dim == 0) 00044 { 00045 grid->SetVerts(cells); 00046 } 00047 else if( dim == 1 ) 00048 { 00049 grid->SetLines(cells); 00050 } 00051 else if( dim == 2) 00052 { 00053 grid->SetPolys(cells); 00054 } 00055 } 00056 }; 00057 00058 } 00059 //basic class that given a range of moab cells will convert them to am unstructured grid. 00060 //holds only references to input cells, so they can't go out of scope 00061 00062 //The Topology tag is used to describe the dimensoniality of the cells we are reading 00063 // 0 == verts, 1 == lines, 2 == tri/quads, 3 = volume elements 00064 class LoadGeometry 00065 { 00066 moab::Interface* Interface; 00067 int TopologyDim; 00068 smoab::detail::MixedCellConnectivity MixConn; 00069 smoab::Range Points; 00070 00071 public: 00072 //---------------------------------------------------------------------------- 00073 //warning cells input to constructor is held by reference 00074 LoadGeometry(const smoab::Range& cells, 00075 int topologyDim, 00076 const smoab::Interface& interface): 00077 Interface(interface.Moab), 00078 TopologyDim(topologyDim), 00079 MixConn(cells,interface.Moab), 00080 Points() 00081 { 00082 } 00083 00084 const smoab::Range& moabPoints() const { return this->Points; } 00085 00086 //---------------------------------------------------------------------------- 00087 //todo: have support for using only a subsection of the input cells 00088 template<typename vtkDataSetType> 00089 void fill(vtkDataSetType *dataSet) 00090 { 00091 //now that mixConn has all the cells properly stored, lets fixup 00092 //the ids so that they start at zero and keep the same logical ordering 00093 //as before. 00094 vtkIdType numCells, connLen; 00095 this->MixConn.compactIds(numCells,connLen); 00096 this->addGridsTopology(dataSet,numCells,connLen); 00097 this->addCoordinates(dataSet); 00098 } 00099 00100 private: 00101 //---------------------------------------------------------------------------- 00102 void addCoordinates(vtkPointSet *grid) 00103 { 00104 //this is sorta of hackish as moabPoints is only valid 00105 //after compactIds has been called 00106 this->MixConn.moabPoints(this->Points); 00107 00108 //since the smoab::range are always unique and sorted 00109 //we can use the more efficient coords_iterate 00110 //call in moab, which returns moab internal allocated memory 00111 vtkNew<vtkPoints> newPoints; 00112 newPoints->SetDataTypeToDouble(); 00113 newPoints->SetNumberOfPoints(this->Points.size()); 00114 00115 //need a pointer to the allocated vtkPoints memory so that we 00116 //don't need to use an extra copy and we can bypass all vtk's check 00117 //on out of bounds 00118 double *rawPoints = static_cast<double*>(newPoints->GetVoidPointer(0)); 00119 this->Interface->get_coords(this->Points,rawPoints); 00120 00121 grid->SetPoints(newPoints.GetPointer()); 00122 } 00123 00124 //---------------------------------------------------------------------------- 00125 template<typename vtkDataSetType> 00126 void addGridsTopology(vtkDataSetType* grid, 00127 vtkIdType numCells, 00128 vtkIdType numConnectivity) const 00129 { 00130 //correct the connectivity size to account for the vtk padding 00131 const vtkIdType vtkConnectivity = numCells + numConnectivity; 00132 00133 vtkNew<vtkIdTypeArray> cellArray; 00134 vtkNew<vtkIdTypeArray> cellLocations; 00135 vtkNew<vtkUnsignedCharArray> cellTypes; 00136 00137 cellArray->SetNumberOfValues(vtkConnectivity); 00138 cellLocations->SetNumberOfValues(numCells); 00139 cellTypes->SetNumberOfValues(numCells); 00140 00141 vtkIdType* rawArray = static_cast<vtkIdType*>(cellArray->GetVoidPointer(0)); 00142 vtkIdType* rawLocations = static_cast<vtkIdType*>(cellLocations->GetVoidPointer(0)); 00143 unsigned char* rawTypes = static_cast<unsigned char*>(cellTypes->GetVoidPointer(0)); 00144 00145 this->MixConn.copyToVtkCellInfo(rawArray,rawLocations,rawTypes); 00146 00147 vtkNew<vtkCellArray> cells; 00148 cells->SetCells(numCells,cellArray.GetPointer()); 00149 00150 00151 setCells<vtkDataSetType>()(grid, 00152 this->TopologyDim, 00153 cellTypes.GetPointer(), 00154 cellLocations.GetPointer(), 00155 cells.GetPointer()); 00156 } 00157 }; 00158 00159 00160 } } 00161 00162 #endif // __smoab_detail_LoadGeometry_h