Mesh Oriented datABase  (version 5.4.1)
Array-based unstructured mesh datastructure
smoab::DataSetConverter Class Reference

#include <DataSetConverter.h>

+ Collaboration diagram for smoab::DataSetConverter:

Public Member Functions

 DataSetConverter (const smoab::Interface &interface, const smoab::Tag *tag)
void readMaterialIds (bool add)
bool readMaterialIds () const
void readProperties (bool readProps)
bool readProperties () const
template<typename VTKGridType >
bool fill (const smoab::Range &entities, VTKGridType *grid) const
template<typename VTKGridType >
bool fill (const smoab::EntityHandle &entity, VTKGridType *grid, const int materialId=0) const

Private Member Functions

void readProperties (smoab::Range const &entities, vtkFieldData *field) const
void readDenseTags (std::vector< moab::Tag > &tags, smoab::Range const &entities, vtkFieldData *field) const

Private Attributes

const smoab::InterfaceInterface
moab::InterfaceMoab
const smoab::TagTag
bool ReadMaterialIds
bool ReadProperties
std::string MaterialName

Detailed Description

Definition at line 20 of file DataSetConverter.h.


Constructor & Destructor Documentation

smoab::DataSetConverter::DataSetConverter ( const smoab::Interface interface,
const smoab::Tag tag 
) [inline]

Definition at line 30 of file DataSetConverter.h.

        : Interface( interface ), Moab( interface.Moab ), Tag( tag ), ReadMaterialIds( false ), ReadProperties( false )
    {
    }

Member Function Documentation

template<typename VTKGridType >
bool smoab::DataSetConverter::fill ( const smoab::Range entities,
VTKGridType *  grid 
) const [inline]

Definition at line 59 of file DataSetConverter.h.

References moab::Range::begin(), dim, moab::Range::end(), smoab::detail::ReadSparseTag::fill(), smoab::detail::LoadGeometry::fill(), smoab::Interface::findEntitiesWithDimension(), smoab::Interface::findHighestDimensionEntities(), smoab::getAllCells(), smoab::Tag::isComparable(), smoab::Tag::name(), readMaterialIds(), Tag, and smoab::Tag::value().

Referenced by vtkMoabReader::CreateSubBlocks().

    {
        //create a helper datastructure which can determines all the unique point ids
        //and converts moab connecitvity info to vtk connectivity

        //get all the cells for each parent entity and create
        // an entity set of those items
        int dim = this->Tag->value();
        typedef smoab::Range::const_iterator iterator;
        smoab::CellSets entitySets;
        for( iterator i = entities.begin(); i != entities.end(); ++i )
        {
            smoab::Range entitiesCells;
            if( this->Tag->isComparable() )
            {
                //if we are comparable only find the cells that match our tags dimension
                entitiesCells = this->Interface.findEntitiesWithDimension( *i, dim, true );
            }
            else
            {
                entitiesCells = this->Interface.findHighestDimensionEntities( *i, true );
            }
            smoab::CellSet set( *i, entitiesCells );
            entitySets.push_back( set );
        }

        moab::Range cells = smoab::getAllCells( entitySets );

        //convert the datastructure from a list of cells to a vtk data set
        detail::LoadGeometry loadGeom( cells, dim, this->Interface );
        loadGeom.fill( grid );

        if( this->readMaterialIds() )
        {

            detail::ReadSparseTag materialTagReading( entitySets, cells, this->Interface );

            smoab::MaterialTag mtag;
            vtkNew< vtkIntArray > materials;
            materials->SetName( mtag.name() );
            materialTagReading.fill( materials.GetPointer(), &mtag );
            grid->GetCellData()->AddArray( materials.GetPointer() );
        }

        //by default we always try to load the default tag
        detail::ReadSparseTag sTagReading( entitySets, cells, this->Interface );

        vtkNew< vtkIntArray > sparseTagData;
        sparseTagData->SetName( this->Tag->name() );
        sTagReading.fill( sparseTagData.GetPointer(), this->Tag );
        grid->GetCellData()->AddArray( sparseTagData.GetPointer() );

        return true;
    }
template<typename VTKGridType >
bool smoab::DataSetConverter::fill ( const smoab::EntityHandle entity,
VTKGridType *  grid,
const int  materialId = 0 
) const [inline]

Definition at line 119 of file DataSetConverter.h.

References dim, smoab::detail::ReadSparseTag::fill(), smoab::detail::LoadGeometry::fill(), smoab::Interface::findEntitiesWithDimension(), smoab::Interface::findHighestDimensionEntities(), smoab::Tag::isComparable(), smoab::detail::LoadGeometry::moabPoints(), smoab::Tag::name(), readMaterialIds(), readProperties(), Tag, and smoab::Tag::value().

    {
        //create a helper datastructure which can determines all the unique point ids
        //and converts moab connecitvity info to vtk connectivity

        smoab::Range cells;
        int dim = this->Tag->value();
        if( this->Tag->isComparable() )
        {
            //if we are comparable only find the cells that match our tags dimension
            cells = this->Interface.findEntitiesWithDimension( entity, dim, true );
        }
        else
        {
            //load subentities
            cells = this->Interface.findHighestDimensionEntities( entity, true );
        }

        //convert the datastructure from a list of cells to a vtk data set
        detail::LoadGeometry loadGeom( cells, dim, this->Interface );
        loadGeom.fill( grid );

        const smoab::Range& points = loadGeom.moabPoints();

        if( this->readProperties() )
        {
            this->readProperties( cells, grid->GetCellData() );
            this->readProperties( points, grid->GetPointData() );
        }

        smoab::CellSets cellSets;
        smoab::CellSet set( entity, cells );
        cellSets.push_back( set );
        if( this->readMaterialIds() )
        {
            smoab::MaterialTag mtag;
            detail::ReadSparseTag materialTagReading( cellSets, cells, this->Interface );

            vtkNew< vtkIntArray > materials;
            materials->SetName( mtag.name() );
            materialTagReading.fill( materials.GetPointer(), &mtag );
            grid->GetCellData()->AddArray( materials.GetPointer() );
        }

        //by default we always try to load the default tag
        detail::ReadSparseTag sTagReading( cellSets, cells, this->Interface );

        vtkNew< vtkIntArray > sparseTagData;
        sparseTagData->SetName( this->Tag->name() );
        sTagReading.fill( sparseTagData.GetPointer(), this->Tag );
        grid->GetCellData()->AddArray( sparseTagData.GetPointer() );

        return true;
    }
void smoab::DataSetConverter::readDenseTags ( std::vector< moab::Tag > &  tags,
smoab::Range const &  entities,
vtkFieldData *  field 
) const [inline, private]

Definition at line 193 of file DataSetConverter.h.

References MB_TAG_DENSE, MB_TYPE_DOUBLE, MB_TYPE_INTEGER, Moab, size, moab::Range::size(), moab::Interface::tag_get_data(), moab::Interface::tag_get_data_type(), moab::Interface::tag_get_length(), moab::Interface::tag_get_name(), moab::Interface::tag_get_type(), and TagType.

Referenced by readProperties().

    {
        typedef std::vector< moab::Tag >::const_iterator iterator;

        for( iterator i = tags.begin(); i != tags.end(); ++i )
        {
            moab::TagType tagType;
            moab::DataType tagDataType;

            this->Moab->tag_get_type( *i, tagType );
            this->Moab->tag_get_data_type( *i, tagDataType );

            //make sure it is only dense
            if( tagType != moab::MB_TAG_DENSE )
            {
                continue;
            }
            //and only integer and double
            if( tagDataType != moab::MB_TYPE_DOUBLE && tagDataType != moab::MB_TYPE_INTEGER )
            {
                //unsupported type, skip to next tag
                continue;
            }

            //read the name of the tag
            std::string name;
            name.reserve( 32 );
            this->Moab->tag_get_name( *i, name );

            //read the number of components of the tag
            int numComps = 1;

            this->Moab->tag_get_length( *i, numComps );

            //read the data if it is one of the two types we support
            int size = entities.size();
            if( tagDataType == moab::MB_TYPE_DOUBLE )
            {
                vtkNew< vtkDoubleArray > array;
                array->SetName( name.c_str() );
                array->SetNumberOfComponents( numComps );
                array->SetNumberOfTuples( size );

                //read directly into the double array
                this->Moab->tag_get_data( *i, entities, array->GetVoidPointer( 0 ) );
                field->AddArray( array.GetPointer() );
            }
            else if( tagDataType == moab::MB_TYPE_INTEGER )
            {
                vtkNew< vtkIntArray > array;
                array->SetName( name.c_str() );
                array->SetNumberOfComponents( numComps );
                array->SetNumberOfTuples( size );

                //read directly into the double array
                this->Moab->tag_get_data( *i, entities, array->GetVoidPointer( 0 ) );
                field->AddArray( array.GetPointer() );
            }
            else
            {
            }
        }
    }
void smoab::DataSetConverter::readMaterialIds ( bool  add) [inline]

Definition at line 35 of file DataSetConverter.h.

References ReadMaterialIds.

Referenced by vtkMoabReader::CreateSubBlocks().

    {
        this->ReadMaterialIds = add;
    }

Definition at line 39 of file DataSetConverter.h.

References ReadMaterialIds.

Referenced by fill().

    {
        return this->ReadMaterialIds;
    }
void smoab::DataSetConverter::readProperties ( bool  readProps) [inline]

Definition at line 44 of file DataSetConverter.h.

References ReadProperties.

Referenced by vtkMoabReader::CreateSubBlocks().

    {
        this->ReadProperties = readProps;
    }
bool smoab::DataSetConverter::readProperties ( ) const [inline]

Definition at line 48 of file DataSetConverter.h.

References ReadProperties.

Referenced by fill().

    {
        return this->ReadProperties;
    }
void smoab::DataSetConverter::readProperties ( smoab::Range const &  entities,
vtkFieldData *  field 
) const [inline, private]

Definition at line 176 of file DataSetConverter.h.

References moab::Range::empty(), moab::Range::front(), Moab, readDenseTags(), and moab::Interface::tag_get_tags_on_entity().

    {
        if( entities.empty() )
        {
            return;
        }

        //so we get all the tags and parse out the sparse and dense tags
        //that we support
        typedef std::vector< moab::Tag >::const_iterator iterator;
        std::vector< moab::Tag > tags;
        this->Moab->tag_get_tags_on_entity( entities.front(), tags );

        this->readDenseTags( tags, entities, field );
    }

Member Data Documentation

Definition at line 22 of file DataSetConverter.h.

Definition at line 27 of file DataSetConverter.h.

Definition at line 23 of file DataSetConverter.h.

Referenced by readDenseTags(), and readProperties().

Definition at line 25 of file DataSetConverter.h.

Referenced by readMaterialIds().

Definition at line 26 of file DataSetConverter.h.

Referenced by readProperties().

Definition at line 24 of file DataSetConverter.h.

Referenced by fill().

List of all members.


The documentation for this class was generated from the following file:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines