MOAB: Mesh Oriented datABase  (version 5.4.1)
ReadVtk.cpp
Go to the documentation of this file.
00001 /**
00002  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
00003  * storing and accessing finite element mesh data.
00004  *
00005  * Copyright 2004 Sandia Corporation.  Under the terms of Contract
00006  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
00007  * retains certain rights in this software.
00008  *
00009  * This library is free software; you can redistribute it and/or
00010  * modify it under the terms of the GNU Lesser General Public
00011  * License as published by the Free Software Foundation; either
00012  * version 2.1 of the License, or (at your option) any later version.
00013  *
00014  */
00015 
00016 /**
00017  * \class ReadVtk
00018  * \brief VTK reader from Mesquite
00019  * \author Jason Kraftcheck
00020  */
00021 
00022 #include <cassert>
00023 #include <cstdlib>
00024 #include <cstring>
00025 
00026 #include "ReadVtk.hpp"
00027 #include "moab/Range.hpp"
00028 #include "Internals.hpp"
00029 #include "moab/Interface.hpp"
00030 #include "moab/ReadUtilIface.hpp"
00031 #include "moab/FileOptions.hpp"
00032 #include "FileTokenizer.hpp"
00033 #include "moab/VtkUtil.hpp"
00034 #include "MBTagConventions.hpp"
00035 
00036 // #define MB_VTK_MATERIAL_SETS
00037 #ifdef MB_VTK_MATERIAL_SETS
00038 #include <map>
00039 #endif
00040 
00041 namespace moab
00042 {
00043 
00044 #ifdef MB_VTK_MATERIAL_SETS
00045 
00046 class Hash
00047 {
00048   public:
00049     unsigned long value;
00050 
00051     Hash()
00052     {
00053         this->value = 5381L;
00054     }
00055 
00056     Hash( const unsigned char* bytes, size_t len )
00057     {  // djb2, a hash by dan bernstein presented on comp.lang.c for hashing strings.
00058         this->value = 5381L;
00059         for( ; len; --len, ++bytes )
00060             this->value = this->value * 33 + ( *bytes );
00061     }
00062 
00063     Hash( bool duh )
00064     {
00065         this->value = duh;  // Hashing a single bit with a long is stupid but easy.
00066     }
00067 
00068     Hash( const Hash& src )
00069     {
00070         this->value = src.value;
00071     }
00072 
00073     Hash& operator=( const Hash& src )
00074     {
00075         this->value = src.value;
00076         return *this;
00077     }
00078 
00079     bool operator<( const Hash& other ) const
00080     {
00081         return this->value < other.value;
00082     }
00083 };
00084 
00085 // Pass this class opaque data + a handle and it adds the handle to a set
00086 // whose opaque data has the same hash. If no set exists for the hash,
00087 // it creates one. Each set is tagged with the opaque data.
00088 // When entities with different opaque data have the same hash, they
00089 // will be placed into the same set.
00090 // There will be no collisions when the opaque data is shorter than an
00091 // unsigned long, and this is really the only case we need to support.
00092 // The rest is bonus. Hash does quite well with strings, even those
00093 // with identical prefixes.
00094 class Modulator : public std::map< Hash, EntityHandle >
00095 {
00096   public:
00097     Modulator( Interface* iface ) : mesh( iface ) {}
00098 
00099     ErrorCode initialize( std::string tag_name, DataType mb_type, size_t sz, size_t per_elem )
00100     {
00101         std::vector< unsigned char > default_val;
00102         default_val.resize( sz * per_elem );
00103         ErrorCode rval = this->mesh->tag_get_handle( tag_name.c_str(), per_elem, mb_type, this->tag,
00104                                                      MB_TAG_SPARSE | MB_TAG_BYTES | MB_TAG_CREAT, &default_val[0] );
00105         return rval;
00106     }
00107 
00108     void add_entity( EntityHandle ent, const unsigned char* bytes, size_t len )
00109     {
00110         Hash h( bytes, len );
00111         EntityHandle mset = this->congruence_class( h, bytes );
00112         ErrorCode rval;
00113         rval = this->mesh->add_entities( mset, &ent, 1 );MB_CHK_SET_ERR_RET( rval, "Failed to add entities to mesh" );
00114     }
00115 
00116     void add_entities( Range& range, const unsigned char* bytes, size_t bytes_per_ent )
00117     {
00118         for( Range::iterator it = range.begin(); it != range.end(); ++it, bytes += bytes_per_ent )
00119         {
00120             Hash h( bytes, bytes_per_ent );
00121             EntityHandle mset = this->congruence_class( h, bytes );
00122             ErrorCode rval;
00123             rval = this->mesh->add_entities( mset, &*it, 1 );MB_CHK_SET_ERR_RET( rval, "Failed to add entities to mesh" );
00124         }
00125     }
00126 
00127     EntityHandle congruence_class( const Hash& h, const void* tag_data )
00128     {
00129         std::map< Hash, EntityHandle >::iterator it = this->find( h );
00130         if( it == this->end() )
00131         {
00132             EntityHandle mset;
00133             Range preexist;
00134             ErrorCode rval;
00135             rval = this->mesh->get_entities_by_type_and_tag( 0, MBENTITYSET, &this->tag, &tag_data, 1, preexist );MB_CHK_SET_ERR_RET_VAL( rval, "Failed to get entities by type and tag", (EntityHandle)0 );
00136             if( preexist.size() )
00137             {
00138                 mset = *preexist.begin();
00139             }
00140             else
00141             {
00142                 rval = this->mesh->create_meshset( MESHSET_SET, mset );MB_CHK_SET_ERR_RET_VAL( rval, "Failed to create mesh set", (EntityHandle)0 );
00143                 rval = this->mesh->tag_set_data( this->tag, &mset, 1, tag_data );MB_CHK_SET_ERR_RET_VAL( rval, "Failed to set tag data", (EntityHandle)0 );
00144             }
00145             ( *this )[h] = mset;
00146             return mset;
00147         }
00148         return it->second;
00149     }
00150 
00151     Interface* mesh;
00152     Tag tag;
00153 };
00154 #endif  // MB_VTK_MATERIAL_SETS
00155 
00156 ReaderIface* ReadVtk::factory( Interface* iface )
00157 {
00158     return new ReadVtk( iface );
00159 }
00160 
00161 ReadVtk::ReadVtk( Interface* impl ) : mdbImpl( impl ), mPartitionTagName( MATERIAL_SET_TAG_NAME )
00162 {
00163     mdbImpl->query_interface( readMeshIface );
00164 }
00165 
00166 ReadVtk::~ReadVtk()
00167 {
00168     if( readMeshIface )
00169     {
00170         mdbImpl->release_interface( readMeshIface );
00171         readMeshIface = 0;
00172     }
00173 }
00174 
00175 const char* const vtk_type_names[] = {
00176     "bit",  "char",          "unsigned_char", "short",  "unsigned_short", "int", "unsigned_int",
00177     "long", "unsigned_long", "float",         "double", "vtkIdType",      0 };
00178 
00179 ErrorCode ReadVtk::read_tag_values( const char* /* file_name */,
00180                                     const char* /* tag_name */,
00181                                     const FileOptions& /* opts */,
00182                                     std::vector< int >& /* tag_values_out */,
00183                                     const SubsetList* /* subset_list */ )
00184 {
00185     return MB_NOT_IMPLEMENTED;
00186 }
00187 
00188 ErrorCode ReadVtk::load_file( const char* filename,
00189                               const EntityHandle* /* file_set */,
00190                               const FileOptions& opts,
00191                               const ReaderIface::SubsetList* subset_list,
00192                               const Tag* file_id_tag )
00193 {
00194     ErrorCode result;
00195 
00196     int major, minor;
00197     char vendor_string[257];
00198     std::vector< Range > element_list;
00199     Range vertices;
00200 
00201     if( subset_list )
00202     {
00203         MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "Reading subset of files not supported for VTK" );
00204     }
00205 
00206     // Does the caller want a field to be used for partitioning the entities?
00207     // If not, we'll assume any scalar integer field named MATERIAL_SET specifies partitions.
00208     std::string partition_tag_name;
00209     result = opts.get_option( "PARTITION", partition_tag_name );
00210     if( result == MB_SUCCESS ) mPartitionTagName = partition_tag_name;
00211 
00212     FILE* file = fopen( filename, "r" );
00213     if( !file ) return MB_FILE_DOES_NOT_EXIST;
00214 
00215     // Read file header
00216 
00217     if( !fgets( vendor_string, sizeof( vendor_string ), file ) )
00218     {
00219         fclose( file );
00220         return MB_FAILURE;
00221     }
00222 
00223     if( !strchr( vendor_string, '\n' ) || 2 != sscanf( vendor_string, "# vtk DataFile Version %d.%d", &major, &minor ) )
00224     {
00225         fclose( file );
00226         return MB_FAILURE;
00227     }
00228 
00229     if( !fgets( vendor_string, sizeof( vendor_string ), file ) )
00230     {
00231         fclose( file );
00232         return MB_FAILURE;
00233     }
00234 
00235     // VTK spec says this should not exceed 256 chars.
00236     if( !strchr( vendor_string, '\n' ) )
00237     {
00238         fclose( file );
00239         MB_SET_ERR( MB_FAILURE, "Vendor string (line 2) exceeds 256 characters" );
00240     }
00241 
00242     // Check file type
00243 
00244     FileTokenizer tokens( file, readMeshIface );
00245     const char* const file_type_names[] = { "ASCII", "BINARY", 0 };
00246     int filetype                        = tokens.match_token( file_type_names );
00247     switch( filetype )
00248     {
00249         case 2:  // BINARY
00250             MB_SET_ERR( MB_FAILURE, "Cannot read BINARY VTK files" );
00251         default:  // ERROR
00252             return MB_FAILURE;
00253         case 1:  // ASCII
00254             break;
00255     }
00256 
00257     // Read the mesh
00258     if( !tokens.match_token( "DATASET" ) ) return MB_FAILURE;
00259     result = vtk_read_dataset( tokens, vertices, element_list );
00260     if( MB_SUCCESS != result ) return result;
00261 
00262     if( file_id_tag )
00263     {
00264         result = store_file_ids( *file_id_tag, vertices, element_list );
00265         if( MB_SUCCESS != result ) return result;
00266     }
00267 
00268     // Count the number of elements read
00269     long elem_count = 0;
00270     for( std::vector< Range >::iterator it = element_list.begin(); it != element_list.end(); ++it )
00271         elem_count += it->size();
00272 
00273     // Read attribute data until end of file.
00274     const char* const block_type_names[] = { "POINT_DATA", "CELL_DATA", 0 };
00275     std::vector< Range > vertex_list( 1 );
00276     vertex_list[0] = vertices;
00277     int blocktype  = 0;
00278     while( !tokens.eof() )
00279     {
00280         // Get POINT_DATA or CELL_DATA
00281         int new_block_type = tokens.match_token( block_type_names, false );
00282         if( tokens.eof() ) break;
00283 
00284         if( !new_block_type )
00285         {
00286             // If next token was neither POINT_DATA nor CELL_DATA,
00287             // then there's another attribute under the current one.
00288             if( blocktype )
00289                 tokens.unget_token();
00290             else
00291                 break;
00292         }
00293         else
00294         {
00295             blocktype = new_block_type;
00296             long count;
00297             if( !tokens.get_long_ints( 1, &count ) ) return MB_FAILURE;
00298 
00299             if( blocktype == 1 && (unsigned long)count != vertices.size() )
00300             {
00301                 MB_SET_ERR( MB_FAILURE, "Count inconsistent with number of vertices at line " << tokens.line_number() );
00302             }
00303             else if( blocktype == 2 && count != elem_count )
00304             {
00305                 MB_SET_ERR( MB_FAILURE, "Count inconsistent with number of elements at line " << tokens.line_number() );
00306             }
00307         }
00308 
00309         if( blocktype == 1 )
00310             result = vtk_read_attrib_data( tokens, vertex_list );
00311         else
00312             result = vtk_read_attrib_data( tokens, element_list );
00313 
00314         if( MB_SUCCESS != result ) return result;
00315     }
00316 
00317     return MB_SUCCESS;
00318 }
00319 
00320 ErrorCode ReadVtk::allocate_vertices( long num_verts,
00321                                       EntityHandle& start_handle_out,
00322                                       double*& x_coord_array_out,
00323                                       double*& y_coord_array_out,
00324                                       double*& z_coord_array_out )
00325 {
00326     ErrorCode result;
00327 
00328     // Create vertices
00329     std::vector< double* > arrays;
00330     start_handle_out = 0;
00331     result           = readMeshIface->get_node_coords( 3, num_verts, MB_START_ID, start_handle_out, arrays );
00332     if( MB_SUCCESS != result ) return result;
00333 
00334     x_coord_array_out = arrays[0];
00335     y_coord_array_out = arrays[1];
00336     z_coord_array_out = arrays[2];
00337 
00338     return MB_SUCCESS;
00339 }
00340 
00341 ErrorCode ReadVtk::read_vertices( FileTokenizer& tokens, long num_verts, EntityHandle& start_handle_out )
00342 {
00343     ErrorCode result;
00344     double *x, *y, *z;
00345 
00346     result = allocate_vertices( num_verts, start_handle_out, x, y, z );
00347     if( MB_SUCCESS != result ) return result;
00348 
00349     // Read vertex coordinates
00350     for( long vtx = 0; vtx < num_verts; ++vtx )
00351     {
00352         if( !tokens.get_doubles( 1, x++ ) || !tokens.get_doubles( 1, y++ ) || !tokens.get_doubles( 1, z++ ) )
00353             return MB_FAILURE;
00354     }
00355 
00356     return MB_SUCCESS;
00357 }
00358 
00359 ErrorCode ReadVtk::allocate_elements( long num_elements,
00360                                       int vert_per_element,
00361                                       EntityType type,
00362                                       EntityHandle& start_handle_out,
00363                                       EntityHandle*& conn_array_out,
00364                                       std::vector< Range >& append_to_this )
00365 {
00366     ErrorCode result;
00367 
00368     start_handle_out = 0;
00369     result = readMeshIface->get_element_connect( num_elements, vert_per_element, type, MB_START_ID, start_handle_out,
00370                                                  conn_array_out );
00371     if( MB_SUCCESS != result ) return result;
00372 
00373     Range range( start_handle_out, start_handle_out + num_elements - 1 );
00374     append_to_this.push_back( range );
00375     return MB_SUCCESS;
00376 }
00377 
00378 ErrorCode ReadVtk::vtk_read_dataset( FileTokenizer& tokens, Range& vertex_list, std::vector< Range >& element_list )
00379 {
00380     const char* const data_type_names[] = {
00381         "STRUCTURED_POINTS", "STRUCTURED_GRID", "UNSTRUCTURED_GRID", "POLYDATA", "RECTILINEAR_GRID", "FIELD", 0 };
00382     int datatype = tokens.match_token( data_type_names );
00383     switch( datatype )
00384     {
00385         case 1:
00386             return vtk_read_structured_points( tokens, vertex_list, element_list );
00387         case 2:
00388             return vtk_read_structured_grid( tokens, vertex_list, element_list );
00389         case 3:
00390             return vtk_read_unstructured_grid( tokens, vertex_list, element_list );
00391         case 4:
00392             return vtk_read_polydata( tokens, vertex_list, element_list );
00393         case 5:
00394             return vtk_read_rectilinear_grid( tokens, vertex_list, element_list );
00395         case 6:
00396             return vtk_read_field( tokens );
00397         default:
00398             return MB_FAILURE;
00399     }
00400 }
00401 
00402 ErrorCode ReadVtk::vtk_read_structured_points( FileTokenizer& tokens,
00403                                                Range& vertex_list,
00404                                                std::vector< Range >& elem_list )
00405 {
00406     long i, j, k;
00407     long dims[3];
00408     double origin[3], space[3];
00409     ErrorCode result;
00410 
00411     if( !tokens.match_token( "DIMENSIONS" ) || !tokens.get_long_ints( 3, dims ) || !tokens.get_newline() )
00412         return MB_FAILURE;
00413 
00414     if( dims[0] < 1 || dims[1] < 1 || dims[2] < 1 )
00415     {
00416         MB_SET_ERR( MB_FAILURE, "Invalid dimension at line " << tokens.line_number() );
00417     }
00418 
00419     if( !tokens.match_token( "ORIGIN" ) || !tokens.get_doubles( 3, origin ) || !tokens.get_newline() )
00420         return MB_FAILURE;
00421 
00422     const char* const spacing_names[] = { "SPACING", "ASPECT_RATIO", 0 };
00423     if( !tokens.match_token( spacing_names ) || !tokens.get_doubles( 3, space ) || !tokens.get_newline() )
00424         return MB_FAILURE;
00425 
00426     // Create vertices
00427     double *x, *y, *z;
00428     EntityHandle start_handle = 0;
00429     long num_verts            = dims[0] * dims[1] * dims[2];
00430     result                    = allocate_vertices( num_verts, start_handle, x, y, z );
00431     if( MB_SUCCESS != result ) return result;
00432     vertex_list.insert( start_handle, start_handle + num_verts - 1 );
00433 
00434     for( k = 0; k < dims[2]; ++k )
00435         for( j = 0; j < dims[1]; ++j )
00436             for( i = 0; i < dims[0]; ++i )
00437             {
00438                 *x = origin[0] + i * space[0];
00439                 ++x;
00440                 *y = origin[1] + j * space[1];
00441                 ++y;
00442                 *z = origin[2] + k * space[2];
00443                 ++z;
00444             }
00445 
00446     return vtk_create_structured_elems( dims, start_handle, elem_list );
00447 }
00448 
00449 ErrorCode ReadVtk::vtk_read_structured_grid( FileTokenizer& tokens,
00450                                              Range& vertex_list,
00451                                              std::vector< Range >& elem_list )
00452 {
00453     long num_verts, dims[3];
00454     ErrorCode result;
00455 
00456     if( !tokens.match_token( "DIMENSIONS" ) || !tokens.get_long_ints( 3, dims ) || !tokens.get_newline() )
00457         return MB_FAILURE;
00458 
00459     if( dims[0] < 1 || dims[1] < 1 || dims[2] < 1 )
00460     {
00461         MB_SET_ERR( MB_FAILURE, "Invalid dimension at line " << tokens.line_number() );
00462     }
00463 
00464     if( !tokens.match_token( "POINTS" ) || !tokens.get_long_ints( 1, &num_verts ) ||
00465         !tokens.match_token( vtk_type_names ) || !tokens.get_newline() )
00466         return MB_FAILURE;
00467 
00468     if( num_verts != ( dims[0] * dims[1] * dims[2] ) )
00469     {
00470         MB_SET_ERR( MB_FAILURE, "Point count not consistent with dimensions at line " << tokens.line_number() );
00471     }
00472 
00473     // Create and read vertices
00474     EntityHandle start_handle = 0;
00475     result                    = read_vertices( tokens, num_verts, start_handle );
00476     if( MB_SUCCESS != result ) return result;
00477     vertex_list.insert( start_handle, start_handle + num_verts - 1 );
00478 
00479     return vtk_create_structured_elems( dims, start_handle, elem_list );
00480 }
00481 
00482 ErrorCode ReadVtk::vtk_read_rectilinear_grid( FileTokenizer& tokens,
00483                                               Range& vertex_list,
00484                                               std::vector< Range >& elem_list )
00485 {
00486     int i, j, k;
00487     long dims[3];
00488     const char* labels[] = { "X_COORDINATES", "Y_COORDINATES", "Z_COORDINATES" };
00489     std::vector< double > coords[3];
00490     ErrorCode result;
00491 
00492     if( !tokens.match_token( "DIMENSIONS" ) || !tokens.get_long_ints( 3, dims ) || !tokens.get_newline() )
00493         return MB_FAILURE;
00494 
00495     if( dims[0] < 1 || dims[1] < 1 || dims[2] < 1 )
00496     {
00497         MB_SET_ERR( MB_FAILURE, "Invalid dimension at line " << tokens.line_number() );
00498     }
00499 
00500     for( i = 0; i < 3; i++ )
00501     {
00502         long count;
00503         if( !tokens.match_token( labels[i] ) || !tokens.get_long_ints( 1, &count ) ||
00504             !tokens.match_token( vtk_type_names ) )
00505             return MB_FAILURE;
00506 
00507         if( count != dims[i] )
00508         {
00509             MB_SET_ERR( MB_FAILURE, "Coordinate count inconsistent with dimensions at line " << tokens.line_number() );
00510         }
00511 
00512         coords[i].resize( count );
00513         if( !tokens.get_doubles( count, &coords[i][0] ) ) return MB_FAILURE;
00514     }
00515 
00516     // Create vertices
00517     double *x, *y, *z;
00518     EntityHandle start_handle = 0;
00519     long num_verts            = dims[0] * dims[1] * dims[2];
00520     result                    = allocate_vertices( num_verts, start_handle, x, y, z );
00521     if( MB_SUCCESS != result ) return result;
00522     vertex_list.insert( start_handle, start_handle + num_verts - 1 );
00523 
00524     // Calculate vertex coordinates
00525     for( k = 0; k < dims[2]; ++k )
00526         for( j = 0; j < dims[1]; ++j )
00527             for( i = 0; i < dims[0]; ++i )
00528             {
00529                 *x = coords[0][i];
00530                 ++x;
00531                 *y = coords[1][j];
00532                 ++y;
00533                 *z = coords[2][k];
00534                 ++z;
00535             }
00536 
00537     return vtk_create_structured_elems( dims, start_handle, elem_list );
00538 }
00539 
00540 ErrorCode ReadVtk::vtk_read_polydata( FileTokenizer& tokens, Range& vertex_list, std::vector< Range >& elem_list )
00541 {
00542     ErrorCode result;
00543     long num_verts;
00544     const char* const poly_data_names[] = { "VERTICES", "LINES", "POLYGONS", "TRIANGLE_STRIPS", 0 };
00545 
00546     if( !tokens.match_token( "POINTS" ) || !tokens.get_long_ints( 1, &num_verts ) ||
00547         !tokens.match_token( vtk_type_names ) || !tokens.get_newline() )
00548         return MB_FAILURE;
00549 
00550     if( num_verts < 1 )
00551     {
00552         MB_SET_ERR( MB_FAILURE, "Invalid point count at line " << tokens.line_number() );
00553     }
00554 
00555     // Create vertices and read coordinates
00556     EntityHandle start_handle = 0;
00557     result                    = read_vertices( tokens, num_verts, start_handle );
00558     if( MB_SUCCESS != result ) return result;
00559     vertex_list.insert( start_handle, start_handle + num_verts - 1 );
00560 
00561     int poly_type = tokens.match_token( poly_data_names );
00562     switch( poly_type )
00563     {
00564         case 0:
00565             result = MB_FAILURE;
00566             break;
00567         case 1:
00568             MB_SET_ERR( MB_FAILURE, "Vertex element type at line " << tokens.line_number() );
00569             break;
00570         case 2:
00571             MB_SET_ERR( MB_FAILURE, "Unsupported type: polylines at line " << tokens.line_number() );
00572             result = MB_FAILURE;
00573             break;
00574         case 3:
00575             result = vtk_read_polygons( tokens, start_handle, elem_list );
00576             break;
00577         case 4:
00578             MB_SET_ERR( MB_FAILURE, "Unsupported type: triangle strips at line " << tokens.line_number() );
00579             result = MB_FAILURE;
00580             break;
00581     }
00582 
00583     return result;
00584 }
00585 
00586 ErrorCode ReadVtk::vtk_read_polygons( FileTokenizer& tokens, EntityHandle first_vtx, std::vector< Range >& elem_list )
00587 {
00588     ErrorCode result;
00589     long size[2];
00590 
00591     if( !tokens.get_long_ints( 2, size ) || !tokens.get_newline() ) return MB_FAILURE;
00592 
00593     const Range empty;
00594     std::vector< EntityHandle > conn_hdl;
00595     std::vector< long > conn_idx;
00596     EntityHandle first = 0, prev = 0, handle;
00597     for( int i = 0; i < size[0]; ++i )
00598     {
00599         long count;
00600         if( !tokens.get_long_ints( 1, &count ) ) return MB_FAILURE;
00601         conn_idx.resize( count );
00602         conn_hdl.resize( count );
00603         if( !tokens.get_long_ints( count, &conn_idx[0] ) ) return MB_FAILURE;
00604 
00605         for( long j = 0; j < count; ++j )
00606             conn_hdl[j] = first_vtx + conn_idx[j];
00607 
00608         result = mdbImpl->create_element( MBPOLYGON, &conn_hdl[0], count, handle );
00609         if( MB_SUCCESS != result ) return result;
00610 
00611         if( prev + 1 != handle )
00612         {
00613             if( first )
00614             {  // True except for first iteration (first == 0)
00615                 if( elem_list.empty() || first < elem_list.back().front() )  // Only need new range if order would get
00616                                                                              // mixed up, or we just began inserting
00617                     elem_list.push_back( empty );
00618                 elem_list.back().insert( first, prev );
00619             }
00620             first = handle;
00621         }
00622         prev = handle;
00623     }
00624     if( first )
00625     {                                                                // True unless no elements (size[0] == 0)
00626         if( elem_list.empty() || first < elem_list.back().front() )  // Only need new range if order would get mixed
00627                                                                      // up, or we just began inserting
00628             elem_list.push_back( empty );
00629         elem_list.back().insert( first, prev );
00630     }
00631 
00632     return MB_SUCCESS;
00633 }
00634 
00635 ErrorCode ReadVtk::vtk_read_unstructured_grid( FileTokenizer& tokens,
00636                                                Range& vertex_list,
00637                                                std::vector< Range >& elem_list )
00638 {
00639     ErrorCode result;
00640     long i, num_verts, num_elems[2];
00641     EntityHandle tmp_conn_list[27];
00642 
00643     // Poorly formatted VTK legacy format document seems to
00644     // lead many to think that a FIELD block can occur within
00645     // an UNSTRUCTURED_GRID dataset rather than as its own data
00646     // set. So allow for field data between other blocks of
00647     // data.
00648 
00649     const char* pts_str[] = { "FIELD", "POINTS", 0 };
00650     while( 1 == ( i = tokens.match_token( pts_str ) ) )
00651     {
00652         result = vtk_read_field( tokens );
00653         if( MB_SUCCESS != result ) return result;
00654     }
00655     if( i != 2 ) return MB_FAILURE;
00656 
00657     if( !tokens.get_long_ints( 1, &num_verts ) || !tokens.match_token( vtk_type_names ) || !tokens.get_newline() )
00658         return MB_FAILURE;
00659 
00660     if( num_verts < 1 )
00661     {
00662         MB_SET_ERR( MB_FAILURE, "Invalid point count at line " << tokens.line_number() );
00663     }
00664 
00665     // Create vertices and read coordinates
00666     EntityHandle first_vertex = 0;
00667     result                    = read_vertices( tokens, num_verts, first_vertex );
00668     if( MB_SUCCESS != result ) return result;
00669     vertex_list.insert( first_vertex, first_vertex + num_verts - 1 );
00670 
00671     const char* cell_str[] = { "FIELD", "CELLS", 0 };
00672     while( 1 == ( i = tokens.match_token( cell_str ) ) )
00673     {
00674         result = vtk_read_field( tokens );
00675         if( MB_SUCCESS != result ) return result;
00676     }
00677     if( i != 2 ) return MB_FAILURE;
00678 
00679     if( !tokens.get_long_ints( 2, num_elems ) || !tokens.get_newline() ) return MB_FAILURE;
00680 
00681     // Read element connectivity for all elements
00682     std::vector< long > connectivity( num_elems[1] );
00683     if( !tokens.get_long_ints( num_elems[1], &connectivity[0] ) ) return MB_FAILURE;
00684 
00685     if( !tokens.match_token( "CELL_TYPES" ) || !tokens.get_long_ints( 1, &num_elems[1] ) || !tokens.get_newline() )
00686         return MB_FAILURE;
00687 
00688     // Read element types
00689     std::vector< long > types( num_elems[0] );
00690     if( !tokens.get_long_ints( num_elems[0], &types[0] ) ) return MB_FAILURE;
00691 
00692     // Create elements in blocks of the same type
00693     // It is important to preserve the order in
00694     // which the elements were read for later reading
00695     // attribute data.
00696     long id                                 = 0;
00697     std::vector< long >::iterator conn_iter = connectivity.begin();
00698     while( id < num_elems[0] )
00699     {
00700         unsigned vtk_type = types[id];
00701         if( vtk_type >= VtkUtil::numVtkElemType ) return MB_FAILURE;
00702 
00703         EntityType type = VtkUtil::vtkElemTypes[vtk_type].mb_type;
00704         if( type == MBMAXTYPE )
00705         {
00706             MB_SET_ERR( MB_FAILURE, "Unsupported VTK element type: " << VtkUtil::vtkElemTypes[vtk_type].name << " ("
00707                                                                      << vtk_type << ")" );
00708         }
00709 
00710         int num_vtx = *conn_iter;
00711         if( type != MBPOLYGON && type != MBPOLYHEDRON && num_vtx != (int)VtkUtil::vtkElemTypes[vtk_type].num_nodes )
00712         {
00713             MB_SET_ERR( MB_FAILURE, "Cell " << id << " is of type '" << VtkUtil::vtkElemTypes[vtk_type].name
00714                                             << "' but has " << num_vtx << " vertices" );
00715         }
00716 
00717         // Find any subsequent elements of the same type
00718         // if polyhedra, need to look at the number of faces to put in the same range
00719         std::vector< long >::iterator conn_iter2 = conn_iter + num_vtx + 1;
00720         long end_id                              = id + 1;
00721         if( MBPOLYHEDRON != type )
00722         {
00723             while( end_id < num_elems[0] && (unsigned)types[end_id] == vtk_type && *conn_iter2 == num_vtx )
00724             {
00725                 ++end_id;
00726                 conn_iter2 += num_vtx + 1;
00727             }
00728         }
00729         else
00730         {
00731             // advance only if next is polyhedron too, and if number of faces is the same
00732             int num_faces = conn_iter[1];
00733             while( end_id < num_elems[0] && (unsigned)types[end_id] == vtk_type && conn_iter2[1] == num_faces )
00734             {
00735                 ++end_id;
00736                 conn_iter2 += num_vtx + 1;
00737             }
00738             // num_vtx becomes in this case num_faces
00739             num_vtx = num_faces;  // for polyhedra, this is what we want
00740             // trigger vertex adjacency call
00741             Range firstFaces;
00742             mdbImpl->get_adjacencies( &first_vertex, 1, 2, false, firstFaces );
00743         }
00744 
00745         // Allocate element block
00746         long num_elem             = end_id - id;
00747         EntityHandle start_handle = 0;
00748         EntityHandle* conn_array;
00749 
00750         // if type is MBVERTEX, skip, we will not create elements with one vertex
00751         if( MBVERTEX == type )
00752         {
00753             id += num_elem;
00754             conn_iter += 2 * num_elem;  // skip 2 * number of 1-vertex elements
00755             continue;
00756         }
00757 
00758         result = allocate_elements( num_elem, num_vtx, type, start_handle, conn_array, elem_list );
00759         if( MB_SUCCESS != result ) return result;
00760 
00761         EntityHandle* conn_sav = conn_array;
00762 
00763         // Store element connectivity
00764         if( type != MBPOLYHEDRON )
00765         {
00766             for( ; id < end_id; ++id )
00767             {
00768                 if( conn_iter == connectivity.end() )
00769                 {
00770                     MB_SET_ERR( MB_FAILURE, "Connectivity data truncated at cell " << id );
00771                 }
00772                 // Make sure connectivity length is correct.
00773                 if( *conn_iter != num_vtx )
00774                 {
00775                     MB_SET_ERR( MB_FAILURE, "Cell " << id << " is of type '" << VtkUtil::vtkElemTypes[vtk_type].name
00776                                                     << "' but has " << num_vtx << " vertices" );
00777                 }
00778                 ++conn_iter;
00779 
00780                 for( i = 0; i < num_vtx; ++i, ++conn_iter )
00781                 {
00782                     if( conn_iter == connectivity.end() )
00783                     {
00784                         MB_SET_ERR( MB_FAILURE, "Connectivity data truncated at cell " << id );
00785                     }
00786 
00787                     conn_array[i] = *conn_iter + first_vertex;
00788                 }
00789 
00790                 const unsigned* order = VtkUtil::vtkElemTypes[vtk_type].node_order;
00791                 if( order )
00792                 {
00793                     assert( num_vtx * sizeof( EntityHandle ) <= sizeof( tmp_conn_list ) );
00794                     memcpy( tmp_conn_list, conn_array, num_vtx * sizeof( EntityHandle ) );
00795                     for( int j = 0; j < num_vtx; ++j )
00796                         conn_array[order[j]] = tmp_conn_list[j];
00797                 }
00798 
00799                 conn_array += num_vtx;
00800             }
00801         }
00802         else  // type == MBPOLYHEDRON
00803         {
00804             // in some cases, we may need to create new elements; will it screw the tags?
00805             // not if the file was not from moab
00806             ErrorCode rv = MB_SUCCESS;
00807             for( ; id < end_id; ++id )
00808             {
00809                 if( conn_iter == connectivity.end() )
00810                 {
00811                     MB_SET_ERR( MB_FAILURE, "Connectivity data truncated at polyhedra cell " << id );
00812                 }
00813                 ++conn_iter;
00814                 // iterator is now at number of faces
00815                 // we should check it is indeed num_vtx
00816                 int num_faces = *conn_iter;
00817                 if( num_faces != num_vtx ) MB_SET_ERR( MB_FAILURE, "Connectivity data wrong at polyhedra cell " << id );
00818 
00819                 EntityHandle connec[20];  // we bet we will have only 20 vertices at most, in a
00820                                           // face in a polyhedra
00821                 for( int j = 0; j < num_faces; j++ )
00822                 {
00823                     conn_iter++;
00824                     int numverticesInFace = (int)*conn_iter;
00825                     if( numverticesInFace > 20 )
00826                         MB_SET_ERR( MB_FAILURE,
00827                                     "too many vertices in face index " << j << " for polyhedra cell " << id );
00828                     // need to find the face, but first fill with vertices
00829                     for( int k = 0; k < numverticesInFace; k++ )
00830                     {
00831                         connec[k] = first_vertex + *( ++conn_iter );  //
00832                     }
00833                     Range adjFaces;
00834                     // find a face with these vertices; if not, we need to create one, on the fly :(
00835                     rv = mdbImpl->get_adjacencies( connec, numverticesInFace, 2, false, adjFaces );MB_CHK_ERR( rv );
00836                     if( adjFaces.size() >= 1 )
00837                     {
00838                         conn_array[j] = adjFaces[0];  // get the first face found
00839                     }
00840                     else
00841                     {
00842                         // create the face; tri, quad or polygon
00843                         EntityType etype = MBTRI;
00844                         if( 4 == numverticesInFace ) etype = MBQUAD;
00845                         if( 4 < numverticesInFace ) etype = MBPOLYGON;
00846 
00847                         rv = mdbImpl->create_element( etype, connec, numverticesInFace, conn_array[j] );MB_CHK_ERR( rv );
00848                     }
00849                 }
00850 
00851                 conn_array += num_vtx;  // advance for next polyhedra
00852                 conn_iter++;            // advance to the next field
00853             }
00854         }
00855 
00856         // Notify MOAB of the new elements
00857         result = readMeshIface->update_adjacencies( start_handle, num_elem, num_vtx, conn_sav );
00858         if( MB_SUCCESS != result ) return result;
00859     }
00860 
00861     return MB_SUCCESS;
00862 }
00863 
00864 ErrorCode ReadVtk::vtk_create_structured_elems( const long* dims,
00865                                                 EntityHandle first_vtx,
00866                                                 std::vector< Range >& elem_list )
00867 {
00868     ErrorCode result;
00869     // int non_zero[3] = {0, 0, 0}; // True if dim > 0 for x, y, z respectively
00870     long elem_dim  = 0;           // Element dimension (2->quad, 3->hex)
00871     long num_elems = 1;           // Total number of elements
00872     long vert_per_elem;           // Element connectivity length
00873     long edims[3] = { 1, 1, 1 };  // Number of elements in each grid direction
00874 
00875     // Populate above data
00876     for( int d = 0; d < 3; d++ )
00877     {
00878         if( dims[d] > 1 )
00879         {
00880             // non_zero[elem_dim] = d;
00881             ++elem_dim;
00882             edims[d] = dims[d] - 1;
00883             num_elems *= edims[d];
00884         }
00885     }
00886     vert_per_elem = 1 << elem_dim;
00887 
00888     // Get element type from element dimension
00889     EntityType type;
00890     switch( elem_dim )
00891     {
00892         case 1:
00893             type = MBEDGE;
00894             break;
00895         case 2:
00896             type = MBQUAD;
00897             break;
00898         case 3:
00899             type = MBHEX;
00900             break;
00901         default:
00902             MB_SET_ERR( MB_FAILURE, "Invalid dimension for structured elements: " << elem_dim );
00903     }
00904 
00905     // Allocate storage for elements
00906     EntityHandle start_handle = 0;
00907     EntityHandle* conn_array;
00908     result = allocate_elements( num_elems, vert_per_elem, type, start_handle, conn_array, elem_list );
00909     if( MB_SUCCESS != result ) return MB_FAILURE;
00910 
00911     EntityHandle* conn_sav = conn_array;
00912 
00913     // Offsets of element vertices in grid relative to corner closest to origin
00914     long k                = dims[0] * dims[1];
00915     const long corners[8] = { 0, 1, 1 + dims[0], dims[0], k, k + 1, k + 1 + dims[0], k + dims[0] };
00916 
00917     // Populate element list
00918     for( long z = 0; z < edims[2]; ++z )
00919         for( long y = 0; y < edims[1]; ++y )
00920             for( long x = 0; x < edims[0]; ++x )
00921             {
00922                 const long index = x + y * dims[0] + z * ( dims[0] * dims[1] );
00923                 for( long j = 0; j < vert_per_elem; ++j, ++conn_array )
00924                     *conn_array = index + corners[j] + first_vtx;
00925             }
00926 
00927     // Notify MOAB of the new elements
00928     result = readMeshIface->update_adjacencies( start_handle, num_elems, vert_per_elem, conn_sav );
00929     if( MB_SUCCESS != result ) return result;
00930 
00931     return MB_SUCCESS;
00932 }
00933 
00934 ErrorCode ReadVtk::vtk_read_field( FileTokenizer& tokens )
00935 {
00936     // This is not supported yet.
00937     // Parse the data but throw it away because
00938     // Mesquite has no internal representation for it.
00939 
00940     // Could save this in tags, but the only useful thing that
00941     // could be done with the data is to write it back out
00942     // with the modified mesh. As there's no way to save the
00943     // type of a tag in Mesquite, it cannot be written back
00944     // out correctly either.
00945     // FIXME: Don't know what to do with this data.
00946     // For now, read it and throw it out.
00947 
00948     long num_arrays;
00949     if( !tokens.get_string() ||  // Name
00950         !tokens.get_long_ints( 1, &num_arrays ) )
00951         return MB_FAILURE;
00952 
00953     for( long i = 0; i < num_arrays; ++i )
00954     {
00955         /*const char* name =*/tokens.get_string();
00956 
00957         long dims[2];
00958         if( !tokens.get_long_ints( 2, dims ) || !tokens.match_token( vtk_type_names ) ) return MB_FAILURE;
00959 
00960         long num_vals = dims[0] * dims[1];
00961 
00962         for( long j = 0; j < num_vals; j++ )
00963         {
00964             double junk;
00965             if( !tokens.get_doubles( 1, &junk ) ) return MB_FAILURE;
00966         }
00967     }
00968 
00969     return MB_SUCCESS;
00970 }
00971 
00972 ErrorCode ReadVtk::vtk_read_attrib_data( FileTokenizer& tokens, std::vector< Range >& entities )
00973 {
00974     const char* const type_names[] = { "SCALARS", "COLOR_SCALARS", "VECTORS", "NORMALS", "TEXTURE_COORDINATES",
00975                                        "TENSORS", "FIELD",         0 };
00976 
00977     int type             = tokens.match_token( type_names );
00978     const char* tmp_name = tokens.get_string();
00979     if( !type || !tmp_name ) return MB_FAILURE;
00980 
00981     std::string name_alloc( tmp_name );
00982     const char* name = name_alloc.c_str();
00983     switch( type )
00984     {
00985         case 1:
00986             return vtk_read_scalar_attrib( tokens, entities, name );
00987         case 2:
00988             return vtk_read_color_attrib( tokens, entities, name );
00989         case 3:
00990             return vtk_read_vector_attrib( tokens, entities, name );
00991         case 4:
00992             return vtk_read_vector_attrib( tokens, entities, name );
00993         case 5:
00994             return vtk_read_texture_attrib( tokens, entities, name );
00995         case 6:
00996             return vtk_read_tensor_attrib( tokens, entities, name );
00997         case 7:
00998             return vtk_read_field_attrib( tokens, entities, name );
00999     }
01000 
01001     return MB_FAILURE;
01002 }
01003 
01004 ErrorCode ReadVtk::vtk_read_tag_data( FileTokenizer& tokens,
01005                                       int type,
01006                                       size_t per_elem,
01007                                       std::vector< Range >& entities,
01008                                       const char* name )
01009 {
01010     ErrorCode result;
01011     DataType mb_type;
01012     if( type == 1 )
01013     {
01014         mb_type = MB_TYPE_BIT;
01015     }
01016     else if( type >= 2 && type <= 9 )
01017     {
01018         mb_type = MB_TYPE_INTEGER;
01019     }
01020     else if( type == 10 || type == 11 )
01021     {
01022         mb_type = MB_TYPE_DOUBLE;
01023     }
01024     else if( type == 12 )
01025     {
01026         mb_type = MB_TYPE_INTEGER;
01027     }
01028     else
01029         return MB_FAILURE;
01030 
01031 #ifdef MB_VTK_MATERIAL_SETS
01032     size_t size;
01033     if( type == 1 )
01034     {
01035         size = sizeof( bool );
01036     }
01037     else if( type >= 2 && type <= 9 )
01038     {
01039         size = sizeof( int );
01040     }
01041     else if( type == 10 || type == 11 )
01042     {
01043         size = sizeof( double );
01044     }
01045     else /* (type == 12) */
01046     {
01047         size = 4;  // Could be 4 or 8, but we don't know. Hope it's 4 because MOAB doesn't support
01048                    // 64-bit ints.
01049     }
01050     Modulator materialMap( this->mdbImpl );
01051     result = materialMap.initialize( this->mPartitionTagName, mb_type, size, per_elem );MB_CHK_SET_ERR( result, "MaterialMap tag (" << this->mPartitionTagName << ") creation failed." );
01052     bool isMaterial = size * per_elem <= 4 &&              // Must have int-sized values (ParallelComm requires it)
01053                       !this->mPartitionTagName.empty() &&  // Must have a non-empty field name...
01054                       !strcmp( name, this->mPartitionTagName.c_str() );  // ... that matches our spec.
01055 #endif                                                                   // MB_VTK_MATERIAL_SETS
01056 
01057     // Get/create tag
01058     Tag handle;
01059     result = mdbImpl->tag_get_handle( name, per_elem, mb_type, handle, MB_TAG_DENSE | MB_TAG_CREAT );MB_CHK_SET_ERR( result, "Tag name conflict for attribute \"" << name << "\" at line " << tokens.line_number() );
01060 
01061     std::vector< Range >::iterator iter;
01062 
01063     if( type == 1 )
01064     {
01065         for( iter = entities.begin(); iter != entities.end(); ++iter )
01066         {
01067             bool* data = new bool[iter->size() * per_elem];
01068             if( !tokens.get_booleans( per_elem * iter->size(), data ) )
01069             {
01070                 delete[] data;
01071                 return MB_FAILURE;
01072             }
01073 
01074             bool* data_iter          = data;
01075             Range::iterator ent_iter = iter->begin();
01076             for( ; ent_iter != iter->end(); ++ent_iter )
01077             {
01078                 unsigned char bits = 0;
01079                 for( unsigned j = 0; j < per_elem; ++j, ++data_iter )
01080                     bits |= (unsigned char)( *data_iter << j );
01081 #ifdef MB_VTK_MATERIAL_SETS
01082                 if( isMaterial ) materialMap.add_entity( *ent_iter, &bits, 1 );
01083 #endif  // MB_VTK_MATERIAL_SETS
01084                 result = mdbImpl->tag_set_data( handle, &*ent_iter, 1, &bits );
01085                 if( MB_SUCCESS != result )
01086                 {
01087                     delete[] data;
01088                     return result;
01089                 }
01090             }
01091             delete[] data;
01092         }
01093     }
01094     else if( ( type >= 2 && type <= 9 ) || type == 12 )
01095     {
01096         std::vector< int > data;
01097         for( iter = entities.begin(); iter != entities.end(); ++iter )
01098         {
01099             data.resize( iter->size() * per_elem );
01100             if( !tokens.get_integers( iter->size() * per_elem, &data[0] ) ) return MB_FAILURE;
01101 #ifdef MB_VTK_MATERIAL_SETS
01102             if( isMaterial ) materialMap.add_entities( *iter, (unsigned char*)&data[0], per_elem * size );
01103 #endif  // MB_VTK_MATERIAL_SETS
01104             result = mdbImpl->tag_set_data( handle, *iter, &data[0] );
01105             if( MB_SUCCESS != result ) return result;
01106         }
01107     }
01108     else if( type == 10 || type == 11 )
01109     {
01110         std::vector< double > data;
01111         for( iter = entities.begin(); iter != entities.end(); ++iter )
01112         {
01113             data.resize( iter->size() * per_elem );
01114             if( !tokens.get_doubles( iter->size() * per_elem, &data[0] ) ) return MB_FAILURE;
01115 #ifdef MB_VTK_MATERIAL_SETS
01116             if( isMaterial ) materialMap.add_entities( *iter, (unsigned char*)&data[0], per_elem * size );
01117 #endif  // MB_VTK_MATERIAL_SETS
01118             result = mdbImpl->tag_set_data( handle, *iter, &data[0] );
01119             if( MB_SUCCESS != result ) return result;
01120         }
01121     }
01122 
01123     return MB_SUCCESS;
01124 }
01125 
01126 ErrorCode ReadVtk::vtk_read_scalar_attrib( FileTokenizer& tokens, std::vector< Range >& entities, const char* name )
01127 {
01128     int type = tokens.match_token( vtk_type_names );
01129     if( !type ) return MB_FAILURE;
01130 
01131     long size;
01132     const char* tok = tokens.get_string();
01133     if( !tok ) return MB_FAILURE;
01134 
01135     const char* end = 0;
01136     size            = strtol( tok, (char**)&end, 0 );
01137     if( *end )
01138     {
01139         size = 1;
01140         tokens.unget_token();
01141     }
01142 
01143     // VTK spec says cannot be greater than 4--do we care?
01144     if( size < 1 )
01145     {  //|| size > 4)
01146         MB_SET_ERR( MB_FAILURE, "Scalar count out of range [1,4] at line " << tokens.line_number() );
01147     }
01148 
01149     if( !tokens.match_token( "LOOKUP_TABLE" ) || !tokens.match_token( "default" ) ) return MB_FAILURE;
01150 
01151     return vtk_read_tag_data( tokens, type, size, entities, name );
01152 }
01153 
01154 ErrorCode ReadVtk::vtk_read_color_attrib( FileTokenizer& tokens, std::vector< Range >& entities, const char* name )
01155 {
01156     long size;
01157     if( !tokens.get_long_ints( 1, &size ) || size < 1 ) return MB_FAILURE;
01158 
01159     return vtk_read_tag_data( tokens, 10, size, entities, name );
01160 }
01161 
01162 ErrorCode ReadVtk::vtk_read_vector_attrib( FileTokenizer& tokens, std::vector< Range >& entities, const char* name )
01163 {
01164     int type = tokens.match_token( vtk_type_names );
01165     if( !type ) return MB_FAILURE;
01166 
01167     return vtk_read_tag_data( tokens, type, 3, entities, name );
01168 }
01169 
01170 ErrorCode ReadVtk::vtk_read_texture_attrib( FileTokenizer& tokens, std::vector< Range >& entities, const char* name )
01171 {
01172     int type, dim;
01173     if( !tokens.get_integers( 1, &dim ) || !( type = tokens.match_token( vtk_type_names ) ) ) return MB_FAILURE;
01174 
01175     if( dim < 1 || dim > 3 )
01176     {
01177         MB_SET_ERR( MB_FAILURE, "Invalid dimension (" << dim << ") at line " << tokens.line_number() );
01178     }
01179 
01180     return vtk_read_tag_data( tokens, type, dim, entities, name );
01181 }
01182 
01183 ErrorCode ReadVtk::vtk_read_tensor_attrib( FileTokenizer& tokens, std::vector< Range >& entities, const char* name )
01184 {
01185     int type = tokens.match_token( vtk_type_names );
01186     if( !type ) return MB_FAILURE;
01187 
01188     return vtk_read_tag_data( tokens, type, 9, entities, name );
01189 }
01190 
01191 ErrorCode ReadVtk::vtk_read_field_attrib( FileTokenizer& tokens, std::vector< Range >& entities, const char* )
01192 {
01193     long num_fields;
01194     if( !tokens.get_long_ints( 1, &num_fields ) ) return MB_FAILURE;
01195 
01196     long i;
01197     for( i = 0; i < num_fields; ++i )
01198     {
01199         const char* tok = tokens.get_string();
01200         if( !tok ) return MB_FAILURE;
01201 
01202         std::string name_alloc( tok );
01203 
01204         long num_comp;
01205         if( !tokens.get_long_ints( 1, &num_comp ) ) return MB_FAILURE;
01206 
01207         long num_tuples;
01208         if( !tokens.get_long_ints( 1, &num_tuples ) ) return MB_FAILURE;
01209 
01210         int type = tokens.match_token( vtk_type_names );
01211         if( !type ) return MB_FAILURE;
01212 
01213         ErrorCode result = vtk_read_tag_data( tokens, type, num_comp, entities, name_alloc.c_str() );MB_CHK_SET_ERR( result, "Error reading data for field \"" << name_alloc << "\" (" << num_comp << " components, "
01214                                                                   << num_tuples << " tuples, type " << type
01215                                                                   << ") at line " << tokens.line_number() );
01216     }
01217 
01218     return MB_SUCCESS;
01219 }
01220 
01221 ErrorCode ReadVtk::store_file_ids( Tag tag, const Range& verts, const std::vector< Range >& elems )
01222 {
01223     ErrorCode rval;
01224 
01225     rval = readMeshIface->assign_ids( tag, verts );
01226     if( MB_SUCCESS != rval ) return rval;
01227 
01228     int id = 0;
01229     for( size_t i = 0; i < elems.size(); ++i )
01230     {
01231         rval = readMeshIface->assign_ids( tag, elems[i], id );
01232         id += elems[i].size();
01233     }
01234 
01235     return MB_SUCCESS;
01236 }
01237 
01238 }  // namespace moab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines