Mesh Oriented datABase  (version 5.5.1)
An array-based unstructured mesh library
ReadRTT.cpp
Go to the documentation of this file.
1 /**
2  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
3  * storing and accessing finite element mesh data.
4  *
5  * Copyright 2004 Sandia Corporation. Under the terms of Contract
6  * DE-AC04-94AL85000 with Sandia Coroporation, the U.S. Government
7  * retains certain rights in this software.
8  *
9  * This library is free software; you can redistribute it and/or
10  * modify it under the terms of the GNU Lesser General Public
11  * License as published by the Free Software Foundation; either
12  * version 2.1 of the License, or (at your option) any later version.
13  *
14  */
15 
16 /**
17  * \class ReadRTT
18  * \brief ReadRTT based on ReadNASTRAN
19  *
20  * See:
21  *
22  * \author Andrew Davis
23  */
24 
25 #include "ReadRTT.hpp"
26 
27 #include <iostream>
28 #include <sstream>
29 #include <fstream>
30 #include <vector>
31 #include <cstdlib>
32 #include <map>
33 #include <cassert>
34 #include <cmath>
35 
36 #include "moab/Interface.hpp"
37 #include "moab/ReadUtilIface.hpp"
38 #include "Internals.hpp" // for MB_START_ID
39 #include "moab/Range.hpp"
40 #include "moab/FileOptions.hpp"
41 #include "FileTokenizer.hpp"
42 #include "MBTagConventions.hpp"
43 #include "moab/CN.hpp"
44 #include "moab/ErrorHandler.hpp"
45 #include "moab/GeomTopoTool.hpp"
46 
47 namespace moab
48 {
49 
51 {
52  return new ReadRTT( iface );
53 }
54 
55 // constructor
57  : MBI( impl ), geom_tag( 0 ), id_tag( 0 ), name_tag( 0 ), category_tag( 0 ), faceting_tol_tag( 0 )
58 {
59  assert( NULL != impl );
60  myGeomTool = new GeomTopoTool( impl );
62  assert( NULL != readMeshIface );
63 
64  // this section copied from ReadCGM initalisation
65  int negone = -1;
66  double zero = 0.;
67  ErrorCode rval;
69  &negone );MB_CHK_ERR_CONT(rval);
70  id_tag = MBI->globalId_tag();
74  rval =
76 }
77 
78 // destructor
80 {
81  if( readMeshIface )
82  {
84  readMeshIface = 0;
85  }
86 
87  delete myGeomTool;
88 }
89 
90 ErrorCode ReadRTT::read_tag_values( const char* /*file_name*/,
91  const char* /*tag_name*/,
92  const FileOptions& /*opts*/,
93  std::vector< int >& /*tag_values_out*/,
94  const SubsetList* /*subset_list*/ )
95 {
96  return MB_NOT_IMPLEMENTED;
97 }
98 
99 // load the file as called by the Interface function
100 ErrorCode ReadRTT::load_file( const char* filename,
101  const EntityHandle*,
102  const FileOptions&,
103  const ReaderIface::SubsetList* subset_list,
104  const Tag* /*file_id_tag*/ )
105 {
106  ErrorCode rval;
107 
108  // at this time there is no support for reading a subset of the file
109  if( subset_list )
110  {
111  std::cout << "Subset reading not supported for RTT meshes" << std::endl;
113  }
114 
115  // test to see if file exists
116  FILE* file = NULL;
117  file = fopen( filename, "r" );
118  if( file == NULL ) return MB_FILE_DOES_NOT_EXIST;
119  // otherwise close the file
120  fclose( file );
121 
122  // read the header
123  rval = ReadRTT::read_header( filename );
124  if( rval != MB_SUCCESS ) return rval;
125 
126  // read the side_flag data
127  std::vector< side > side_data;
128  rval = ReadRTT::read_sides( filename, side_data );
129  if( rval != MB_SUCCESS ) return rval;
130 
131  // read the cell data
132  std::vector< cell > cell_data;
133  rval = ReadRTT::read_cells( filename, cell_data );
134  if( rval != MB_SUCCESS ) return rval;
135 
136  // read the node data
137  std::vector< node > node_data;
138  rval = ReadRTT::read_nodes( filename, node_data );
139  if( rval != MB_SUCCESS ) return rval;
140 
141  // read the facet data
142  std::vector< facet > facet_data;
143  rval = ReadRTT::read_facets( filename, facet_data );
144  if( rval != MB_SUCCESS ) return rval;
145 
146  // read the tetrahedra data
147  std::vector< tet > tet_data;
148  rval = ReadRTT::read_tets( filename, tet_data );
149  if( rval != MB_SUCCESS ) return rval;
150 
151  // make the map of surface number in the rttmesh to the surface meshset
152  std::map< int, EntityHandle > surface_map; // corrsespondance of surface number to entity handle
153  rval = ReadRTT::generate_topology( side_data, cell_data, surface_map );
154  if( rval != MB_SUCCESS ) return rval;
155 
156  // generate the rest of the database, triangles to surface meshsets etc
157  rval = ReadRTT::build_moab( node_data, facet_data, tet_data, surface_map );
158  if( rval != MB_SUCCESS ) return rval;
159 
160  return MB_SUCCESS;
161 }
162 
163 /*
164  * builds the topology of the problem
165  */
166 ErrorCode ReadRTT::generate_topology( std::vector< side > side_data,
167  std::vector< cell > cell_data,
168  std::map< int, EntityHandle >& surface_map )
169 {
170 
171  ErrorCode rval;
172  std::vector< EntityHandle > entmap[4];
173  int num_ents[4]; // number of entities in each dimension
174 
175  const char geom_categories[][CATEGORY_TAG_SIZE] = { "Vertex\0", "Curve\0", "Surface\0", "Volume\0", "Group\0" };
176 
177  std::vector< int > surface_numbers; // the surface numbers in the problem
178 
179  // corresponds to number of cad like surfaces and cad like volumes
180  num_ents[2] = side_data.size();
181  num_ents[3] = cell_data.size();
182 
183  // loop over surfaces & volumes
184  for( int dim = 2; dim <= 3; dim++ )
185  {
186  for( int i = 0; i != num_ents[dim]; i++ )
187  {
188  EntityHandle handle;
189  // create a meshset for each entity surface/volume
190  rval = MBI->create_meshset( dim == 1 ? MESHSET_ORDERED : MESHSET_SET, handle );
191  // if failure
192  if( rval != MB_SUCCESS ) return rval;
193 
194  // collect the entity handles into an
195  entmap[dim].push_back( handle );
196 
197  // set the dimension tag
198  rval = MBI->tag_set_data( geom_tag, &handle, 1, &dim );
199  // if fail
200  if( MB_SUCCESS != rval ) return rval;
201  // if we are a surface
202  if( dim == 2 )
203  {
204  // tag the id onto the surface meshset
205  rval = MBI->tag_set_data( id_tag, &handle, 1, &side_data[i].id );
206  // inesert entity into the map
207  surface_map[side_data[i].id] = handle;
208  }
209  else
210  {
211  // otherwise we set the volume tag data, loop is only 2 & 3 dim
212  rval = MBI->tag_set_data( id_tag, &handle, 1, &cell_data[i].id );
213  }
214  // if fail
215  if( MB_SUCCESS != rval ) return rval;
216  // set the category tag
217  rval = MBI->tag_set_data( category_tag, &handle, 1, &geom_categories[dim] );
218  if( MB_SUCCESS != rval ) return rval;
219  }
220  }
221 
222  // generate parent child links
223  // best to loop over the surfaces and assign them to volumes, we can then assign facets to
224  // to each surface
225  generate_parent_child_links( num_ents, entmap, side_data, cell_data );
226 
227  // set the surface senses
228  set_surface_senses( num_ents, entmap, side_data, cell_data );
229 
230  // set the group data
231  rval = setup_group_data( entmap );
232 
233  return MB_SUCCESS;
234 }
235 
236 /*
237  * builds the moab representation of the mesh
238  */
239 ErrorCode ReadRTT::build_moab( std::vector< node > node_data,
240  std::vector< facet > facet_data,
241  std::vector< tet > tet_data,
242  std::map< int, EntityHandle > surface_map )
243 {
244 
245  ErrorCode rval; // reusable return value
246  EntityHandle file_set; // the file handle
247  // create the file set
248  rval = MBI->create_meshset( MESHSET_SET, file_set );
249  if( MB_SUCCESS != rval ) return rval;
250 
251  // create the vertices
252  EntityHandle handle;
253  std::vector< node >::iterator it; // iterate over the nodes
254  Range mb_coords; // range of coordinates
255  for( it = node_data.begin(); it != node_data.end(); ++it )
256  {
257  node tmp = *it;
258  double coords[3] = { tmp.x, tmp.y, tmp.z };
259  rval = MBI->create_vertex( coords, handle );
260  if( MB_SUCCESS != rval ) return rval;
261  mb_coords.insert( handle ); // inesert handle into the coordinate range
262  }
263 
264  // add verts to set
265  rval = MBI->add_entities( file_set, mb_coords );
266 
267  // create sense tag
268  Tag side_id_tag, surface_number_tag;
269  // int zero = 0;
270  rval = MBI->tag_get_handle( "SIDEID_TAG", 1, MB_TYPE_INTEGER, side_id_tag, MB_TAG_SPARSE | MB_TAG_CREAT );
271  rval =
272  MBI->tag_get_handle( "SURFACE_NUMBER", 1, MB_TYPE_INTEGER, surface_number_tag, MB_TAG_SPARSE | MB_TAG_CREAT );
273 
274  // create the facets
275  EntityHandle triangle;
276  std::vector< facet >::iterator it_f;
277  // range of triangles
278  Range mb_tris;
279  // loop over the facet data
280  for( it_f = facet_data.begin(); it_f != facet_data.end(); ++it_f )
281  {
282  facet tmp = *it_f;
283  // get the nodes for the triangle
284  EntityHandle tri_nodes[3] = { mb_coords[tmp.connectivity[0] - 1], mb_coords[tmp.connectivity[1] - 1],
285  mb_coords[tmp.connectivity[2] - 1] };
286  // create a triangle element
287  rval = MBI->create_element( MBTRI, tri_nodes, 3, triangle );
288  // tag in side id on the triangle
289  rval = MBI->tag_set_data( side_id_tag, &triangle, 1, &tmp.side_id );
290  // tag the surface number on the triangle
291  rval = MBI->tag_set_data( surface_number_tag, &triangle, 1, &tmp.surface_number );
292  // insert vertices and triangles into the appropriate surface meshset
293  EntityHandle meshset_handle = surface_map[tmp.surface_number];
294  // also set surface tag
295  rval = MBI->tag_set_data( side_id_tag, &meshset_handle, 1, &tmp.side_id );
296  rval = MBI->tag_set_data( surface_number_tag, &meshset_handle, 1, &tmp.surface_number );
297  // add vertices to the mesh
298  rval = MBI->add_entities( meshset_handle, &( *tri_nodes ), 3 );
299  // add triangles to the meshset
300  rval = MBI->add_entities( meshset_handle, &triangle, 1 );
301  // ineter triangles into large run
302  mb_tris.insert( triangle );
303  }
304  // add tris to set to fileset
305  rval = MBI->add_entities( file_set, mb_tris );
306 
307  // create material number tag
308  Tag mat_num_tag;
309  // int zero = 0;
310  rval = MBI->tag_get_handle( "MATERIAL_NUMBER", 1, MB_TYPE_INTEGER, mat_num_tag, MB_TAG_SPARSE | MB_TAG_CREAT );
311 
312  // create the tets
313  EntityHandle tetra; // handle for a specific tet
314  std::vector< tet >::iterator it_t;
315  Range mb_tets;
316  // loop over all tets
317  for( it_t = tet_data.begin(); it_t != tet_data.end(); ++it_t )
318  {
319  tet tmp = *it_t;
320  // get the handles for the tet
321  EntityHandle tet_nodes[4] = { mb_coords[tmp.connectivity[0] - 1], mb_coords[tmp.connectivity[1] - 1],
322  mb_coords[tmp.connectivity[2] - 1], mb_coords[tmp.connectivity[3] - 1] };
323  // create the tet
324  rval = MBI->create_element( MBTET, tet_nodes, 4, tetra );
325  int mat_number = tmp.material_number;
326  // tag the tet with the material number
327  rval = MBI->tag_set_data( mat_num_tag, &tetra, 1, &mat_number );
328  // set the tag data
329  mb_tets.insert( tetra );
330  }
331  // add tris to set
332  rval = MBI->add_entities( file_set, mb_tets );
333 
334  return MB_SUCCESS;
335 }
336 
337 /*
338  * read the header data from the filename pointed to
339  */
340 ErrorCode ReadRTT::read_header( const char* filename )
341 {
342  std::ifstream input_file( filename ); // filename for rtt file
343  // file ok?
344  if( !input_file.good() )
345  {
346  std::cout << "Problems reading file = " << filename << std::endl;
347  return MB_FAILURE;
348  }
349 
350  // if it works
351  std::string line;
352  moab::ErrorCode rval = MB_FAILURE;
353  if( input_file.is_open() )
354  {
355  while( std::getline( input_file, line ) )
356  {
357  if( line.compare( "header" ) == 0 )
358  {
359  rval = get_header_data( input_file );
360  }
361  }
362  input_file.close();
363  }
364  return rval;
365 }
366 
367 /*
368  * reads the side data from the filename pointed to
369  */
370 ErrorCode ReadRTT::read_sides( const char* filename, std::vector< side >& side_data )
371 {
372  std::string line; // the current line being read
373  std::ifstream input_file( filename ); // filestream for rttfile
374  // file ok?
375  if( !input_file.good() )
376  {
377  std::cout << "Problems reading file = " << filename << std::endl;
378  return MB_FAILURE;
379  }
380  // if it works
381  if( input_file.is_open() )
382  {
383  while( std::getline( input_file, line ) )
384  {
385  if( line.compare( " 2 FACES\0" ) == 0 )
386  {
387  // read lines until find end nodes
388  while( std::getline( input_file, line ) )
389  {
390  if( line.compare( "end_side_flags\0" ) == 0 ) break;
391  side data = ReadRTT::get_side_data( line );
392  side_data.push_back( data );
393  }
394  }
395  }
396  input_file.close();
397  }
398  if( side_data.size() == 0 ) return MB_FAILURE;
399  return MB_SUCCESS;
400 }
401 
402 /*
403  * reads the cell data from the filename pointed to
404  */
405 ErrorCode ReadRTT::read_cells( const char* filename, std::vector< cell >& cell_data )
406 {
407  std::string line; // the current line being read
408  std::ifstream input_file( filename ); // filestream for rttfile
409  // file ok?
410  if( !input_file.good() )
411  {
412  std::cout << "Problems reading file = " << filename << std::endl;
413  return MB_FAILURE;
414  }
415  // if it works
416  if( input_file.is_open() )
417  {
418  while( std::getline( input_file, line ) )
419  {
420  if( line.compare( " 1 REGIONS\0" ) == 0 )
421  {
422  // read lines until find end nodes
423  while( std::getline( input_file, line ) )
424  {
425  if( line.compare( "end_cell_flags\0" ) == 0 ) break;
426  cell data = ReadRTT::get_cell_data( line );
427  cell_data.push_back( data );
428  }
429  }
430  }
431  input_file.close();
432  }
433  if( cell_data.size() == 0 ) return MB_FAILURE;
434  return MB_SUCCESS;
435 }
436 
437 /*
438  * Reads the node data fromt the filename pointed to
439  */
440 ErrorCode ReadRTT::read_nodes( const char* filename, std::vector< node >& node_data )
441 {
442  std::string line; // the current line being read
443  std::ifstream input_file( filename ); // filestream for rttfile
444  // file ok?
445  if( !input_file.good() )
446  {
447  std::cout << "Problems reading file = " << filename << std::endl;
448  return MB_FAILURE;
449  }
450 
451  // if it works
452  if( input_file.is_open() )
453  {
454  while( std::getline( input_file, line ) )
455  {
456  if( line.compare( "nodes\0" ) == 0 )
457  {
458  // read lines until find end nodes
459  while( std::getline( input_file, line ) )
460  {
461  if( line.compare( "end_nodes\0" ) == 0 ) break;
462  node data = ReadRTT::get_node_data( line );
463  node_data.push_back( data );
464  }
465  }
466  }
467  input_file.close();
468  }
469  if( node_data.size() == 0 ) return MB_FAILURE;
470  return MB_SUCCESS;
471 }
472 
473 /*
474  * Reads the facet data fromt the filename pointed to
475  */
476 ErrorCode ReadRTT::read_facets( const char* filename, std::vector< facet >& facet_data )
477 {
478  std::string line; // the current line being read
479  std::ifstream input_file( filename ); // filestream for rttfile
480  // file ok?
481  if( !input_file.good() )
482  {
483  std::cout << "Problems reading file = " << filename << std::endl;
484  return MB_FAILURE;
485  }
486 
487  // if it works
488  if( input_file.is_open() )
489  {
490  while( std::getline( input_file, line ) )
491  {
492  if( line.compare( "sides\0" ) == 0 )
493  {
494  // read lines until find end nodes
495  while( std::getline( input_file, line ) )
496  {
497  if( line.compare( "end_sides\0" ) == 0 ) break;
498  facet data = ReadRTT::get_facet_data( line );
499  facet_data.push_back( data );
500  }
501  }
502  }
503  input_file.close();
504  }
505  if( facet_data.size() == 0 ) return MB_FAILURE;
506  return MB_SUCCESS;
507 }
508 
509 /*
510  * Reads the facet data fromt the filename pointed to
511  */
512 ErrorCode ReadRTT::read_tets( const char* filename, std::vector< tet >& tet_data )
513 {
514  std::string line; // the current line being read
515  std::ifstream input_file( filename ); // filestream for rttfile
516  // file ok?
517  if( !input_file.good() )
518  {
519  std::cout << "Problems reading file = " << filename << std::endl;
520  return MB_FAILURE;
521  }
522  // if it works
523  if( input_file.is_open() )
524  {
525  while( std::getline( input_file, line ) )
526  {
527  if( line.compare( "cells\0" ) == 0 )
528  {
529  // read lines until find end nodes
530  while( std::getline( input_file, line ) )
531  {
532  if( line.compare( "end_cells\0" ) == 0 ) break;
533  tet data = ReadRTT::get_tet_data( line );
534  tet_data.push_back( data );
535  }
536  }
537  }
538  input_file.close();
539  }
540  if( tet_data.size() == 0 ) return MB_FAILURE;
541  return MB_SUCCESS;
542 }
543 
544 /*
545  * given the open file handle read until we find
546  */
548 {
549  std::string line;
550  while( std::getline( input_file, line ) )
551  {
552 
553  // tokenize the line
554  std::istringstream iss( line );
555  std::vector< std::string > split_string;
556  do
557  {
558  std::string sub_string;
559  iss >> sub_string;
560  split_string.push_back( sub_string );
561  } while( iss );
562 
563  // if we find version
564  if( line.find( "version" ) != std::string::npos )
565  {
566  if( split_string[1].find( "v" ) != std::string::npos &&
567  split_string[0].find( "version" ) != std::string::npos )
568  {
570  }
571  }
572 
573  if( line.find( "title" ) != std::string::npos )
574  {
576  }
577  if( line.find( "date" ) != std::string::npos )
578  {
580  }
581  if( line.find( "end_header" ) != std::string::npos )
582  {
583  return MB_SUCCESS;
584  }
585  }
586 
587  // otherwise we never found the end_header keyword
588  return MB_FAILURE;
589 }
590 
591 /*
592  * given the string sidedata, get the id number, senses and names of the sides
593  */
594 ReadRTT::side ReadRTT::get_side_data( std::string sidedata )
595 {
596  side new_side;
597  std::vector< std::string > tokens;
598  tokens = ReadRTT::split_string( sidedata, ' ' );
599 
600  // set the side id
601  if( tokens.size() != 2 )
602  {
603  MB_SET_ERR_RET_VAL( "Error, too many tokens found from side_data", new_side );
604  }
605  // create the new side
606  new_side.id = std::atoi( tokens[0].c_str() );
607 
608  std::vector< std::string > cell_names = ReadRTT::split_string( tokens[1], '/' );
609  // get the boundary
610  boundary new_bnd = ReadRTT::split_name( cell_names[0] );
611  // set the surface sense and name
612  new_side.senses[0] = new_bnd.sense;
613  new_side.names[0] = new_bnd.name;
614  //
615  if( cell_names.size() > 1 )
616  {
617  boundary bnd = ReadRTT::split_name( cell_names[1] );
618  new_side.senses[1] = bnd.sense;
619  new_side.names[1] = bnd.name;
620  }
621  else
622  {
623  new_side.senses[1] = 0;
624  new_side.names[1] = "\0";
625  }
626 
627  return new_side;
628 }
629 
630 /*
631  * given the string celldata, get the id number and name of each cell
632  */
633 ReadRTT::cell ReadRTT::get_cell_data( std::string celldata )
634 {
635  cell new_cell;
636  std::vector< std::string > tokens;
637  tokens = ReadRTT::split_string( celldata, ' ' );
638 
639  // set the side id
640  if( tokens.size() != 2 )
641  {
642  MB_SET_ERR_RET_VAL( "Error, too many tokens found from cell_data", new_cell );
643  }
644  // create the new side
645  new_cell.id = std::atoi( tokens[0].c_str() );
646  new_cell.name = tokens[1];
647 
648  return new_cell;
649 }
650 
651 /*
652  * given the string nodedata, get the id number and coordinates of the node
653  */
654 ReadRTT::node ReadRTT::get_node_data( std::string nodedata )
655 {
656  node new_node;
657  std::vector< std::string > tokens;
658  tokens = ReadRTT::split_string( nodedata, ' ' );
659 
660  // set the side id
661  if( tokens.size() != 5 )
662  {
663  MB_SET_ERR_RET_VAL( "Error, too many tokens found from get_node_data", new_node );
664  }
665  new_node.id = std::atoi( tokens[0].c_str() );
666  new_node.x = std::atof( tokens[1].c_str() );
667  new_node.y = std::atof( tokens[2].c_str() );
668  new_node.z = std::atof( tokens[3].c_str() );
669  return new_node;
670 }
671 
672 /*
673  * given the string nodedata, get the id number, connectivity and sense data
674  */
675 ReadRTT::facet ReadRTT::get_facet_data( std::string facetdata )
676 {
677  facet new_facet;
678  std::vector< std::string > tokens;
679  tokens = ReadRTT::split_string( facetdata, ' ' );
680 
681  // set the side id
682  if( tokens.size() != 7 )
683  {
684  MB_SET_ERR_RET_VAL( "Error, too many tokens found from get_facet_data", new_facet );
685  }
686 
687  new_facet.id = std::atoi( tokens[0].c_str() );
688  // branch on the rtt version number
689  if( header_data.version == "v1.0.0" )
690  {
691  new_facet.connectivity[0] = std::atoi( tokens[1].c_str() );
692  new_facet.connectivity[1] = std::atoi( tokens[2].c_str() );
693  new_facet.connectivity[2] = std::atoi( tokens[3].c_str() );
694  new_facet.side_id = std::atoi( tokens[4].c_str() );
695  new_facet.surface_number = std::atoi( tokens[5].c_str() );
696  }
697  else if( header_data.version == "v1.0.1" )
698  {
699  new_facet.connectivity[0] = std::atoi( tokens[2].c_str() );
700  new_facet.connectivity[1] = std::atoi( tokens[3].c_str() );
701  new_facet.connectivity[2] = std::atoi( tokens[4].c_str() );
702  new_facet.side_id = std::atoi( tokens[5].c_str() );
703  new_facet.surface_number = std::atoi( tokens[6].c_str() );
704  }
705  else
706  {
707  MB_SET_ERR_RET_VAL( "Error, version number not understood", new_facet );
708  }
709 
710  return new_facet;
711 }
712 
713 /*
714  * given the string tetdata, get the id number, connectivity and mat num of the tet
715  */
716 ReadRTT::tet ReadRTT::get_tet_data( std::string tetdata )
717 {
718  tet new_tet;
719  std::vector< std::string > tokens;
720  tokens = ReadRTT::split_string( tetdata, ' ' );
721 
722  // set the side id
723  if( tokens.size() != 7 )
724  {
725  MB_SET_ERR_RET_VAL( "Error, too many tokens found from get_tet_data", new_tet );
726  }
727  new_tet.id = std::atoi( tokens[0].c_str() );
728  // branch on the version number
729  if( header_data.version == "v1.0.0" )
730  {
731  new_tet.connectivity[0] = std::atoi( tokens[1].c_str() );
732  new_tet.connectivity[1] = std::atoi( tokens[2].c_str() );
733  new_tet.connectivity[2] = std::atoi( tokens[3].c_str() );
734  new_tet.connectivity[3] = std::atoi( tokens[4].c_str() );
735  new_tet.material_number = std::atoi( tokens[5].c_str() );
736  }
737  else if( header_data.version == "v1.0.1" )
738  {
739  new_tet.connectivity[0] = std::atoi( tokens[2].c_str() );
740  new_tet.connectivity[1] = std::atoi( tokens[3].c_str() );
741  new_tet.connectivity[2] = std::atoi( tokens[4].c_str() );
742  new_tet.connectivity[3] = std::atoi( tokens[5].c_str() );
743  new_tet.material_number = std::atoi( tokens[6].c_str() );
744  }
745  else
746  {
747  MB_SET_ERR_RET_VAL( "Error, version number not supported", new_tet );
748  }
749 
750  return new_tet;
751 }
752 
753 /*
754  * splits string into sense and name, to later facilitate the building
755  * of sense data, strips off the tailing @ if it exists
756  */
757 ReadRTT::boundary ReadRTT::split_name( std::string atilla_cellname )
758 {
759  boundary new_boundary;
760  // default initialisation
761  new_boundary.sense = 0;
762  new_boundary.name = "\0";
763  // +ve sense
764  if( atilla_cellname.find( "+" ) != std::string::npos )
765  {
766  new_boundary.sense = 1;
767  // look for the @# we do not want it
768  std::size_t found = atilla_cellname.find( "@" );
769  if( found != std::string::npos )
770  new_boundary.name = atilla_cellname.substr( 3, found );
771  else
772  new_boundary.name = atilla_cellname.substr( 3, atilla_cellname.length() );
773  }
774  else if( atilla_cellname.find( "-" ) != std::string::npos )
775  {
776  // negative sense
777  new_boundary.sense = -1;
778  new_boundary.name = atilla_cellname.substr( 3, atilla_cellname.length() );
779  }
780  return new_boundary;
781 }
782 
783 /*
784  * splits a string in a vector of strings split by spaces
785  */
786 std::vector< std::string > ReadRTT::split_string( std::string string_to_split, char split_char )
787 {
788  std::istringstream ss( string_to_split );
789  std::vector< std::string > tokens;
790  while( !ss.eof() )
791  {
792  std::string x; // here's a nice, empty string
793  std::getline( ss, x, split_char ); // try to read the next field into it
794  tokens.push_back( x );
795  }
796 
797  // remove empty tokens
798  std::vector< std::string >::iterator it;
799  for( it = tokens.begin(); it != tokens.end(); )
800  {
801  std::string string = *it;
802  if( string.compare( "\0" ) == 0 )
803  it = tokens.erase( it );
804  else
805  ++it;
806  }
807  return tokens;
808 }
809 
810 /*
811  * Generate the parent-child links bwetween the cell and surface meshsets
812  */
814  std::vector< EntityHandle > entity_map[4],
815  std::vector< side > side_data,
816  std::vector< cell > cell_data )
817 {
818  ErrorCode rval; // return value
819  // loop over the number of surfaces
820  for( int i = 0; i < num_ents[2]; i++ )
821  {
822  // get the surface handle
823  EntityHandle surf_handle = entity_map[2][i];
824  // there are volumes that share this face
825  for( unsigned int shared = 0; shared <= 1; shared++ )
826  {
827  std::string parent_name = side_data[i].names[shared];
828  // find the @ sign
829  unsigned pos = parent_name.find( "@" );
830  parent_name = parent_name.substr( 0, pos );
831 
832  // loop over tets looking for matching name
833  for( int j = 0; j < num_ents[3]; j++ )
834  {
835  // if match found
836  if( cell_data[j].name.compare( parent_name ) == 0 )
837  {
838  EntityHandle cell_handle = entity_map[3][j];
839  // parent
840  rval = MBI->add_parent_child( cell_handle, surf_handle );
841  if( rval != MB_SUCCESS )
842  {
843  std::cerr << "Failed to add parent child relationship" << std::endl;
844  }
845  }
846  }
847  }
848  }
849  return;
850 }
851 
852 /*
853  * sets the sense of the surfaces wrt to volumes using geom topo tool
854  */
855 void ReadRTT::set_surface_senses( int num_ents[4],
856  std::vector< EntityHandle > entity_map[4],
857  std::vector< side > side_data,
858  std::vector< cell > cell_data )
859 {
860 
861  ErrorCode rval; // return value
862  // loop over the number of surfaces
863  for( int i = 0; i < num_ents[2]; i++ )
864  {
865  EntityHandle surf_handle = entity_map[2][i];
866  // there are 2 volumes that share this face
867  for( unsigned int shared = 0; shared <= 1; shared++ )
868  {
869  std::string parent_name = side_data[i].names[shared];
870  unsigned pos = parent_name.find( "@" );
871  parent_name = parent_name.substr( 0, pos );
872  // loop over tets looking for matching name
873  for( int j = 0; j < num_ents[3]; j++ )
874  {
875  // if match found
876  if( cell_data[j].name.compare( parent_name ) == 0 )
877  {
878  EntityHandle cell_handle = entity_map[3][j];
879  // in rtt mesh +represents the inside and -represents outside
880  // in moab reverse is outside and forward is inside
881  if( side_data[i].senses[shared] == 1 )
882  rval = myGeomTool->set_sense( surf_handle, cell_handle, SENSE_FORWARD );
883  else if( side_data[i].senses[shared] == -1 )
884  rval = myGeomTool->set_sense( surf_handle, cell_handle, SENSE_REVERSE );
885  else
886  rval = myGeomTool->set_sense( surf_handle, 0, SENSE_REVERSE );
887 
888  if( rval != MB_SUCCESS )
889  {
890  std::cerr << "Failed to set sense appropriately" << std::endl;
891  }
892  }
893  }
894  }
895  }
896  return;
897 }
898 
899 /*
900  * Add all entities that are to be part of the graveyard
901  */
902 ErrorCode ReadRTT::setup_group_data( std::vector< EntityHandle > entity_map[4] )
903 {
904  ErrorCode rval; // error codes
905  EntityHandle handle;
906  handle = create_group( "graveyard_comp", 1 );
907 
908  // add any volume to group graveyard, it is ignored by dag
909  EntityHandle vol_handle = entity_map[3][0];
910  rval = MBI->add_entities( handle, &vol_handle, 1 );
911  return rval;
912 }
913 
914 /*
915  * create a new group of with the name group name and id
916  */
917 EntityHandle ReadRTT::create_group( std::string group_name, int id )
918 {
919  ErrorCode rval;
920  // category tags
921  const char geom_categories[][CATEGORY_TAG_SIZE] = { "Vertex\0", "Curve\0", "Surface\0", "Volume\0", "Group\0" };
922 
923  EntityHandle handle;
924  rval = MBI->create_meshset( MESHSET_SET, handle );
925  if( MB_SUCCESS != rval ) return rval;
926 
927  rval = MBI->tag_set_data( name_tag, &handle, 1, group_name.c_str() );
928  if( MB_SUCCESS != rval ) return MB_FAILURE;
929 
930  rval = MBI->tag_set_data( id_tag, &handle, 1, &id );
931  if( MB_SUCCESS != rval ) return MB_FAILURE;
932 
933  rval = MBI->tag_set_data( category_tag, &handle, 1, &geom_categories[4] );
934  if( MB_SUCCESS != rval ) return MB_FAILURE;
935 
936  return handle;
937 }
938 
939 } // namespace moab