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