MOAB: Mesh Oriented datABase  (version 5.2.1)
DataSetConverter.h
Go to the documentation of this file.
00001 #ifndef __smoab_DataSetConverter_h
00002 #define __smoab_DataSetConverter_h
00003 
00004 #include "SimpleMoab.h"
00005 #include "CellSets.h"
00006 #include "detail/LoadGeometry.h"
00007 #include "detail/ReadSparseTag.h"
00008 
00009 #include <vtkCellData.h>
00010 #include <vtkDoubleArray.h>
00011 #include <vtkFieldData.h>
00012 #include <vtkIntArray.h>
00013 #include <vtkIdTypeArray.h>
00014 #include <vtkNew.h>
00015 
00016 #include <algorithm>
00017 
00018 namespace smoab{
00019 class DataSetConverter
00020 {
00021   const smoab::Interface& Interface;
00022   moab::Interface* Moab;
00023   const smoab::Tag *Tag;
00024   bool ReadMaterialIds;
00025   bool ReadProperties;
00026   std::string MaterialName;
00027 
00028 public:
00029   DataSetConverter(const smoab::Interface& interface, const smoab::Tag* tag):
00030     Interface(interface),
00031     Moab(interface.Moab),
00032     Tag(tag),
00033     ReadMaterialIds(false),
00034     ReadProperties(false)
00035     {
00036     }
00037 
00038   void readMaterialIds(bool add) { this->ReadMaterialIds = add; }
00039   bool readMaterialIds() const { return this->ReadMaterialIds; }
00040 
00041   void readProperties(bool readProps) { this->ReadProperties = readProps; }
00042   bool readProperties() const { return this->ReadProperties; }
00043 
00044   //----------------------------------------------------------------------------
00045   //given a range of entity handles merge them all into a single unstructured
00046   //grid. Currently doesn't support reading properties.
00047   //Will read in material ids,  if no material id is assigned to an entity,
00048   //its cells will be given an unique id
00049   template<typename VTKGridType>
00050   bool fill(const smoab::Range& entities,
00051             VTKGridType* grid) const
00052     {
00053     //create a helper datastructure which can determines all the unique point ids
00054     //and converts moab connecitvity info to vtk connectivity
00055 
00056 
00057     //get all the cells for each parent entity and create
00058     // an entity set of those items
00059     int dim = this->Tag->value();
00060     typedef smoab::Range::const_iterator iterator;
00061     smoab::CellSets entitySets;
00062     for(iterator i=entities.begin(); i!= entities.end(); ++i)
00063       {
00064       smoab::Range entitiesCells;
00065       if(this->Tag->isComparable())
00066         {
00067         //if we are comparable only find the cells that match our tags dimension
00068         entitiesCells = this->Interface.findEntitiesWithDimension(*i,dim,true);
00069         }
00070       else
00071         {
00072         entitiesCells = this->Interface.findHighestDimensionEntities(*i,true);
00073         }
00074       smoab::CellSet set(*i,entitiesCells);
00075       entitySets.push_back(set);
00076       }
00077 
00078     moab::Range cells = smoab::getAllCells(entitySets);
00079 
00080 
00081     //convert the datastructure from a list of cells to a vtk data set
00082     detail::LoadGeometry loadGeom(cells,dim,this->Interface);
00083     loadGeom.fill(grid);
00084 
00085     if(this->readMaterialIds())
00086       {
00087 
00088       detail::ReadSparseTag materialTagReading(entitySets,
00089                                                cells,
00090                                                this->Interface);
00091 
00092       smoab::MaterialTag mtag;
00093       vtkNew<vtkIntArray> materials;
00094       materials->SetName(mtag.name());
00095       materialTagReading.fill(materials.GetPointer(),&mtag);
00096       grid->GetCellData()->AddArray(materials.GetPointer());
00097       }
00098 
00099     //by default we always try to load the default tag
00100     detail::ReadSparseTag sTagReading(entitySets,
00101                                       cells,
00102                                       this->Interface);
00103 
00104     vtkNew<vtkIntArray> sparseTagData;
00105     sparseTagData->SetName(this->Tag->name());
00106     sTagReading.fill(sparseTagData.GetPointer(),this->Tag);
00107     grid->GetCellData()->AddArray(sparseTagData.GetPointer());
00108 
00109     return true;
00110     }
00111 
00112   //----------------------------------------------------------------------------
00113   //given a single entity handle create a unstructured grid from it.
00114   //optional third parameter is the material id to use if readMaterialIds
00115   //is on, and no material sparse tag is found for this entity
00116   template<typename VTKGridType>
00117   bool fill(const smoab::EntityHandle& entity,
00118             VTKGridType* grid,
00119             const int materialId=0) const
00120     {
00121     //create a helper datastructure which can determines all the unique point ids
00122     //and converts moab connecitvity info to vtk connectivity
00123 
00124     smoab::Range cells;
00125     int dim = this->Tag->value();
00126     if(this->Tag->isComparable())
00127       {
00128       //if we are comparable only find the cells that match our tags dimension
00129       cells = this->Interface.findEntitiesWithDimension(entity,dim,true);
00130       }
00131     else
00132       {
00133       //load subentities
00134       cells = this->Interface.findHighestDimensionEntities(entity,true);
00135       }
00136 
00137     //convert the datastructure from a list of cells to a vtk data set
00138     detail::LoadGeometry loadGeom(cells,dim,this->Interface);
00139     loadGeom.fill(grid);
00140 
00141     const smoab::Range& points = loadGeom.moabPoints();
00142 
00143     if(this->readProperties())
00144       {
00145       this->readProperties(cells,grid->GetCellData());
00146       this->readProperties(points,grid->GetPointData());
00147       }
00148 
00149     smoab::CellSets cellSets;
00150     smoab::CellSet set(entity,cells);
00151     cellSets.push_back(set);
00152     if(this->readMaterialIds())
00153       {
00154       smoab::MaterialTag mtag;
00155       detail::ReadSparseTag materialTagReading(cellSets,
00156                                                cells,
00157                                                this->Interface);
00158 
00159       vtkNew<vtkIntArray> materials;
00160       materials->SetName(mtag.name());
00161       materialTagReading.fill(materials.GetPointer(),&mtag);
00162       grid->GetCellData()->AddArray(materials.GetPointer());
00163 
00164       }
00165 
00166     //by default we always try to load the default tag
00167     detail::ReadSparseTag sTagReading(cellSets,
00168                                       cells,
00169                                       this->Interface);
00170 
00171     vtkNew<vtkIntArray> sparseTagData;
00172     sparseTagData->SetName(this->Tag->name());
00173     sTagReading.fill(sparseTagData.GetPointer(),this->Tag);
00174     grid->GetCellData()->AddArray(sparseTagData.GetPointer());
00175 
00176     return true;
00177     }
00178 
00179 private:
00180   //----------------------------------------------------------------------------
00181   void readProperties(smoab::Range const& entities,
00182                            vtkFieldData* field) const
00183     {
00184     if(entities.empty()) { return; }
00185 
00186     //so we get all the tags and parse out the sparse and dense tags
00187     //that we support
00188     typedef std::vector<moab::Tag>::const_iterator iterator;
00189     std::vector<moab::Tag> tags;
00190     this->Moab->tag_get_tags_on_entity(entities.front(),tags);
00191 
00192     this->readDenseTags(tags,entities,field);
00193     }
00194 
00195   //----------------------------------------------------------------------------
00196   void readDenseTags(std::vector<moab::Tag> &tags,
00197                      smoab::Range const& entities,
00198                      vtkFieldData* field) const
00199     {
00200     typedef std::vector<moab::Tag>::const_iterator iterator;
00201 
00202     for(iterator i=tags.begin();i!=tags.end();++i)
00203       {
00204       moab::TagType tagType;
00205       moab::DataType tagDataType;
00206 
00207       this->Moab->tag_get_type(*i,tagType);
00208       this->Moab->tag_get_data_type(*i,tagDataType);
00209 
00210       //make sure it is only dense
00211       if(tagType != moab::MB_TAG_DENSE)
00212         {
00213         continue;
00214         }
00215       //and only integer and double
00216       if(tagDataType != moab::MB_TYPE_DOUBLE &&
00217          tagDataType != moab::MB_TYPE_INTEGER)
00218         {
00219         //unsupported type, skip to next tag
00220         continue;
00221         }
00222 
00223       //read the name of the tag
00224       std::string name;
00225       name.reserve(32);
00226       this->Moab->tag_get_name(*i,name);
00227 
00228       //read the number of components of the tag
00229       int numComps = 1;
00230 
00231       this->Moab->tag_get_length(*i,numComps);
00232 
00233       //read the data if it is one of the two types we support
00234       int size = entities.size();
00235       if(tagDataType == moab::MB_TYPE_DOUBLE)
00236         {
00237         vtkNew<vtkDoubleArray> array;
00238         array->SetName(name.c_str());
00239         array->SetNumberOfComponents(numComps);
00240         array->SetNumberOfTuples(size);
00241 
00242         //read directly into the double array
00243         this->Moab->tag_get_data(*i,entities,
00244                                  array->GetVoidPointer(0));
00245         field->AddArray(array.GetPointer());
00246         }
00247       else if(tagDataType == moab::MB_TYPE_INTEGER)
00248         {
00249         vtkNew<vtkIntArray> array;
00250         array->SetName(name.c_str());
00251         array->SetNumberOfComponents(numComps);
00252         array->SetNumberOfTuples(size);
00253 
00254         //read directly into the double array
00255         this->Moab->tag_get_data(*i,entities,
00256                                  array->GetVoidPointer(0));
00257         field->AddArray(array.GetPointer());
00258         }
00259       else
00260         {
00261         }
00262       }
00263     }
00264 
00265 
00266 };
00267 
00268 }
00269 #endif // __smoab_DataSetConverter_h
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines