MOAB: Mesh Oriented datABase
(version 5.2.1)
|
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 #ifdef WIN32 00017 #ifdef _DEBUG 00018 // turn off warnings that say they debugging identifier has been truncated 00019 // this warning comes up when using some STL containers 00020 #pragma warning( disable : 4786 ) 00021 #endif 00022 #endif 00023 00024 #include "WriteVtk.hpp" 00025 #include "moab/VtkUtil.hpp" 00026 #include "SysUtil.hpp" 00027 00028 #include <fstream> 00029 #include <iostream> 00030 #include <cstdio> 00031 #include <cassert> 00032 #include <vector> 00033 #include <set> 00034 #include <map> 00035 #include <iterator> 00036 00037 #include "moab/Interface.hpp" 00038 #include "moab/Range.hpp" 00039 #include "moab/CN.hpp" 00040 #include "MBTagConventions.hpp" 00041 #include "moab/WriteUtilIface.hpp" 00042 #include "Internals.hpp" 00043 #include "moab/FileOptions.hpp" 00044 00045 #define INS_ID( stringvar, prefix, id ) sprintf( stringvar, prefix, id ) 00046 00047 namespace moab 00048 { 00049 00050 const int DEFAULT_PRECISION = 10; 00051 const bool DEFAULT_STRICT = true; 00052 00053 WriterIface* WriteVtk::factory( Interface* iface ) 00054 { 00055 return new WriteVtk( iface ); 00056 } 00057 00058 WriteVtk::WriteVtk( Interface* impl ) 00059 : mbImpl( impl ), writeTool( 0 ), mStrict( DEFAULT_STRICT ), freeNodes( 0 ), createOneNodeCells( false ) 00060 { 00061 assert( impl != NULL ); 00062 impl->query_interface( writeTool ); 00063 } 00064 00065 WriteVtk::~WriteVtk() 00066 { 00067 mbImpl->release_interface( writeTool ); 00068 } 00069 00070 ErrorCode WriteVtk::write_file( const char* file_name, const bool overwrite, const FileOptions& opts, 00071 const EntityHandle* output_list, const int num_sets, 00072 const std::vector< std::string >& /* qa_list */, const Tag* tag_list, int num_tags, 00073 int /* export_dimension */ ) 00074 { 00075 ErrorCode rval; 00076 00077 // Get precision for node coordinates 00078 int precision; 00079 if( MB_SUCCESS != opts.get_int_option( "PRECISION", precision ) ) precision = DEFAULT_PRECISION; 00080 00081 if( MB_SUCCESS == opts.get_null_option( "STRICT" ) ) 00082 mStrict = true; 00083 else if( MB_SUCCESS == opts.get_null_option( "RELAXED" ) ) 00084 mStrict = false; 00085 else 00086 mStrict = DEFAULT_STRICT; 00087 00088 if( MB_SUCCESS == opts.get_null_option( "CREATE_ONE_NODE_CELLS" ) ) createOneNodeCells = true; 00089 00090 // Get entities to write 00091 Range nodes, elems; 00092 rval = gather_mesh( output_list, num_sets, nodes, elems ); 00093 if( MB_SUCCESS != rval ) return rval; 00094 00095 // Honor overwrite flag 00096 if( !overwrite ) 00097 { 00098 rval = writeTool->check_doesnt_exist( file_name ); 00099 if( MB_SUCCESS != rval ) return rval; 00100 } 00101 00102 // Create file 00103 std::ofstream file( file_name ); 00104 if( !file ) { MB_SET_ERR( MB_FILE_WRITE_ERROR, "Could not open file: " << file_name ); } 00105 file.precision( precision ); 00106 00107 // Write file 00108 if( ( rval = write_header( file ) ) != MB_SUCCESS || ( rval = write_nodes( file, nodes ) ) != MB_SUCCESS || 00109 ( rval = write_elems( file, nodes, elems ) ) != MB_SUCCESS || 00110 ( rval = write_tags( file, true, nodes, tag_list, num_tags ) ) != MB_SUCCESS || 00111 ( rval = write_tags( file, false, elems, tag_list, num_tags ) ) != MB_SUCCESS ) 00112 { 00113 file.close(); 00114 remove( file_name ); 00115 return rval; 00116 } 00117 00118 return MB_SUCCESS; 00119 } 00120 00121 ErrorCode WriteVtk::gather_mesh( const EntityHandle* set_list, int num_sets, Range& nodes, Range& elems ) 00122 { 00123 ErrorCode rval; 00124 int e; 00125 00126 if( !set_list || !num_sets ) 00127 { 00128 Range a; 00129 rval = mbImpl->get_entities_by_handle( 0, a ); 00130 if( MB_SUCCESS != rval ) return rval; 00131 00132 Range::const_iterator node_i, elem_i, set_i; 00133 node_i = a.lower_bound( a.begin(), a.end(), CREATE_HANDLE( MBVERTEX, 0, e ) ); 00134 elem_i = a.lower_bound( node_i, a.end(), CREATE_HANDLE( MBEDGE, 0, e ) ); 00135 set_i = a.lower_bound( elem_i, a.end(), CREATE_HANDLE( MBENTITYSET, 0, e ) ); 00136 nodes.merge( node_i, elem_i ); 00137 elems.merge( elem_i, set_i ); 00138 00139 // Filter out unsupported element types 00140 EntityType et = MBEDGE; 00141 for( et++; et < MBENTITYSET; et++ ) 00142 { 00143 if( VtkUtil::get_vtk_type( et, CN::VerticesPerEntity( et ) ) ) continue; 00144 Range::iterator eit = elems.lower_bound( elems.begin(), elems.end(), CREATE_HANDLE( et, 0, e ) ), 00145 ep1it = elems.lower_bound( elems.begin(), elems.end(), CREATE_HANDLE( et + 1, 0, e ) ); 00146 elems.erase( eit, ep1it ); 00147 } 00148 } 00149 else 00150 { 00151 std::set< EntityHandle > visited; 00152 std::vector< EntityHandle > sets; 00153 sets.reserve( num_sets ); 00154 std::copy( set_list, set_list + num_sets, std::back_inserter( sets ) ); 00155 while( !sets.empty() ) 00156 { 00157 // Get next set 00158 EntityHandle set = sets.back(); 00159 sets.pop_back(); 00160 // Skip sets we've already done 00161 if( !visited.insert( set ).second ) continue; 00162 00163 Range a; 00164 rval = mbImpl->get_entities_by_handle( set, a ); 00165 if( MB_SUCCESS != rval ) return rval; 00166 00167 Range::const_iterator node_i, elem_i, set_i; 00168 node_i = a.lower_bound( a.begin(), a.end(), CREATE_HANDLE( MBVERTEX, 0, e ) ); 00169 elem_i = a.lower_bound( node_i, a.end(), CREATE_HANDLE( MBEDGE, 0, e ) ); 00170 set_i = a.lower_bound( elem_i, a.end(), CREATE_HANDLE( MBENTITYSET, 0, e ) ); 00171 nodes.merge( node_i, elem_i ); 00172 elems.merge( elem_i, set_i ); 00173 std::copy( set_i, a.end(), std::back_inserter( sets ) ); 00174 00175 a.clear(); 00176 rval = mbImpl->get_child_meshsets( set, a ); 00177 std::copy( a.begin(), a.end(), std::back_inserter( sets ) ); 00178 } 00179 00180 for( Range::const_iterator ei = elems.begin(); ei != elems.end(); ++ei ) 00181 { 00182 std::vector< EntityHandle > connect; 00183 rval = mbImpl->get_connectivity( &( *ei ), 1, connect ); 00184 if( MB_SUCCESS != rval ) return rval; 00185 00186 for( unsigned int i = 0; i < connect.size(); ++i ) 00187 nodes.insert( connect[i] ); 00188 } 00189 } 00190 00191 if( nodes.empty() ) { MB_SET_ERR( MB_ENTITY_NOT_FOUND, "Nothing to write" ); } 00192 00193 return MB_SUCCESS; 00194 } 00195 00196 ErrorCode WriteVtk::write_header( std::ostream& stream ) 00197 { 00198 stream << "# vtk DataFile Version 3.0" << std::endl; 00199 stream << MOAB_VERSION_STRING << std::endl; 00200 stream << "ASCII" << std::endl; 00201 stream << "DATASET UNSTRUCTURED_GRID" << std::endl; 00202 return MB_SUCCESS; 00203 } 00204 00205 ErrorCode WriteVtk::write_nodes( std::ostream& stream, const Range& nodes ) 00206 { 00207 ErrorCode rval; 00208 00209 stream << "POINTS " << nodes.size() << " double" << std::endl; 00210 00211 double coords[3]; 00212 for( Range::const_iterator i = nodes.begin(); i != nodes.end(); ++i ) 00213 { 00214 coords[1] = coords[2] = 0.0; 00215 rval = mbImpl->get_coords( &( *i ), 1, coords ); 00216 if( MB_SUCCESS != rval ) return rval; 00217 stream << coords[0] << ' ' << coords[1] << ' ' << coords[2] << std::endl; 00218 } 00219 00220 return MB_SUCCESS; 00221 } 00222 00223 ErrorCode WriteVtk::write_elems( std::ostream& stream, const Range& nodes, const Range& elems ) 00224 { 00225 ErrorCode rval; 00226 00227 Range connectivity; // because we now support polyhedra, it could contain faces 00228 rval = mbImpl->get_connectivity( elems, connectivity );MB_CHK_ERR( rval ); 00229 00230 Range nodes_from_connectivity = connectivity.subset_by_type( MBVERTEX ); 00231 Range faces_from_connectivity = 00232 subtract( connectivity, nodes_from_connectivity ); // these could be faces of polyhedra 00233 00234 Range connected_nodes; 00235 rval = mbImpl->get_connectivity( faces_from_connectivity, connected_nodes );MB_CHK_ERR( rval ); 00236 connected_nodes.merge( nodes_from_connectivity ); 00237 00238 Range free_nodes = subtract( nodes, connected_nodes ); 00239 00240 // Get and write counts 00241 unsigned long num_elems, num_uses; 00242 num_elems = num_uses = elems.size(); 00243 00244 std::map< EntityHandle, int > sizeFieldsPolyhedra; 00245 00246 for( Range::const_iterator i = elems.begin(); i != elems.end(); ++i ) 00247 { 00248 EntityType type = mbImpl->type_from_handle( *i ); 00249 if( !VtkUtil::get_vtk_type( type, CN::VerticesPerEntity( type ) ) ) continue; 00250 00251 EntityHandle elem = *i; 00252 const EntityHandle* connect = NULL; 00253 int conn_len = 0; 00254 // Dummy storage vector for structured mesh "get_connectivity" function 00255 std::vector< EntityHandle > storage; 00256 rval = mbImpl->get_connectivity( elem, connect, conn_len, false, &storage );MB_CHK_ERR( rval ); 00257 00258 num_uses += conn_len; 00259 // if polyhedra, we will count the number of nodes in each face too 00260 if( TYPE_FROM_HANDLE( elem ) == MBPOLYHEDRON ) 00261 { 00262 int numFields = 1; // there will be one for number of faces; forgot about this one 00263 for( int j = 0; j < conn_len; j++ ) 00264 { 00265 const EntityHandle* conn = NULL; 00266 int num_nd = 0; 00267 rval = mbImpl->get_connectivity( connect[j], conn, num_nd );MB_CHK_ERR( rval ); 00268 numFields += num_nd + 1; 00269 } 00270 sizeFieldsPolyhedra[elem] = numFields; // will be used later, at writing 00271 num_uses += ( numFields - conn_len ); 00272 } 00273 } 00274 freeNodes = (int)free_nodes.size(); 00275 if( !createOneNodeCells ) freeNodes = 0; // do not create one node cells 00276 stream << "CELLS " << num_elems + freeNodes << ' ' << num_uses + 2 * freeNodes << std::endl; 00277 00278 // Write element connectivity 00279 std::vector< int > conn_data; 00280 std::vector< unsigned > vtk_types( elems.size() + freeNodes ); 00281 std::vector< unsigned >::iterator t = vtk_types.begin(); 00282 for( Range::const_iterator i = elems.begin(); i != elems.end(); ++i ) 00283 { 00284 // Get type information for element 00285 EntityHandle elem = *i; 00286 EntityType type = TYPE_FROM_HANDLE( elem ); 00287 00288 // Get element connectivity 00289 const EntityHandle* connect = NULL; 00290 int conn_len = 0; 00291 // Dummy storage vector for structured mesh "get_connectivity" function 00292 std::vector< EntityHandle > storage; 00293 rval = mbImpl->get_connectivity( elem, connect, conn_len, false, &storage );MB_CHK_ERR( rval ); 00294 00295 // Get VTK type 00296 const VtkElemType* vtk_type = VtkUtil::get_vtk_type( type, conn_len ); 00297 if( !vtk_type ) 00298 { 00299 // Try connectivity with 1 fewer node 00300 vtk_type = VtkUtil::get_vtk_type( type, conn_len - 1 ); 00301 if( vtk_type ) 00302 conn_len--; 00303 else 00304 { 00305 MB_SET_ERR( MB_FAILURE, "Vtk file format does not support elements of type " 00306 << CN::EntityTypeName( type ) << " (" << (int)type << ") with " << conn_len 00307 << " nodes" ); 00308 } 00309 } 00310 00311 // Save VTK type index for later 00312 *t = vtk_type->vtk_type; 00313 ++t; 00314 00315 if( type != MBPOLYHEDRON ) 00316 { 00317 // Get IDs from vertex handles 00318 assert( conn_len > 0 ); 00319 conn_data.resize( conn_len ); 00320 for( int j = 0; j < conn_len; ++j ) 00321 conn_data[j] = nodes.index( connect[j] ); 00322 00323 // Write connectivity list 00324 stream << conn_len; 00325 if( vtk_type->node_order ) 00326 for( int k = 0; k < conn_len; ++k ) 00327 stream << ' ' << conn_data[vtk_type->node_order[k]]; 00328 else 00329 for( int k = 0; k < conn_len; ++k ) 00330 stream << ' ' << conn_data[k]; 00331 stream << std::endl; 00332 } 00333 else 00334 { 00335 // POLYHEDRON needs a special case, loop over faces to get nodes 00336 stream << sizeFieldsPolyhedra[elem] << " " << conn_len; 00337 for( int k = 0; k < conn_len; k++ ) 00338 { 00339 EntityHandle face = connect[k]; 00340 const EntityHandle* conn = NULL; 00341 int num_nodes = 0; 00342 rval = mbImpl->get_connectivity( face, conn, num_nodes );MB_CHK_ERR( rval ); 00343 // num_uses += num_nd + 1; // 1 for number of vertices in face 00344 conn_data.resize( num_nodes ); 00345 for( int j = 0; j < num_nodes; ++j ) 00346 conn_data[j] = nodes.index( conn[j] ); 00347 00348 stream << ' ' << num_nodes; 00349 00350 for( int j = 0; j < num_nodes; ++j ) 00351 stream << ' ' << conn_data[j]; 00352 } 00353 stream << std::endl; 00354 } 00355 } 00356 00357 if( createOneNodeCells ) 00358 for( Range::const_iterator v = free_nodes.begin(); v != free_nodes.end(); ++v, ++t ) 00359 { 00360 EntityHandle node = *v; 00361 stream << "1 " << nodes.index( node ) << std::endl; 00362 *t = 1; 00363 } 00364 00365 // Write element types 00366 stream << "CELL_TYPES " << vtk_types.size() << std::endl; 00367 for( std::vector< unsigned >::const_iterator i = vtk_types.begin(); i != vtk_types.end(); ++i ) 00368 stream << *i << std::endl; 00369 00370 return MB_SUCCESS; 00371 } 00372 00373 ErrorCode WriteVtk::write_tags( std::ostream& stream, bool nodes, const Range& entities, const Tag* tag_list, 00374 int num_tags ) 00375 { 00376 ErrorCode rval; 00377 00378 // The #$%@#$% MOAB interface does not have a function to retrieve 00379 // all entities with a tag, only all entities with a specified type 00380 // and tag. Define types to loop over depending on the if vertex 00381 // or element tag data is being written. It seems horribly inefficient 00382 // to have the implementation subdivide the results by type and have 00383 // to call that function once for each type just to recombine the results. 00384 // Unfortunately, there doesn't seem to be any other way. 00385 EntityType low_type, high_type; 00386 if( nodes ) 00387 { 00388 low_type = MBVERTEX; 00389 high_type = MBEDGE; 00390 } 00391 else 00392 { 00393 low_type = MBEDGE; 00394 high_type = MBENTITYSET; 00395 } 00396 00397 // Get all defined tags 00398 std::vector< Tag > tags; 00399 std::vector< Tag >::iterator i; 00400 rval = writeTool->get_tag_list( tags, tag_list, num_tags, false ); 00401 if( MB_SUCCESS != rval ) return rval; 00402 00403 // For each tag... 00404 bool entities_have_tags = false; 00405 for( i = tags.begin(); i != tags.end(); ++i ) 00406 { 00407 // Skip tags holding entity handles -- no way to save them 00408 DataType dtype; 00409 rval = mbImpl->tag_get_data_type( *i, dtype ); 00410 if( MB_SUCCESS != rval ) return rval; 00411 if( dtype == MB_TYPE_HANDLE ) continue; 00412 00413 // If in strict mode, don't write tags that do not fit in any 00414 // attribute type (SCALAR : 1 to 4 values, VECTOR : 3 values, TENSOR : 9 values) 00415 if( mStrict ) 00416 { 00417 int count; 00418 rval = mbImpl->tag_get_length( *i, count ); 00419 if( MB_SUCCESS != rval ) return rval; 00420 if( count < 1 || ( count > 4 && count != 9 ) ) continue; 00421 } 00422 00423 // Get subset of input entities that have the tag set 00424 Range tagged; 00425 for( EntityType type = low_type; type < high_type; ++type ) 00426 { 00427 Range tmp_tagged; 00428 rval = mbImpl->get_entities_by_type_and_tag( 0, type, &( *i ), 0, 1, tmp_tagged ); 00429 if( MB_SUCCESS != rval ) return rval; 00430 tmp_tagged = intersect( tmp_tagged, entities ); 00431 tagged.merge( tmp_tagged ); 00432 } 00433 00434 // If any entities were tagged 00435 if( !tagged.empty() ) 00436 { 00437 // If this is the first tag being written for the 00438 // entity type, write the label marking the beginning 00439 // of the tag data. 00440 if( !entities_have_tags ) 00441 { 00442 entities_have_tags = true; 00443 if( nodes ) 00444 stream << "POINT_DATA " << entities.size() << std::endl; 00445 else 00446 stream << "CELL_DATA " << entities.size() + freeNodes << std::endl; 00447 } 00448 00449 // Write the tag 00450 rval = write_tag( stream, *i, entities, tagged ); 00451 if( MB_SUCCESS != rval ) return rval; 00452 } 00453 } 00454 00455 return MB_SUCCESS; 00456 } 00457 00458 template < typename T > 00459 void WriteVtk::write_data( std::ostream& stream, const std::vector< T >& data, unsigned vals_per_tag ) 00460 { 00461 typename std::vector< T >::const_iterator d = data.begin(); 00462 const unsigned long n = data.size() / vals_per_tag; 00463 00464 for( unsigned long i = 0; i < n; ++i ) 00465 { 00466 for( unsigned j = 0; j < vals_per_tag; ++j, ++d ) 00467 { 00468 if( sizeof( T ) == 1 ) 00469 stream << (unsigned int)*d << ' '; 00470 else 00471 stream << *d << ' '; 00472 } 00473 stream << std::endl; 00474 } 00475 } 00476 00477 // template <> 00478 // void WriteVtk::write_data<unsigned char>(std::ostream& stream, 00479 // const std::vector<unsigned char>& data, 00480 // unsigned vals_per_tag) 00481 //{ 00482 // std::vector<unsigned char>::const_iterator d = data.begin(); 00483 // const unsigned long n = data.size() / vals_per_tag; 00484 // 00485 // for (unsigned long i = 0; i < n; ++i) { 00486 // for (unsigned j = 0; j < vals_per_tag; ++j, ++d) 00487 // stream << (unsigned int)*d << ' '; 00488 // stream << std::endl; 00489 // } 00490 //} 00491 00492 template < typename T > 00493 ErrorCode WriteVtk::write_tag( std::ostream& stream, Tag tag, const Range& entities, const Range& tagged, const int ) 00494 { 00495 ErrorCode rval; 00496 int addFreeNodes = 0; 00497 if( TYPE_FROM_HANDLE( entities[0] ) > MBVERTEX ) addFreeNodes = freeNodes; 00498 // we created freeNodes 1-node cells, so we have to augment cell data too 00499 // we know that the 1 node cells are added at the end, after all other cells; 00500 // so the default values will be set to those extra , artificial cells 00501 const unsigned long n = entities.size() + addFreeNodes; 00502 00503 // Get tag properties 00504 00505 std::string name; 00506 int vals_per_tag; 00507 if( MB_SUCCESS != mbImpl->tag_get_name( tag, name ) || MB_SUCCESS != mbImpl->tag_get_length( tag, vals_per_tag ) ) 00508 return MB_FAILURE; 00509 00510 // Get a tag value for each entity. Do this by initializing the 00511 // "data" vector with zero, and then filling in the values for 00512 // the entities that actually have the tag set. 00513 std::vector< T > data; 00514 data.resize( n * vals_per_tag, 0 ); 00515 // If there is a default value for the tag, set the actual default value 00516 std::vector< T > def_value( vals_per_tag ); 00517 rval = mbImpl->tag_get_default_value( tag, &( def_value[0] ) ); 00518 if( MB_SUCCESS == rval ) SysUtil::setmem( &( data[0] ), &( def_value[0] ), vals_per_tag * sizeof( T ), n ); 00519 00520 Range::const_iterator t = tagged.begin(); 00521 typename std::vector< T >::iterator d = data.begin(); 00522 for( Range::const_iterator i = entities.begin(); i != entities.end() && t != tagged.end(); ++i, d += vals_per_tag ) 00523 { 00524 if( *i == *t ) 00525 { 00526 ++t; 00527 rval = mbImpl->tag_get_data( tag, &( *i ), 1, &( *d ) ); 00528 if( MB_SUCCESS != rval ) return rval; 00529 } 00530 } 00531 00532 // Write the tag values, one entity per line. 00533 write_data( stream, data, vals_per_tag ); 00534 00535 return MB_SUCCESS; 00536 } 00537 00538 ErrorCode WriteVtk::write_bit_tag( std::ostream& stream, Tag tag, const Range& entities, const Range& tagged ) 00539 { 00540 ErrorCode rval; 00541 const unsigned long n = entities.size(); 00542 00543 // Get tag properties 00544 00545 std::string name; 00546 int vals_per_tag; 00547 if( MB_SUCCESS != mbImpl->tag_get_name( tag, name ) || MB_SUCCESS != mbImpl->tag_get_length( tag, vals_per_tag ) ) 00548 return MB_FAILURE; 00549 00550 if( vals_per_tag > 8 ) { MB_SET_ERR( MB_FAILURE, "Invalid tag size for bit tag \"" << name << "\"" ); } 00551 00552 // Get a tag value for each entity. 00553 // Get bits for each entity and unpack into 00554 // one integer in the 'data' array for each bit. 00555 // Initialize 'data' to zero because we will skip 00556 // those entities for which the tag is not set. 00557 std::vector< unsigned short > data; 00558 data.resize( n * vals_per_tag, 0 ); 00559 Range::const_iterator t = tagged.begin(); 00560 std::vector< unsigned short >::iterator d = data.begin(); 00561 for( Range::const_iterator i = entities.begin(); i != entities.end() && t != tagged.end(); ++i ) 00562 { 00563 if( *i == *t ) 00564 { 00565 ++t; 00566 unsigned char value; 00567 rval = mbImpl->tag_get_data( tag, &( *i ), 1, &value ); 00568 for( int j = 0; j < vals_per_tag; ++j, ++d ) 00569 *d = (unsigned short)( value & ( 1 << j ) ? 1 : 0 ); 00570 if( MB_SUCCESS != rval ) return rval; 00571 } 00572 else 00573 { 00574 // If tag is not set for entity, skip values in array 00575 d += vals_per_tag; 00576 } 00577 } 00578 00579 // Write the tag values, one entity per line. 00580 write_data( stream, data, vals_per_tag ); 00581 00582 return MB_SUCCESS; 00583 } 00584 00585 ErrorCode WriteVtk::write_tag( std::ostream& s, Tag tag, const Range& entities, const Range& tagged ) 00586 { 00587 // Get tag properties 00588 std::string name; 00589 DataType type; 00590 int vals_per_tag; 00591 if( MB_SUCCESS != mbImpl->tag_get_name( tag, name ) || MB_SUCCESS != mbImpl->tag_get_length( tag, vals_per_tag ) || 00592 MB_SUCCESS != mbImpl->tag_get_data_type( tag, type ) ) 00593 return MB_FAILURE; 00594 00595 // Skip tags of type ENTITY_HANDLE 00596 if( MB_TYPE_HANDLE == type ) return MB_FAILURE; 00597 00598 // Now that we're past the point where the name would be used in 00599 // an error message, remove any spaces to conform to VTK file. 00600 for( std::string::iterator i = name.begin(); i != name.end(); ++i ) 00601 { 00602 if( isspace( *i ) || iscntrl( *i ) ) *i = '_'; 00603 } 00604 00605 // Write the tag description 00606 if( 3 == vals_per_tag && MB_TYPE_DOUBLE == type ) 00607 s << "VECTORS " << name << ' ' << VtkUtil::vtkTypeNames[type] << std::endl; 00608 else if( 9 == vals_per_tag ) 00609 s << "TENSORS " << name << ' ' << VtkUtil::vtkTypeNames[type] << std::endl; 00610 else 00611 s << "SCALARS " << name << ' ' << VtkUtil::vtkTypeNames[type] << ' ' << vals_per_tag << std::endl 00612 << "LOOKUP_TABLE default" << std::endl; 00613 00614 // Write the tag data 00615 switch( type ) 00616 { 00617 case MB_TYPE_OPAQUE: 00618 return write_tag< unsigned char >( s, tag, entities, tagged, 0 ); 00619 case MB_TYPE_INTEGER: 00620 return write_tag< int >( s, tag, entities, tagged, 0 ); 00621 case MB_TYPE_DOUBLE: 00622 return write_tag< double >( s, tag, entities, tagged, 0 ); 00623 case MB_TYPE_BIT: 00624 return write_bit_tag( s, tag, entities, tagged ); 00625 default: 00626 return MB_FAILURE; 00627 } 00628 } 00629 00630 } // namespace moab