moab
|
#include <ReadVtk.hpp>
Public Member Functions | |
ErrorCode | load_file (const char *file_name, const EntityHandle *file_set, const FileOptions &opts, const SubsetList *subset_list=0, const Tag *file_id_tag=0) |
load a file | |
ErrorCode | read_tag_values (const char *file_name, const char *tag_name, const FileOptions &opts, std::vector< int > &tag_values_out, const SubsetList *subset_list=0) |
Read tag values from a file. | |
ReadVtk (Interface *impl=NULL) | |
Constructor. | |
virtual | ~ReadVtk () |
Destructor. | |
Static Public Member Functions | |
static ReaderIface * | factory (Interface *) |
Protected Member Functions | |
ErrorCode | allocate_vertices (long num_vtx, EntityHandle &start_handle_out, double *&x_coord_array_out, double *&y_coord_array_out, double *&z_coord_array_out) |
ErrorCode | read_vertices (FileTokenizer &tokens, long num_verts, EntityHandle &start_handle_out) |
ErrorCode | allocate_elements (long num_elements, int vert_per_element, EntityType type, EntityHandle &start_handle_out, EntityHandle *&conn_array_out, std::vector< Range > &append_to_this) |
ErrorCode | vtk_read_dataset (FileTokenizer &tokens, Range &vertex_list, std::vector< Range > &element_list) |
ErrorCode | vtk_read_structured_points (FileTokenizer &tokens, Range &vertex_list, std::vector< Range > &elem_list) |
ErrorCode | vtk_read_structured_grid (FileTokenizer &tokens, Range &vertex_list, std::vector< Range > &elem_list) |
ErrorCode | vtk_read_rectilinear_grid (FileTokenizer &tokens, Range &vertex_list, std::vector< Range > &elem_list) |
ErrorCode | vtk_read_polydata (FileTokenizer &tokens, Range &vertex_list, std::vector< Range > &elem_list) |
ErrorCode | vtk_read_polygons (FileTokenizer &tokens, EntityHandle first_vtx, std::vector< Range > &elem_list) |
ErrorCode | vtk_read_unstructured_grid (FileTokenizer &tokens, Range &vertex_list, std::vector< Range > &elem_list) |
ErrorCode | vtk_create_structured_elems (const long *dims, EntityHandle first_vtx, std::vector< Range > &elem_list) |
ErrorCode | vtk_read_field (FileTokenizer &tokens) |
ErrorCode | vtk_read_attrib_data (FileTokenizer &tokens, std::vector< Range > &entities) |
ErrorCode | vtk_read_tag_data (FileTokenizer &tokens, int type, size_t per_elem, std::vector< Range > &entities, const char *name) |
ErrorCode | vtk_read_scalar_attrib (FileTokenizer &tokens, std::vector< Range > &entities, const char *name) |
ErrorCode | vtk_read_color_attrib (FileTokenizer &tokens, std::vector< Range > &entities, const char *name) |
ErrorCode | vtk_read_vector_attrib (FileTokenizer &tokens, std::vector< Range > &entities, const char *name) |
ErrorCode | vtk_read_texture_attrib (FileTokenizer &tokens, std::vector< Range > &entities, const char *name) |
ErrorCode | vtk_read_tensor_attrib (FileTokenizer &tokens, std::vector< Range > &entities, const char *name) |
ErrorCode | vtk_read_field_attrib (FileTokenizer &tokens, std::vector< Range > &entities, const char *name) |
ErrorCode | store_file_ids (Tag tag, const Range &vertices, const std::vector< Range > &elements) |
Private Attributes | |
ReadUtilIface * | readMeshIface |
Interface * | mdbImpl |
interface instance | |
EntityHandle | mCurrentMeshHandle |
Meshset Handle for the mesh that is currently being read. | |
std::string | mPartitionTagName |
A field which, if present and having a single integer for storage, should be used to partition the mesh by range. Defaults to MATERIAL_SET_TAG_NAME. |
Definition at line 29 of file ReadVtk.hpp.
ReadVtk::ReadVtk | ( | Interface * | impl = NULL | ) |
Constructor.
Definition at line 143 of file ReadVtk.cpp.
: mdbImpl(impl), mPartitionTagName(MATERIAL_SET_TAG_NAME) { mdbImpl->query_interface(readMeshIface); }
ReadVtk::~ReadVtk | ( | ) | [virtual] |
Destructor.
Definition at line 149 of file ReadVtk.cpp.
{ if (readMeshIface) { mdbImpl->release_interface(readMeshIface); readMeshIface = 0; } }
ErrorCode ReadVtk::allocate_elements | ( | long | num_elements, |
int | vert_per_element, | ||
EntityType | type, | ||
EntityHandle & | start_handle_out, | ||
EntityHandle *& | conn_array_out, | ||
std::vector< Range > & | append_to_this | ||
) | [protected] |
Definition at line 378 of file ReadVtk.cpp.
{ ErrorCode result; start_handle_out = 0; result = readMeshIface->get_element_connect( num_elements, vert_per_element, type, MB_START_ID, start_handle_out, conn_array_out ); if (MB_SUCCESS != result) return result; Range range(start_handle_out, start_handle_out+num_elements-1); append_to_this.push_back( range ); return MB_SUCCESS; }
ErrorCode ReadVtk::allocate_vertices | ( | long | num_vtx, |
EntityHandle & | start_handle_out, | ||
double *& | x_coord_array_out, | ||
double *& | y_coord_array_out, | ||
double *& | z_coord_array_out | ||
) | [protected] |
Definition at line 332 of file ReadVtk.cpp.
{ ErrorCode result; // Create vertices std::vector<double*> arrays; start_handle_out = 0; result = readMeshIface->get_node_coords( 3, num_verts, MB_START_ID, start_handle_out, arrays ); if (MB_SUCCESS != result) return result; x_coord_array_out = arrays[0]; y_coord_array_out = arrays[1]; z_coord_array_out = arrays[2]; return MB_SUCCESS; }
ReaderIface * ReadVtk::factory | ( | Interface * | iface | ) | [static] |
Definition at line 140 of file ReadVtk.cpp.
ErrorCode ReadVtk::load_file | ( | const char * | file_name, |
const EntityHandle * | file_set, | ||
const FileOptions & | opts, | ||
const SubsetList * | subset_list = 0 , |
||
const Tag * | file_id_tag = 0 |
||
) | [virtual] |
load a file
Implements moab::ReaderIface.
Definition at line 182 of file ReadVtk.cpp.
{ ErrorCode result; int major, minor; char vendor_string[257]; std::vector<Range> element_list; Range vertices; if (subset_list) { readMeshIface->report_error( "Reading subset of files not supported for VTK." ); return MB_UNSUPPORTED_OPERATION; } // Does the caller want a field to be used for partitioning the entities? // If not, we'll assume any scalar integer field named MATERIAL_SET specifies partitions. std::string partition_tag_name; result = opts.get_option( "PARTITION", partition_tag_name ); if ( result == MB_SUCCESS ) mPartitionTagName = partition_tag_name; FILE* file = fopen( filename, "r" ); if (!file) { return MB_FILE_DOES_NOT_EXIST; } // Read file header if (!fgets( vendor_string, sizeof(vendor_string), file )) { fclose( file ); return MB_FAILURE; } if (!strchr( vendor_string, '\n' ) || 2 != sscanf( vendor_string, "# vtk DataFile Version %d.%d", &major, &minor )) { fclose( file ); return MB_FAILURE; } if (!fgets( vendor_string, sizeof(vendor_string), file )) { fclose( file ); return MB_FAILURE; } // VTK spec says this should not exceed 256 chars. if (!strchr( vendor_string, '\n' )) { readMeshIface->report_error( "Vendor string (line 2) exceeds 256 characters." ); fclose( file ); return MB_FAILURE; } // Check file type FileTokenizer tokens( file, readMeshIface ); const char* const file_type_names[] = { "ASCII", "BINARY", 0 }; int filetype = tokens.match_token( file_type_names ); switch (filetype) { case 2: // BINARY readMeshIface->report_error( "Cannot read BINARY VTK files" ); default: // ERROR return MB_FAILURE; case 1: // ASCII break; } // Read the mesh if (!tokens.match_token( "DATASET" )) return MB_FAILURE; result = vtk_read_dataset( tokens, vertices, element_list ); if (MB_SUCCESS != result) return result; if (file_id_tag) { result = store_file_ids( *file_id_tag, vertices, element_list ); if (MB_SUCCESS != result) return result; } // Count the number of elements read long elem_count = 0; for (std::vector<Range>::iterator it = element_list.begin(); it != element_list.end(); ++it ) elem_count += it->size(); // Read attribute data until end of file. const char* const block_type_names[] = { "POINT_DATA", "CELL_DATA", 0 }; std::vector<Range> vertex_list(1); vertex_list[0] = vertices; int blocktype = 0; while (!tokens.eof()) { // get POINT_DATA or CELL_DATA int new_block_type = tokens.match_token( block_type_names, false ); if (tokens.eof()) break; if (!new_block_type) { // If next token was neither POINT_DATA nor CELL_DATA, // then there's another attribute under the current one. if (blocktype) tokens.unget_token(); else break; } else { blocktype = new_block_type; long count; if (!tokens.get_long_ints( 1, &count )) return MB_FAILURE; if (blocktype == 1 && (unsigned long)count != vertices.size()) { readMeshIface->report_error( "Count inconsistent with number of vertices at line %d.", tokens.line_number()); return MB_FAILURE; } else if (blocktype == 2 && count != elem_count) { readMeshIface->report_error( "Count inconsistent with number of elements at line %d.", tokens.line_number()); return MB_FAILURE; } } if (blocktype == 1) result = vtk_read_attrib_data( tokens, vertex_list ); else result = vtk_read_attrib_data( tokens, element_list ); if (MB_SUCCESS != result) return result; } return MB_SUCCESS; }
ErrorCode ReadVtk::read_tag_values | ( | const char * | file_name, |
const char * | tag_name, | ||
const FileOptions & | opts, | ||
std::vector< int > & | tag_values_out, | ||
const SubsetList * | subset_list = 0 |
||
) | [virtual] |
Read tag values from a file.
Read the list if all integer tag values from the file for a tag that is a single integer value per entity.
file_name | The file to read. |
tag_name | The tag for which to read values |
tag_values_out | Output: The list of tag values. |
subset_list | An array of tag name and value sets specifying the subset of the file to read. If multiple tags are specified, the sets that match all tags (intersection) should be read. |
subset_list_length | The length of the 'subset_list' array. |
Implements moab::ReaderIface.
Definition at line 172 of file ReadVtk.cpp.
{ return MB_NOT_IMPLEMENTED; }
ErrorCode ReadVtk::read_vertices | ( | FileTokenizer & | tokens, |
long | num_verts, | ||
EntityHandle & | start_handle_out | ||
) | [protected] |
Definition at line 355 of file ReadVtk.cpp.
{ ErrorCode result; double *x, *y, *z; result = allocate_vertices( num_verts, start_handle_out, x, y, z ); if (MB_SUCCESS != result) return result; // Read vertex coordinates for (long vtx = 0; vtx < num_verts; ++vtx) { if (!tokens.get_doubles( 1, x++ ) || !tokens.get_doubles( 1, y++ ) || !tokens.get_doubles( 1, z++ )) return MB_FAILURE; } return MB_SUCCESS; }
ErrorCode ReadVtk::store_file_ids | ( | Tag | tag, |
const Range & | vertices, | ||
const std::vector< Range > & | elements | ||
) | [protected] |
Definition at line 1284 of file ReadVtk.cpp.
{ ErrorCode rval; rval = readMeshIface->assign_ids( tag, verts ); if (MB_SUCCESS != rval) return rval; int id = 0; for (size_t i = 0; i < elems.size(); ++i) { rval = readMeshIface->assign_ids( tag, elems[i], id ); id += elems[i].size(); } return MB_SUCCESS; }
ErrorCode ReadVtk::vtk_create_structured_elems | ( | const long * | dims, |
EntityHandle | first_vtx, | ||
std::vector< Range > & | elem_list | ||
) | [protected] |
Definition at line 878 of file ReadVtk.cpp.
{ ErrorCode result; //int non_zero[3] = {0,0,0}; // True if dim > 0 for x, y, z respectively long elem_dim = 0; // Element dimension (2->quad, 3->hex) long num_elems = 1; // Total number of elements long vert_per_elem; // Element connectivity length long edims[3] = { 1, 1, 1 };// Number of elements in each grid direction // Populate above data for (int d = 0; d < 3; d++) if (dims[d] > 1) { //non_zero[elem_dim] = d; ++elem_dim; edims[d] = dims[d] - 1; num_elems *= edims[d]; } vert_per_elem = 1 << elem_dim; // Get element type from element dimension EntityType type; switch( elem_dim ) { case 1: type = MBEDGE; break; case 2: type = MBQUAD; break; case 3: type = MBHEX; break; default: readMeshIface->report_error( "Invalid dimension for structured elements: %ld\n", elem_dim ); return MB_FAILURE; } // Allocate storage for elements EntityHandle start_handle = 0; EntityHandle* conn_array; result = allocate_elements( num_elems, vert_per_elem, type, start_handle, conn_array, elem_list ); if (MB_SUCCESS != result) return MB_FAILURE; EntityHandle *conn_sav = conn_array; // Offsets of element vertices in grid relative to corner closest to origin long k = dims[0]*dims[1]; const long corners[8] = { 0, 1, 1+dims[0], dims[0], k, k+1, k+1+dims[0], k+dims[0] }; // Populate element list for (long z = 0; z < edims[2]; ++z) for (long y = 0; y < edims[1]; ++y) for (long x = 0; x < edims[0]; ++x) { const long index = x + y*dims[0] + z*(dims[0]*dims[1]); for (long j = 0; j < vert_per_elem; ++j, ++conn_array) *conn_array = index + corners[j] + first_vtx; } // notify MOAB of the new elements result = readMeshIface->update_adjacencies(start_handle, num_elems, vert_per_elem, conn_sav); if (MB_SUCCESS != result) return result; return MB_SUCCESS; }
ErrorCode ReadVtk::vtk_read_attrib_data | ( | FileTokenizer & | tokens, |
std::vector< Range > & | entities | ||
) | [protected] |
Definition at line 986 of file ReadVtk.cpp.
{ const char* const type_names[] = { "SCALARS", "COLOR_SCALARS", "VECTORS", "NORMALS", "TEXTURE_COORDINATES", "TENSORS", "FIELD", 0 }; int type = tokens.match_token( type_names ); const char* tmp_name = tokens.get_string( ); if (!type || !tmp_name) return MB_FAILURE; std::string name_alloc( tmp_name ); const char* name = name_alloc.c_str(); switch( type ) { case 1: return vtk_read_scalar_attrib ( tokens, entities, name ); case 2: return vtk_read_color_attrib ( tokens, entities, name ); case 3: return vtk_read_vector_attrib ( tokens, entities, name ); case 4: return vtk_read_vector_attrib ( tokens, entities, name ); case 5: return vtk_read_texture_attrib( tokens, entities, name ); case 6: return vtk_read_tensor_attrib ( tokens, entities, name ); case 7: return vtk_read_field_attrib ( tokens, entities, name ); } return MB_FAILURE; }
ErrorCode ReadVtk::vtk_read_color_attrib | ( | FileTokenizer & | tokens, |
std::vector< Range > & | entities, | ||
const char * | name | ||
) | [protected] |
Definition at line 1190 of file ReadVtk.cpp.
{ long size; if (!tokens.get_long_ints( 1, &size ) || size < 1) return MB_FAILURE; return vtk_read_tag_data( tokens, 10, size, entities, name ); }
ErrorCode ReadVtk::vtk_read_dataset | ( | FileTokenizer & | tokens, |
Range & | vertex_list, | ||
std::vector< Range > & | element_list | ||
) | [protected] |
Definition at line 402 of file ReadVtk.cpp.
{ const char* const data_type_names[] = { "STRUCTURED_POINTS", "STRUCTURED_GRID", "UNSTRUCTURED_GRID", "POLYDATA", "RECTILINEAR_GRID", "FIELD", 0 }; int datatype = tokens.match_token( data_type_names ); switch( datatype ) { case 1: return vtk_read_structured_points( tokens, vertex_list, element_list ); case 2: return vtk_read_structured_grid ( tokens, vertex_list, element_list ); case 3: return vtk_read_unstructured_grid( tokens, vertex_list, element_list ); case 4: return vtk_read_polydata ( tokens, vertex_list, element_list ); case 5: return vtk_read_rectilinear_grid ( tokens, vertex_list, element_list ); case 6: return vtk_read_field ( tokens ); default: return MB_FAILURE; } }
ErrorCode ReadVtk::vtk_read_field | ( | FileTokenizer & | tokens | ) | [protected] |
Definition at line 945 of file ReadVtk.cpp.
{ // This is not supported yet. // Parse the data but throw it away because // Mesquite has no internal representation for it. // Could save this in tags, but the only useful thing that // could be done with the data is to write it back out // with the modified mesh. As there's no way to save the // type of a tag in Mesquite, it cannot be written back // out correctly either. // FIXME: Don't know what to do with this data. // For now, read it and throw it out. long num_arrays; if (!tokens.get_string() || // name !tokens.get_long_ints( 1, &num_arrays )) return MB_FAILURE; for (long i = 0; i < num_arrays; ++i) { /*const char* name =*/ tokens.get_string( ); long dims[2]; if (!tokens.get_long_ints( 2, dims ) || !tokens.match_token( vtk_type_names )) return MB_FAILURE; long num_vals = dims[0] * dims[1]; for (long j = 0; j < num_vals; j++) { double junk; if (!tokens.get_doubles( 1, &junk )) return MB_FAILURE; } } return MB_SUCCESS; }
ErrorCode ReadVtk::vtk_read_field_attrib | ( | FileTokenizer & | tokens, |
std::vector< Range > & | entities, | ||
const char * | name | ||
) | [protected] |
Definition at line 1242 of file ReadVtk.cpp.
{ long num_fields; if (!tokens.get_long_ints(1, &num_fields)) return MB_FAILURE; long i; for ( i = 0; i < num_fields; ++ i ) { const char* tok = tokens.get_string( ); if ( ! tok ) return MB_FAILURE; std::string name_alloc( tok ); long num_comp; if (!tokens.get_long_ints( 1, &num_comp)) return MB_FAILURE; long num_tuples; if (!tokens.get_long_ints( 1, &num_tuples)) return MB_FAILURE; int type = tokens.match_token( vtk_type_names ); if (!type) return MB_FAILURE; ErrorCode result; if ( ( result = vtk_read_tag_data( tokens, type, num_comp, entities, name_alloc.c_str() ) ) != MB_SUCCESS ) { readMeshIface->report_error( "Error reading data for field \"%s\" (%ld components, %ld tuples, type %d) at line %d", name_alloc.c_str(), num_comp, num_tuples, (int)type, tokens.line_number()); return result; } } return MB_SUCCESS; }
ErrorCode ReadVtk::vtk_read_polydata | ( | FileTokenizer & | tokens, |
Range & | vertex_list, | ||
std::vector< Range > & | elem_list | ||
) | [protected] |
Definition at line 590 of file ReadVtk.cpp.
{ ErrorCode result; long num_verts; std::vector<int> connectivity; const char* const poly_data_names[] = { "VERTICES", "LINES", "POLYGONS", "TRIANGLE_STRIPS", 0 }; if (!tokens.match_token( "POINTS" ) || !tokens.get_long_ints( 1, &num_verts ) || !tokens.match_token( vtk_type_names ) || !tokens.get_newline( )) return MB_FAILURE; if (num_verts < 1) { readMeshIface->report_error( "Invalid point count at line %d", tokens.line_number()); return MB_FAILURE; } // Create vertices and read coordinates EntityHandle start_handle = 0; result = read_vertices( tokens, num_verts, start_handle ); if (MB_SUCCESS != result) return result; vertex_list.insert( start_handle, start_handle + num_verts - 1 ); int poly_type = tokens.match_token( poly_data_names ); switch (poly_type) { case 0: result = MB_FAILURE; break; case 1: readMeshIface->report_error( "Vertex element type at line %d", tokens.line_number() ); result = MB_FAILURE; break; case 2: readMeshIface->report_error( "Unsupported type: polylines at line %d", tokens.line_number() ); result = MB_FAILURE; break; case 3: result = vtk_read_polygons( tokens, start_handle, elem_list ); break; case 4: readMeshIface->report_error( "Unsupported type: triangle strips at line %d", tokens.line_number() ); result = MB_FAILURE; break; } return result; }
ErrorCode ReadVtk::vtk_read_polygons | ( | FileTokenizer & | tokens, |
EntityHandle | first_vtx, | ||
std::vector< Range > & | elem_list | ||
) | [protected] |
Definition at line 655 of file ReadVtk.cpp.
{ ErrorCode result; long size[2]; if (!tokens.get_long_ints( 2, size ) || !tokens.get_newline( )) return MB_FAILURE; const Range empty; std::vector<EntityHandle> conn_hdl; std::vector<long> conn_idx; EntityHandle first = 0, prev = 0, handle; for (int i = 0; i < size[0]; ++i) { long count; if (!tokens.get_long_ints( 1, &count )) return MB_FAILURE; conn_idx.resize(count); conn_hdl.resize(count); if (!tokens.get_long_ints( count, &conn_idx[0])) return MB_FAILURE; for (long j = 0; j < count; ++j) conn_hdl[j] = first_vtx + conn_idx[j]; result = mdbImpl->create_element( MBPOLYGON, &conn_hdl[0], count, handle ); if (MB_SUCCESS != result) return result; if (prev +1 != handle) { if (first) { // true except for first iteration (first == 0) if (first < elem_list.back().front()) // only need new range if order would get mixed up elem_list.push_back( empty ); elem_list.back().insert( first, prev ); } first = handle; } prev = handle; } if (first) { // true unless no elements (size[0] == 0) if (first < elem_list.back().front()) // only need new range if order would get mixed up elem_list.push_back( empty ); elem_list.back().insert( first, prev ); } return MB_SUCCESS; }
ErrorCode ReadVtk::vtk_read_rectilinear_grid | ( | FileTokenizer & | tokens, |
Range & | vertex_list, | ||
std::vector< Range > & | elem_list | ||
) | [protected] |
Definition at line 524 of file ReadVtk.cpp.
{ int i, j, k; long dims[3]; const char* labels[] = { "X_COORDINATES", "Y_COORDINATES", "Z_COORDINATES" }; std::vector<double> coords[3]; ErrorCode result; if (!tokens.match_token( "DIMENSIONS" )|| !tokens.get_long_ints( 3, dims ) || !tokens.get_newline( )) return MB_FAILURE; if (dims[0] < 1 || dims[1] < 1 || dims[2] < 1) { readMeshIface->report_error( "Invalid dimension at line %d", tokens.line_number()); return MB_FAILURE; } for (i = 0; i < 3; i++) { long count; if (!tokens.match_token( labels[i] ) || !tokens.get_long_ints( 1, &count ) || !tokens.match_token( vtk_type_names )) return MB_FAILURE; if (count != dims[i]) { readMeshIface->report_error( "Coordinate count inconsistent with dimensions at line %d", tokens.line_number() ); return MB_FAILURE; } coords[i].resize(count); if (!tokens.get_doubles( count, &coords[i][0] )) return MB_FAILURE; } // Create vertices double *x, *y, *z; EntityHandle start_handle = 0; long num_verts = dims[0]*dims[1]*dims[2]; result = allocate_vertices( num_verts, start_handle, x, y, z ); if (MB_SUCCESS != result) return result; vertex_list.insert( start_handle, start_handle + num_verts - 1 ); // Calculate vertex coordinates for (k = 0; k < dims[2]; ++k) for (j = 0; j < dims[1]; ++j) for (i = 0; i < dims[0]; ++i) { *x = coords[0][i]; ++x; *y = coords[1][j]; ++y; *z = coords[2][k]; ++z; } return vtk_create_structured_elems( dims, start_handle, elem_list ); }
ErrorCode ReadVtk::vtk_read_scalar_attrib | ( | FileTokenizer & | tokens, |
std::vector< Range > & | entities, | ||
const char * | name | ||
) | [protected] |
Definition at line 1152 of file ReadVtk.cpp.
{ int type = tokens.match_token( vtk_type_names ); if (!type) return MB_FAILURE; long size; const char* tok = tokens.get_string( ); if (!tok) return MB_FAILURE; const char* end = 0; size = strtol( tok, (char**)&end, 0 ); if (*end) { size = 1; tokens.unget_token(); } // VTK spec says cannot be greater than 4--do we care? if (size < 1) //|| size > 4) { readMeshIface->report_error( "Scalar count out of range [1,4]" " at line %d", tokens.line_number()); return MB_FAILURE; } if (!tokens.match_token("LOOKUP_TABLE") || !tokens.match_token("default")) return MB_FAILURE; return vtk_read_tag_data( tokens, type, size, entities, name ); }
ErrorCode ReadVtk::vtk_read_structured_grid | ( | FileTokenizer & | tokens, |
Range & | vertex_list, | ||
std::vector< Range > & | elem_list | ||
) | [protected] |
Definition at line 481 of file ReadVtk.cpp.
{ long num_verts, dims[3]; ErrorCode result; if (!tokens.match_token( "DIMENSIONS" ) || !tokens.get_long_ints( 3, dims ) || !tokens.get_newline( ) ) return MB_FAILURE; if (dims[0] < 1 || dims[1] < 1 || dims[2] < 1) { readMeshIface->report_error( "Invalid dimension at line %d", tokens.line_number()); return MB_FAILURE; } if (!tokens.match_token( "POINTS" ) || !tokens.get_long_ints( 1, &num_verts ) || !tokens.match_token( vtk_type_names ) || !tokens.get_newline()) return MB_FAILURE; if (num_verts != (dims[0] * dims[1] * dims[2])) { readMeshIface->report_error( "Point count not consistent with dimensions at line %d", tokens.line_number() ); return MB_FAILURE; } // Create and read vertices EntityHandle start_handle = 0; result = read_vertices( tokens, num_verts, start_handle ); if (MB_SUCCESS != result) return result; vertex_list.insert( start_handle, start_handle + num_verts - 1 ); return vtk_create_structured_elems( dims, start_handle, elem_list ); }
ErrorCode ReadVtk::vtk_read_structured_points | ( | FileTokenizer & | tokens, |
Range & | vertex_list, | ||
std::vector< Range > & | elem_list | ||
) | [protected] |
Definition at line 428 of file ReadVtk.cpp.
{ long i, j, k; long dims[3]; double origin[3], space[3]; ErrorCode result; if (!tokens.match_token( "DIMENSIONS" ) || !tokens.get_long_ints( 3, dims ) || !tokens.get_newline( )) return MB_FAILURE; if (dims[0] < 1 || dims[1] < 1 || dims[2] < 1) { readMeshIface->report_error( "Invalid dimension at line %d", tokens.line_number()); return MB_FAILURE; } if (!tokens.match_token( "ORIGIN" ) || !tokens.get_doubles( 3, origin ) || !tokens.get_newline( )) return MB_FAILURE; const char* const spacing_names[] = { "SPACING", "ASPECT_RATIO", 0 }; if (!tokens.match_token( spacing_names ) || !tokens.get_doubles( 3, space ) || !tokens.get_newline()) return MB_FAILURE; // Create vertices double *x, *y, *z; EntityHandle start_handle = 0; long num_verts = dims[0]*dims[1]*dims[2]; result = allocate_vertices( num_verts, start_handle, x, y, z ); if (MB_SUCCESS != result) return result; vertex_list.insert( start_handle, start_handle + num_verts - 1 ); for (k = 0; k < dims[2]; ++k) for (j = 0; j < dims[1]; ++j) for (i = 0; i < dims[0]; ++i) { *x = origin[0] + i*space[0]; ++x; *y = origin[1] + j*space[1]; ++y; *z = origin[2] + k*space[2]; ++z; } return vtk_create_structured_elems( dims, start_handle, elem_list ); }
ErrorCode ReadVtk::vtk_read_tag_data | ( | FileTokenizer & | tokens, |
int | type, | ||
size_t | per_elem, | ||
std::vector< Range > & | entities, | ||
const char * | name | ||
) | [protected] |
Definition at line 1019 of file ReadVtk.cpp.
{ ErrorCode result; DataType mb_type; size_t size; if (type == 1) { mb_type = MB_TYPE_BIT; size = sizeof(bool); } else if (type >= 2 && type <= 9) { mb_type = MB_TYPE_INTEGER; size = sizeof(int); } else if (type == 10 || type == 11) { mb_type = MB_TYPE_DOUBLE; size = sizeof(double); } else if ( type == 12 ) { mb_type = MB_TYPE_INTEGER; size = 4; // could be 4 or 8, but we don't know. Hope it's 4 because MOAB doesn't support 64-bit ints. } else { return MB_FAILURE; } #ifdef MB_VTK_MATERIAL_SETS Modulator materialMap( this->mdbImpl, this->mPartitionTagName, mb_type, size, per_elem ); bool isMaterial = size * per_elem <= 4 && // must have int-sized values (ParallelComm requires it) ! this->mPartitionTagName.empty() && // must have a non-empty field name... ! strcmp( name, this->mPartitionTagName.c_str() ); // ... that matches our spec. #endif // MB_VTK_MATERIAL_SETS // get/create tag Tag handle; result = mdbImpl->tag_get_handle( name, per_elem, mb_type, handle, MB_TAG_DENSE|MB_TAG_CREAT ); if (MB_SUCCESS != result) { readMeshIface->report_error( "Tag name conflict for attribute \"%s\" at line %d\n", name, tokens.line_number() ); return result; } // Count number of entities long count = 0; std::vector<Range>::iterator iter; for (iter = entities.begin(); iter != entities.end(); ++iter) count += iter->size(); if (type == 1) { for (iter = entities.begin(); iter != entities.end(); ++iter) { bool *data = new bool[iter->size() * per_elem]; if (!tokens.get_booleans( per_elem*iter->size(), data )) { delete [] data; return MB_FAILURE; } bool* data_iter = data; Range::iterator ent_iter = iter->begin(); for ( ; ent_iter != iter->end(); ++ent_iter) { unsigned char bits = 0; for (unsigned j = 0; j < per_elem; ++j, ++data_iter) bits |= (unsigned char)(*data_iter << j); #ifdef MB_VTK_MATERIAL_SETS if ( isMaterial ) materialMap.add_entity( *ent_iter, &bits, 1 ); #endif // MB_VTK_MATERIAL_SETS result = mdbImpl->tag_set_data( handle, &*ent_iter, 1, &bits ); if (MB_SUCCESS != result) { delete [] data; return result; } } delete [] data; } } else if ((type >= 2 && type <= 9) || type == 12) { std::vector<int> data; for (iter = entities.begin(); iter != entities.end(); ++iter) { data.resize( iter->size() * per_elem ); if (!tokens.get_integers( iter->size() * per_elem, &data[0] )) return MB_FAILURE; #ifdef MB_VTK_MATERIAL_SETS if ( isMaterial ) materialMap.add_entities( *iter, (unsigned char*) &data[0], per_elem * size ); #endif // MB_VTK_MATERIAL_SETS result = mdbImpl->tag_set_data( handle, *iter, &data[0] ); if (MB_SUCCESS != result) return result; } } else if (type == 10 || type == 11) { std::vector<double> data; for (iter = entities.begin(); iter != entities.end(); ++iter) { data.resize( iter->size() * per_elem ); if (!tokens.get_doubles( iter->size() * per_elem, &data[0] )) return MB_FAILURE; #ifdef MB_VTK_MATERIAL_SETS if ( isMaterial ) materialMap.add_entities( *iter, (unsigned char*) &data[0], per_elem * size ); #endif // MB_VTK_MATERIAL_SETS result = mdbImpl->tag_set_data( handle, *iter, &data[0] ); if (MB_SUCCESS != result) return result; } } else { return MB_FAILURE; } return MB_SUCCESS; }
ErrorCode ReadVtk::vtk_read_tensor_attrib | ( | FileTokenizer & | tokens, |
std::vector< Range > & | entities, | ||
const char * | name | ||
) | [protected] |
Definition at line 1231 of file ReadVtk.cpp.
{ int type = tokens.match_token( vtk_type_names ); if (!type) return MB_FAILURE; return vtk_read_tag_data( tokens, type, 9, entities, name ); }
ErrorCode ReadVtk::vtk_read_texture_attrib | ( | FileTokenizer & | tokens, |
std::vector< Range > & | entities, | ||
const char * | name | ||
) | [protected] |
Definition at line 1212 of file ReadVtk.cpp.
{ int type, dim; if (!tokens.get_integers( 1, &dim ) || !(type = tokens.match_token( vtk_type_names ))) return MB_FAILURE; if (dim < 1 || dim > 3) { readMeshIface->report_error( "Invalid dimension (%d) at line %d.", dim, tokens.line_number() ); return MB_FAILURE; } return vtk_read_tag_data( tokens, type, dim, entities, name ); }
ErrorCode ReadVtk::vtk_read_unstructured_grid | ( | FileTokenizer & | tokens, |
Range & | vertex_list, | ||
std::vector< Range > & | elem_list | ||
) | [protected] |
Definition at line 707 of file ReadVtk.cpp.
{ ErrorCode result; long i, num_verts, num_elems[2]; EntityHandle tmp_conn_list[27]; // Poorly formatted VTK legacy format document seems to // lead many to think that a FIELD block can occur within // an UNSTRUCTURED_GRID dataset rather than as its own data // set. So allow for field data between other blocks of // data. const char* pts_str[] = { "FIELD", "POINTS", 0 }; while (1 == (i = tokens.match_token(pts_str))) { result = vtk_read_field(tokens); if (MB_SUCCESS != result) return result; } if (i != 2) return MB_FAILURE; if (!tokens.get_long_ints( 1, &num_verts ) || !tokens.match_token( vtk_type_names) || !tokens.get_newline( )) return MB_FAILURE; if (num_verts < 1) { readMeshIface->report_error( "Invalid point count at line %d", tokens.line_number()); return MB_FAILURE; } // Create vertices and read coordinates EntityHandle first_vertex = 0; result = read_vertices( tokens, num_verts, first_vertex ); if (MB_SUCCESS != result) return result; vertex_list.insert( first_vertex, first_vertex + num_verts - 1 ); const char* cell_str[] = { "FIELD", "CELLS", 0 }; while (1 == (i = tokens.match_token(cell_str))) { result = vtk_read_field(tokens); if (MB_SUCCESS != result) return result; } if (i != 2) return MB_FAILURE; if (!tokens.get_long_ints( 2, num_elems ) || !tokens.get_newline( )) return MB_FAILURE; // Read element connectivity for all elements std::vector<long> connectivity( num_elems[1] ); if (!tokens.get_long_ints( num_elems[1], &connectivity[0] )) return MB_FAILURE; if (!tokens.match_token( "CELL_TYPES" ) || !tokens.get_long_ints( 1, &num_elems[1]) || !tokens.get_newline( )) return MB_FAILURE; // Read element types std::vector<long> types( num_elems[0] ); if (!tokens.get_long_ints( num_elems[0], &types[0] )) return MB_FAILURE; // Create elements in blocks of the same type // It is important to preserve the order in // which the elements were read for later reading // attribute data. long id = 0; std::vector<long>::iterator conn_iter = connectivity.begin(); while (id < num_elems[0]) { unsigned vtk_type = types[id]; if (vtk_type >= VtkUtil::numVtkElemType) return MB_FAILURE; EntityType type = VtkUtil::vtkElemTypes[vtk_type].mb_type; if (type == MBMAXTYPE) { readMeshIface->report_error( "Unsupported VTK element type: %s (%d)\n", VtkUtil::vtkElemTypes[vtk_type].name, vtk_type ); return MB_FAILURE; } int num_vtx = *conn_iter; if (type != MBPOLYGON && num_vtx != (int)VtkUtil::vtkElemTypes[vtk_type].num_nodes) { readMeshIface->report_error( "Cell %ld is of type '%s' but has %u vertices.\n", id, VtkUtil::vtkElemTypes[vtk_type].name, num_vtx ); return MB_FAILURE; } // Find any subsequent elements of the same type std::vector<long>::iterator conn_iter2 = conn_iter + num_vtx + 1; long end_id = id + 1; while ( end_id < num_elems[0] && (unsigned)types[end_id] == vtk_type && *conn_iter2 == num_vtx) { ++end_id; conn_iter2 += num_vtx + 1; } // Allocate element block long num_elem = end_id - id; EntityHandle start_handle = 0; EntityHandle* conn_array; result = allocate_elements( num_elem, num_vtx, type, start_handle, conn_array, elem_list ); if (MB_SUCCESS != result) return result; EntityHandle *conn_sav = conn_array; // Store element connectivity for ( ; id < end_id; ++id) { if (conn_iter == connectivity.end()) { readMeshIface->report_error( "Connectivity data truncated at cell %ld\n", id ); return MB_FAILURE; } // make sure connectivity length is correct. if (*conn_iter != num_vtx) { readMeshIface->report_error( "Cell %ld is of type '%s' but has %u vertices.\n", id, VtkUtil::vtkElemTypes[vtk_type].name, num_vtx ); return MB_FAILURE; } ++conn_iter; for (i = 0; i < num_vtx; ++i, ++conn_iter) { if (conn_iter == connectivity.end()) { readMeshIface->report_error( "Connectivity data truncated at cell %ld\n", id ); return MB_FAILURE; } conn_array[i] = *conn_iter + first_vertex; } const unsigned* order = VtkUtil::vtkElemTypes[vtk_type].node_order; if ( order ) { assert( num_vtx * sizeof(EntityHandle) <= sizeof(tmp_conn_list) ); memcpy( tmp_conn_list, conn_array, num_vtx * sizeof(EntityHandle) ); for (int j = 0; j < num_vtx; ++j) conn_array[order[j]] = tmp_conn_list[j]; } conn_array += num_vtx; } // notify MOAB of the new elements result = readMeshIface->update_adjacencies(start_handle, num_elem, num_vtx, conn_sav); if (MB_SUCCESS != result) return result; } return MB_SUCCESS; }
ErrorCode ReadVtk::vtk_read_vector_attrib | ( | FileTokenizer & | tokens, |
std::vector< Range > & | entities, | ||
const char * | name | ||
) | [protected] |
Definition at line 1201 of file ReadVtk.cpp.
{ int type = tokens.match_token( vtk_type_names ); if (!type) return MB_FAILURE; return vtk_read_tag_data( tokens, type, 3, entities, name ); }
Meshset Handle for the mesh that is currently being read.
Definition at line 155 of file ReadVtk.hpp.
Interface* moab::ReadVtk::mdbImpl [private] |
interface instance
Definition at line 152 of file ReadVtk.hpp.
std::string moab::ReadVtk::mPartitionTagName [private] |
A field which, if present and having a single integer for storage, should be used to partition the mesh by range. Defaults to MATERIAL_SET_TAG_NAME.
Definition at line 158 of file ReadVtk.hpp.
ReadUtilIface* moab::ReadVtk::readMeshIface [private] |
Definition at line 147 of file ReadVtk.hpp.