Mesh Oriented datABase  (version 5.4.1)
Array-based unstructured mesh datastructure
vtkMoabReader.cxx
Go to the documentation of this file.
00001 #include "vtkMoabReader.h"
00002 
00003 #include "SimpleMoab.h"
00004 #include "DataSetConverter.h"
00005 
00006 
00007 #include "vtkNew.h"
00008 #include "vtkObjectFactory.h"
00009 #include "vtkInformation.h"
00010 #include "vtkInformationVector.h"
00011 #include "vtkMultiBlockDataSet.h"
00012 #include "vtkUnstructuredGrid.h"
00013 
00014 #include "vtkPolyData.h"
00015 
00016 
00017 vtkStandardNewMacro(vtkMoabReader)
00018 //------------------------------------------------------------------------------
00019 vtkMoabReader::vtkMoabReader()
00020   {
00021   this->SetNumberOfInputPorts(0);
00022   this->FileName = NULL;
00023   }
00024 
00025 //------------------------------------------------------------------------------
00026 vtkMoabReader::~vtkMoabReader()
00027   {
00028   }
00029 
00030 //------------------------------------------------------------------------------
00031 int vtkMoabReader::RequestInformation(vtkInformation *request,
00032                        vtkInformationVector **inputVector,
00033                        vtkInformationVector *outputVector)
00034 {
00035 
00036   //todo. Walk the file and display all the 2d and 3d elements that the users
00037   //could possibly want to load
00038   return this->Superclass::RequestInformation(request,inputVector,outputVector);
00039 }
00040 
00041 //------------------------------------------------------------------------------
00042 int vtkMoabReader::RequestData(vtkInformation *vtkNotUsed(request),
00043                 vtkInformationVector **vtkNotUsed(inputVector),
00044                 vtkInformationVector *outputVector)
00045 {
00046   //First pass is lets load in all 3d elements in a block called Volumes,
00047   //and load all 2d elements in a block called Surfaces
00048 
00049   vtkInformation* outInfo = outputVector->GetInformationObject(0);
00050   vtkMultiBlockDataSet *output =
00051     vtkMultiBlockDataSet::SafeDownCast(outInfo->Get(vtkDataObject::DATA_OBJECT()));
00052 
00053   vtkNew<vtkMultiBlockDataSet> volumeRoot;
00054   vtkNew<vtkMultiBlockDataSet> boundaryRoot;
00055   vtkNew<vtkMultiBlockDataSet> surfaceRoot;
00056   vtkNew<vtkMultiBlockDataSet> materialRoot;
00057   vtkNew<vtkMultiBlockDataSet> neumannRoot;
00058   vtkNew<vtkMultiBlockDataSet> dirichletRoot;
00059 
00060   const int blockIndex = output->GetNumberOfBlocks();
00061   output->SetBlock(blockIndex,volumeRoot.GetPointer());
00062   output->SetBlock(blockIndex+1,boundaryRoot.GetPointer());
00063   output->SetBlock(blockIndex+2,surfaceRoot.GetPointer());
00064   output->SetBlock(blockIndex+3,materialRoot.GetPointer());
00065   output->SetBlock(blockIndex+4,neumannRoot.GetPointer());
00066   output->SetBlock(blockIndex+5,dirichletRoot.GetPointer());
00067 
00068   //boring work, set the names of the blocks
00069   output->GetMetaData(blockIndex)->Set(vtkCompositeDataSet::NAME(), "Volumes");
00070   output->GetMetaData(blockIndex+1)->Set(vtkCompositeDataSet::NAME(), "Boundary");
00071   output->GetMetaData(blockIndex+2)->Set(vtkCompositeDataSet::NAME(), "Surfaces");
00072   output->GetMetaData(blockIndex+3)->Set(vtkCompositeDataSet::NAME(), "Material");
00073   output->GetMetaData(blockIndex+4)->Set(vtkCompositeDataSet::NAME(), "Neumann Sets");
00074   output->GetMetaData(blockIndex+5)->Set(vtkCompositeDataSet::NAME(), "Dirichlet Sets");
00075 
00076   smoab::GeomTag geom3Tag(3);
00077   smoab::GeomTag geom2Tag(2);
00078   smoab::GeomTag geom1Tag(1);
00079   smoab::MaterialTag matTag;
00080   smoab::NeumannTag neTag;
00081   smoab::DirichletTag diTag;
00082 
00083   smoab::Interface interface(this->FileName);
00084   this->CreateSubBlocks(volumeRoot, &interface, &geom3Tag);
00085   this->CreateSubBlocks(boundaryRoot, &interface, &geom3Tag, &geom2Tag);
00086 
00087   this->CreateSubBlocks(surfaceRoot, &interface, &geom2Tag);
00088   this->CreateSubBlocks(materialRoot, &interface, &matTag, &geom3Tag);
00089   this->CreateSubBlocks(neumannRoot, &interface, &neTag);
00090   this->CreateSubBlocks(dirichletRoot, &interface, &diTag);
00091 
00092   return 1;
00093 }
00094 
00095 
00096 //------------------------------------------------------------------------------
00097 void vtkMoabReader::CreateSubBlocks(vtkNew<vtkMultiBlockDataSet> & root,
00098                                     smoab::Interface* interface,
00099                                     smoab::Tag const* parentTag,
00100                                     smoab::Tag const* extractTag)
00101 {
00102   if(!extractTag)
00103     {
00104     extractTag = parentTag;
00105     }
00106   //basic premise: query the database for all tagged elements and create a new
00107   //multiblock elemenent for each
00108   smoab::DataSetConverter converter(*interface,extractTag);
00109   converter.readMaterialIds(true);
00110   converter.readProperties(true);
00111 
00112   smoab::EntityHandle rootHandle = interface->getRoot();
00113   smoab::Range parents = interface->findEntityRootParents(rootHandle);
00114   smoab::Range dimEnts = interface->findEntitiesWithTag(*parentTag,
00115                                                        rootHandle);
00116 
00117   smoab::Range geomParents = smoab::intersect(parents,dimEnts);
00118 
00119   parents.clear(); //remove this range as it is unneeded
00120   dimEnts.clear();
00121 
00122   //now each item in range can be extracted into a different grid
00123   typedef smoab::Range::iterator iterator;
00124   vtkIdType index = 0;
00125   for(iterator i=geomParents.begin(); i != geomParents.end(); ++i)
00126     {
00127     vtkNew<vtkUnstructuredGrid> block;
00128     //fill the dataset with geometry and properties
00129     converter.fill(*i, block.GetPointer(),index);
00130 
00131     //only add it if we have cells found
00132     if(block->GetNumberOfCells() > 0)
00133       {
00134       root->SetBlock(index,block.GetPointer());
00135       std::string name = interface->name(*i);
00136       if(name.size() > 0)
00137         {
00138         root->GetMetaData(index)->Set(vtkCompositeDataSet::NAME(), name.c_str());
00139         }
00140       ++index;
00141       }
00142     }
00143 }
00144 
00145 //------------------------------------------------------------------------------
00146 void vtkMoabReader::PrintSelf(ostream& os, vtkIndent indent)
00147 {
00148   this->Superclass::PrintSelf(os,indent);
00149 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines