Mesh Oriented datABase  (version 5.4.1)
Array-based unstructured mesh datastructure
LoadPoly.h
Go to the documentation of this file.
00001 #ifndef __smoab_detail_LoadPoly_h
00002 #define __smoab_detail_LoadPoly_h
00003 
00004 #include "SimpleMoab.h"
00005 #include "LinearCellConnectivity.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 
00015 namespace smoab
00016 {
00017 namespace detail
00018 {
00019 
00020     //basic class that given a range of moab cells will convert them to a
00021     //vtk poly data. If given any high order cells will do a simple linearization
00022     //of the cells
00023     class LoadPoly
00024     {
00025         moab::Interface* Interface;
00026         smoab::detail::LinearCellConnectivity CellConn;
00027         smoab::Range Points;
00028 
00029       public:
00030         //----------------------------------------------------------------------------
00031         //warning cells input to constructor is held by reference
00032         LoadPoly( const smoab::Range& cells, const smoab::Interface& interface )
00033             : Interface( interface.Moab ), CellConn( cells, interface.Moab ), Points()
00034         {
00035         }
00036 
00037         const smoab::Range& moabPoints() const
00038         {
00039             return this->Points;
00040         }
00041 
00042         //----------------------------------------------------------------------------
00043         //todo: have support for using only a subsection of the input cells
00044         void fill( vtkPolyData* dataSet )
00045         {
00046             //now that CellConn has all the cells properly stored, lets fixup
00047             //the ids so that they start at zero and keep the same logical ordering
00048             //as before.
00049             vtkIdType numCells, connLen;
00050             this->CellConn.compactIds( numCells, connLen );
00051             this->addGridsTopology( dataSet, numCells, connLen );
00052             this->addCoordinates( dataSet );
00053         }
00054 
00055       private:
00056         //----------------------------------------------------------------------------
00057         void addCoordinates( vtkPointSet* grid )
00058         {
00059             //this is sorta of hackish as moabPoints is only valid
00060             //after compactIds has been called
00061             this->CellConn.moabPoints( this->Points );
00062 
00063             //since the smoab::range are always unique and sorted
00064             //we can use the more efficient coords_iterate
00065             //call in moab, which returns moab internal allocated memory
00066             vtkNew< vtkPoints > newPoints;
00067             newPoints->SetDataTypeToDouble();
00068             newPoints->SetNumberOfPoints( this->Points.size() );
00069 
00070             //need a pointer to the allocated vtkPoints memory so that we
00071             //don't need to use an extra copy and we can bypass all vtk's check
00072             //on out of bounds
00073             double* rawPoints = static_cast< double* >( newPoints->GetVoidPointer( 0 ) );
00074             this->Interface->get_coords( this->Points, rawPoints );
00075 
00076             grid->SetPoints( newPoints.GetPointer() );
00077         }
00078 
00079         //----------------------------------------------------------------------------
00080         void addGridsTopology( vtkPolyData* data, vtkIdType numCells, vtkIdType numConnectivity ) const
00081         {
00082             //correct the connectivity size to account for the vtk padding
00083             const vtkIdType vtkConnectivity = numCells + numConnectivity;
00084 
00085             vtkNew< vtkIdTypeArray > cellArray;
00086             cellArray->SetNumberOfValues( vtkConnectivity );
00087 
00088             vtkIdType* rawArray = static_cast< vtkIdType* >( cellArray->GetVoidPointer( 0 ) );
00089             this->CellConn.copyToVtkCellInfo( rawArray );
00090 
00091             vtkNew< vtkCellArray > cells;
00092             cells->SetCells( numCells, cellArray.GetPointer() );
00093 
00094             data->SetPolys( cells.GetPointer() );
00095         }
00096     };
00097 
00098 }  // namespace detail
00099 }  // namespace smoab
00100 
00101 #endif  // __smoab_detail_LoadPoly_h
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines